Lernziele in R

Explorative Datenanalyse
R Grafiken mit base
kleine Random-Effects-Metaanalyse mit Odds Ratios

Seit 2016 führen wir in der Vorlesung im 1. Semester Methodenlehre eine kleine Datenerhebung mit Studierenden durch. Eine der Fragen bezieht sich auf die Ernährungsform:

Sind Sie Vegetarier? (Ja / Nein)

Wir wollen die Entwicklung der vegetarischen Ernährung über die letzten Jahre betrachten. Für das Jahr 2021 liegen keine Daten vor, da wir im Rahmen eines Forschungsprojekts eine andere Umfrage durchgeführt haben.

1. Anteil an Vegetarier:innen

# Datenimport
library(readr)
Veg <- read_delim("https://www-user.tu-chemnitz.de/~burma/blog_data/Vegetarier.csv", 
 delim = ";", escape_double = FALSE, locale = locale(encoding = "ISO-8859-1"), 
 trim_ws = TRUE)

# Daten-Split nach Jahren
Veg_split1 <- split(Veg, f = Veg$Zeitstempel)
Veg_ja <- unlist(
  lapply(Veg_split1, function(x) (table(x$Veg)/sum(table(x$Veg)))[1])
)
N <- unlist(lapply(Veg_split1, function(x) sum(table(x$Veg))))   

# Konfidenzintervalle 
library(Hmisc)
CI <- matrix(NA, nrow = length(N), ncol = 2)
for (i in 1:7){
  CI[i, ] <- binconf(x = Veg_ja[i] * N[[i]], n = N[[i]], alpha = .05)[2:3]
  }
Veg_table <- cbind(Veg_ja, N, CI)
colnames(Veg_table)[3:4] <- c("lower", "upper")

			95 % - KI
	Anteil Vegetarier	N	UG	OG
2016	0.08	89	0.04	0.15
2017	0.15	110	0.10	0.23
2018	0.31	127	0.23	0.39
2019	0.28	110	0.21	0.37
2020	0.40	87	0.31	0.51
2022	0.38	55	0.27	0.51
2023	0.40	90	0.30	0.50

# Plot
plot(1:7, Veg_table[,1], ylim=c(0,1), ylab= "Anteil der Vegetarier im 1. Semester",
     xlab = "Jahr", xaxt = "n")
axis(1, at = 1:7, labels = c(2016:2020, 2022:2023))
for(i in 1:7){
  arrows(x0 = i, CI[i,1], i, CI[i,2], angle = 90, code = 3, length = .1)
}

Offensichtlich hat der Anteil an Vegetarier:innen in den letzten Jahren stark zugenommen auf nun mehr ca. 40 % (stabil in den letzten 3 Jahren). Laut einer repräsentativen Forsa-Umfrage von 2023 entspricht der Vegetarier:innen-Anteil in Deutschland 12 %. Unsere Stichprobe der Erstsemester (Studiengänge Psychologie und Sensorik und kognitive Psychologie) liegt deutlich darüber. Dieser Befund deckt sich mit früheren Befragungen (Mensink, Barbosa & Brettschneider 2016). Das Interessante ist für uns die Veränderung innerhalb weniger Jahre!

Ursachen den Anstieg vegetarischer Ernährung

Hopwood et al. (2020) konnten 3 von einander unabhängige Motive für Vegetarismus empirisch abgrenzen:

Gesundheitsbewusstsein
Umwelt- und Klimaschutz
Tierschutz und Tierrechte

In meinen eigenen Lehrveranstaltungen und Diskussionen mit Studierenden bemerke ich großes Bewusstsein für Umwelt- und Klimaschutz. Das Aufkommen der “Fridays for Future”-Bewegung seit 2019 legt den Schluss nahe, dass Klimaschutz eine wesentliche Motivation für das veränderte Ernährungsverhalten sein könnte. Für das kommende Jahr werden wir die Motivation der Studierenden erfragen, um diese Annahme empirisch zu prüfen.

2. Geschlechtereffekte

In der oben zitierten Umfrage “… ist der Anteil der Frauen, die sich vegetarisch ernähren, doppelt so hoch wie bei den Männern.” Auch in unseren Kohorten ernähren sich mehr Frauen vegetarisch. Wir können das oben generierte Script einfach erneut nutzen und erweitern in split das Argument f.

# Daten-Split nach Jahren und Geschlecht
Veg_split2 <- split(Veg, f = list(Veg$Geschlecht, Veg$Zeitstempel))
Veg_ja <- unlist(
  lapply(Veg_split2, function(x) (table(x$Veg)/sum(table(x$Veg)))[1])
)
N <- unlist(lapply(Veg_split2, function(x) sum(table(x$Veg)))   )

N <- unlist(lapply(Veg_split2, function(x) sum(table(x$Veg)))   )

# Konfidenzintervalle 
CI <- matrix(NA, nrow = length(N), ncol = 2)
for (i in 1:length(N)){
    ifelse(N[i] == 0, CI[i, ] <- c(NA, NA), 
         CI[i, ] <- binconf(x = Veg_ja[i] * N[[i]], n = N[[i]], alpha = .05)[2:3]
    )
  }
Veg_table_sex <- cbind(Veg_ja, N, CI)
colnames(Veg_table_sex)[3:4] <- c("lower", "upper")

Nun ist noch eine kleine Nachbereitung der Tabelle notwendig. Aufgrund der geringen Fallzahlen entfernen wir die Personen aus dem Datensatz, die sich als “divers” bezeichnen. Diese Schätzungen sind wenig vertrauenswürdig. Außerdem ersetzen wir die Zeile 2 durch 0 % Vegetarier (anstelle von 100 % Nicht-Vegetarier.) Nun zeichnen wir auch hier die Anteile mit den 95 %- Konfidenzintervallen.

# Nun ist noch eine Nachbereitung der Tabelle notwendig:
Veg_table_sex <- Veg_table_sex[-seq(1,21,3), ]
Veg_table_sex[1, c(1, 3, 4)] <- c(binconf(0,24))
rownames(Veg_table_sex)[1] <- "männlich.2016.ja"

# Grafik
xp <- c(1:21)[-seq(3,21,3)]
plot(xp,
  Veg_table_sex[,1], ylim = c(0,1), 
  ylab = "Anteil der Vegetarier im 1. Semester",
  xlab = "Jahr", xaxt = "n")
axis(1, at = xp[seq(1, 14, 2)] + .5, 
     labels = c(2016:2020, 2022:2023))
for (i in 1:14){
  arrows(x0 = xp[i], Veg_table_sex[i,3], 
         xp[i], Veg_table_sex[i,4], 
         angle = 90, code = 3, length = .1, 
         col = rep(c("blue", "darkred"), 7)[i]
  )
legend("topleft", legend = c("M", "W"), 
       fill = c("blue", "darkred"))
}

Zunächst lasst sich der bereits bekannte Anstieg im Verlauf der Zeit wieder finden. Außerdem lässt sich ein tendenziell höherer Vegetarieranteil bei den Frauen beobachten. Wir berechnen das mittlere \(OR\) aus den Daten. Dabei lassen wir 2016 (\(p_{Veg,männlich} = 0\)) und 2022 (\(N_{Veg,männlich} = 2\)) außen vor. Da es sich offensichtlich um unterschiedliche Populationseffekte handelt, ist eine Random-Effects-Metaanalyse das Mittel der Wahl.

library(metafor)
ORs <- NULL
for (i in seq(1,14,2)){
  ORs[[i]] <- unlist(escalc(
    measure = "OR",
    # Veg. männlich
    ai = Veg_table_sex[i, "Veg_ja"] * Veg_table_sex[i, "N"],
    # Nicht Veg. männlich
    bi =  (1 - Veg_table_sex[i, "Veg_ja"]) * Veg_table_sex[i, "N"],
    # Veg. weiblich
    ci = Veg_table_sex[i + 1 , "Veg_ja"] * Veg_table_sex[i + 1, "N"],
    # Nicht-Veg. weiblich
    di = (1 - Veg_table_sex[i + 1 , "Veg_ja"]) * 
      Veg_table_sex[i + 1, "N"]
    ))[c("yi", "vi")]
}
ORs <- do.call(rbind, ORs)[c(2,3,4,5,7), ]   
meta_ORs <- rma(yi = ORs[, 1], vi = ORs[ , 2], measure = "OR")
meta_ORs


Random-Effects Model (k = 5; tau^2 estimator: REML)

tau^2 (estimated amount of total heterogeneity): 0.2804 (SE = 0.4697)
tau (square root of estimated tau^2 value):      0.5295
I^2 (total heterogeneity / total variability):   42.56%
H^2 (total variability / sampling variability):  1.74

Test for Heterogeneity:
Q(df = 4) = 6.7249, p-val = 0.1512

Model Results:

estimate      se     zval    pval    ci.lb    ci.ub    
 -0.9177  0.3656  -2.5099  0.0121  -1.6344  -0.2011  * 

---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Die Heterogenität der Effekte ist nicht überraschend. Zum einen ist der Anteil an männlichen Befragten generell kleiner und somit unsicherer, zum anderen haben wir bereits die Veränderung des Ernährungsverhaltens in den letzten Jahren ausgemacht.

predict(meta_ORs, transf = exp, digits = 2)


 pred ci.lb ci.ub pi.lb pi.ub 
 0.40  0.20  0.82  0.11  1.41

Schließlich zeigt sich ein mittleres Odds Ratio von \(OR = 0.4\). Die Chance, unter den Vegetarier:innen eine männliche Person anzutreffen, ist lediglich 0.4-mal so hoch wie bei den Frauen.

Literatur

Mensink, G., Barbosa, C. L., & Brettschneider, A. K. (2016). Verbreitung der vegetarischen Ernährungsweise in Deutschland.

Hopwood CJ, Bleidorn W, Schwaba T, Chen S (2020) Health, environmental, and animal rights motives for vegetarian eating. PLoS ONE 15(4): e0230609. https://doi.org/10.1371/journal.pone.0230609