Lernziele in R

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

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

Sind Sie Vegetarier? (Ja / Nein)

Wir wollen die Entwicklung des Anteils vegetarisch lebender Studierender ü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 nunmehr rund 40 % (stabil in den letzten drei 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). Für uns ist vor allem die Veränderung innerhalb weniger Jahre interessant.

Ursachen des Anstiegs vegetarischer Ernährung

Hopwood et al. (2020) konnten drei voneinander unabhängige Motive für Vegetarismus empirisch abgrenzen:

Gesundheitsbewusstsein
Umwelt- und Klimaschutz
Tierschutz und Tierrechte

In meinen eigenen Lehrveranstaltungen und Diskussionen mit Studierenden beobachte ich ein 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 heißt es, dass der Anteil der Frauen, die sich vegetarisch ernähren, etwa doppelt so hoch ist 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)))   )

# 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 sehr geringen Fallzahlen für die Kategorie divers können wir hier keine belastbaren Schätzungen vornehmen und schließen diese Kategorie in den folgenden Analysen aus. Dies ist eine rein statistische Entscheidung und sagt nichts über die Relevanz dieser Gruppe aus. Außerdem setzen wir die Zeile 2 auf 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 lässt sich der bereits bekannte Anstieg im Verlauf der Zeit wiederfinden. 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