21 Verallgemeinerte lineare Modelle
21.1 Folien
21.2 Daten zur heutigen Sitzung
21.3 Code und Ausgaben aus der Vorlesung
Laden der relevanten Pakete
library(MASS) # Für Negativ-binomiale Regression
library(marginaleffects) # Modellvorhersagen und Vergleiche
library(performance) # Modellbeurteilung (hier vor allem R2 und Dispersion)
library(parameters) # Modellschätzungen extrahieren und darstellen
library(report) # Einfaches Erstellen von statistischen Berichten
library(tidyverse) # Datenmanagement und Visualisierung: https://www.tidyverse.org/── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr 1.2.0 ✔ readr 2.1.6
✔ forcats 1.0.1 ✔ stringr 1.6.0
✔ ggplot2 4.0.2 ✔ tibble 3.3.1
✔ lubridate 1.9.4 ✔ tidyr 1.3.2
✔ purrr 1.2.1
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
✖ dplyr::select() masks MASS::select()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorstheme_set(theme_classic(base_size = 12)) # Voreinstellung für PlotsLesen und Aufbereiten des Datensatz von Van Erkel & Van Aelst
d <- haven::read_stata(here::here("data/Vanerkel_Vanaelst_2021.dta")) |>
rename(
Political_knowledge = PK,
Personalized_news = personalized_news,
Radio = News_channels_w4_1,
Television = News_channels_w4_2,
Newspapers = News_channels_w4_3,
Online_news_sites = News_channels_w4_4,
Twitter = News_channels_w4_5,
Facebook = News_channels_w4_6,
) |>
mutate(
Gender = as_factor(Gender),
Education = as_factor(Education),
Age10 = Age / 10
)Beschreibung einzelen Fragen
d |>
select(starts_with("PK")) |>
mutate(across(everything(), ~ factor(., labels = c("incorrect", "correct")))) |>
report_sample()# Descriptive Statistics
Variable | Summary
--------------------------
PK1 [correct], % | 83.1
PK2 [correct], % | 68.7
PK3 [correct], % | 40.7
PK4 [correct], % | 87.9
PK5 [correct], % | 62.8
PK6 [correct], % | 24.1Beschreibung Wissensindex
d |>
ggplot(aes(Political_knowledge)) +
geom_bar()
Daten aus Faehnrich et al. laden
d2 <- read_rds(here::here("data/Faehnrich_2020.rds"))Auszug aus dem Datensatz von Faehnrich et al. - Zaehlvariablen
d2 |>
select(id, likes_count, comments_count, shares_count) |>
slice_sample(n = 10)# A tibble: 10 × 4
id likes_count comments_count shares_count
<chr> <dbl> <dbl> <dbl>
1 102471088936_10151782113633937 32 0 2
2 103256838688_10152076126808689 248 4 22
3 21489041474_10151617226346475 70 0 3
4 16686610106_10154986047745107 513 22 29
5 140105122708_10152984518697709 336 6 98
6 245640871929_10152948041316930 775 17 91
7 49702985881_10151464322755882 34 0 9
8 374809010319_10152312330145320 86 2 10
9 10111634660_706131796084572 114 1 15
10 20835777216_10152626212717217 133 0 22Plot Lineare Regression mit PK3 - Simuliertes Beispiel
set.seed(1)
d_fake <- d |>
mutate(PK3_fake = rbinom(n = n(), size = 1, prob = plogis(-10 + 2 * Age10)))
plt_bin_fake <- d_fake |>
ggplot(aes(Age10, PK3_fake)) +
geom_point(position = position_jitter(width = 0, height = 0.03, seed = 1), shape = 1) +
labs(x = "Alter in Jahrzehnten (Age / 10)", y = "Frage 3 (0 = falsch, 1 = richtig)") +
scale_y_continuous(limits = c(-0.2, 1.2), breaks = (0:4) / 4)
plt_bin_fake +
geom_smooth(method = "lm", linewidth = 2, se = F)`geom_smooth()` using formula = 'y ~ x'Warning: Removed 6 rows containing missing values or values outside the scale range
(`geom_smooth()`).
Plot Logistische Regression mit PK3 - Simuliertes Beispiel
plt_bin_fake +
geom_smooth(
method = "glm", method.args = list(family = "binomial"),
se = F, linewidth = 2
)`geom_smooth()` using formula = 'y ~ x'
Plot Lineare Regression mit PK3
plt_bin <- d |>
ggplot(aes(Age10, PK3)) +
geom_point(position = position_jitter(width = 0, height = 0.03, seed = 1), shape = 1) +
labs(x = "Alter in Jahrzehnten (Age / 10)", y = "Frage 3 (0 = falsch, 1 = richtig)") +
scale_y_continuous(limits = c(-0.2, 1.2), breaks = (0:4) / 4)
plt_bin +
geom_smooth(method = "lm", linewidth = 2, se = F)`geom_smooth()` using formula = 'y ~ x'
Plot Logistische Regression mit PK3
plt_bin +
geom_smooth(
method = "glm", method.args = list(family = "binomial"),
se = F, linewidth = 2
)`geom_smooth()` using formula = 'y ~ x'
Null-Modell Logistische Regression mit PK3
m0 <- glm(
PK3 ~ 1, # Null-Modell nur mit Intercept (Konstante)
family = binomial(link = "logit"), # aV ist binär, logit-Link = Logistische Regression
data = d
)Null-Modell Logistische Regression mit PK3 Tabelle
m0 |>
report_table(metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | p | Fit
-----------------------------------------------------------------
(Intercept) | -0.38 | [-0.50, -0.25] | -5.84 | < .001 |
| | | | |
Tjur's R2 | | | | | 0Invers logit Formel
1 / (1 + exp(0.38))[1] 0.4061269Invers logit Funktion
plogis(-0.38)[1] 0.4061269PK 3 deskriptiv
d |>
select(PK3) |>
mutate(across(everything(), ~ factor(., labels = c("incorrect", "correct")))) |>
report_sample()# Descriptive Statistics
Variable | Summary
--------------------------
PK3 [correct], % | 40.7Bivariate Logistische Regression mit PK3 und Alter
d <- d |>
mutate(
Age10_c = Age10 - mean(Age10)
)
m1 <- glm(
PK3 ~ Age10_c, # Regressionsgleichung
family = binomial(link = "logit"), # aV ist binär, logit-Link = Logistische Regression
data = d
)Bivariate Logistische Regression mit PK3 und Alter Tabelle
m1 |>
report_table(metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | p | Fit
------------------------------------------------------------------
(Intercept) | -0.39 | [-0.52, -0.26] | -5.97 | < .001 |
Age10 c | 0.28 | [ 0.18, 0.38] | 5.60 | < .001 |
| | | | |
Tjur's R2 | | | | | 0.03Bivariate Logistische Regression mit PK3 und Alter Baseline Wahrscheinlichkeit
plogis(-0.39)[1] 0.4037173Bivariate Logistische Regression mit PK3 und Alter Average counterfactual comparison M+1
m1 |>
avg_comparisons(
variables = list(Age10_c = 1), # Welcher Vergleich? Hier: Mittelwert und M + 1
type = "response" # In welcher Metrik? Hier: Wahrscheinlichkeit
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.066 0.0112 5.91 <0.001 28.2 0.0442 0.0879
Term: Age10_c
Type: response
Comparison: +1Bivariate Logistische Regression mit PK3 und Alter Average counterfactual comparison M-1
m1 |>
avg_comparisons(
variables = list(Age10_c = -1), # Welcher Vergleich? Hier: Mittelwert und M - 1
type = "response" # In welcher Metrik? Hier: Wahrscheinlichkeit
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
-0.063 0.0103 -6.14 <0.001 30.2 -0.0831 -0.0429
Term: Age10_c
Type: response
Comparison: +-1Bivariate Logistische Regression mit PK3 und Alter Average counterfactual comparison Junge
m1 |>
avg_comparisons(
variables = list(Age10_c = c(min(d$Age10_c), min(d$Age10_c) + 1)), # Welcher Vergleich? Hier: minimum + 10 Jahre
type = "response" # In welcher Metrik? Hier: Wahrscheinlichkeit
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.0494 0.00516 9.59 <0.001 69.9 0.0393 0.0596
Term: Age10_c
Type: response
Comparison: -2.397583081571 - -3.397583081571Bivariate Logistische Regression mit PK3 und Alter Average counterfactual comparison Alte
m1 |>
avg_comparisons(
variables = list(Age10_c = c(max(d$Age10_c) - 1, max(d$Age10_c))), # Welcher Vergleich? Hier: maximum - 10 Jahre
type = "response" # In welcher Metrik? Hier: Wahrscheinlichkeit
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.0693 0.0124 5.59 <0.001 25.4 0.045 0.0936
Term: Age10_c
Type: response
Comparison: 1.802416918429 - 0.802416918429003Bivariate Logistische Regression mit PK3 und Alter Odds Ratios
m1 |>
report_table(metrics = "R2", include_effectsize = FALSE, exponentiate = TRUE)Parameter | Coefficient | 95% CI | z | p | Fit
----------------------------------------------------------------
(Intercept) | 0.67 | [0.59, 0.77] | -5.97 | < .001 |
Age10 c | 1.32 | [1.20, 1.46] | 5.60 | < .001 |
| | | | |
Tjur's R2 | | | | | 0.03Multiple Logistische Regression mit PK3
d <- d |>
mutate(
Online_news_sites_c = Online_news_sites - mean(Online_news_sites),
Twitter_c = Twitter - mean(Twitter),
Facebook_c = Facebook - mean(Facebook),
Age10_c = Age10 - mean(Age10)
)
m2 <- glm(
PK3 ~ Online_news_sites_c + Twitter_c + Facebook_c + Gender + Age10_c, # Regressionsgleichung
family = binomial(link = "logit"), # aV ist binär, logit-Link = Logistische Regression
data = d
)Multiple Logistische Regression mit PK3 Tabelle
m2 |>
report_table(metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | p | Fit
--------------------------------------------------------------------------
(Intercept) | 0.13 | [-0.05, 0.31] | 1.42 | 0.157 |
Online news sites c | 0.19 | [ 0.11, 0.28] | 4.50 | < .001 |
Twitter c | -0.03 | [-0.18, 0.12] | -0.42 | 0.675 |
Facebook c | -0.11 | [-0.19, -0.04] | -2.82 | 0.005 |
Gender [female] | -1.20 | [-1.48, -0.92] | -8.35 | < .001 |
Age10 c | 0.19 | [ 0.09, 0.30] | 3.55 | < .001 |
| | | | |
Tjur's R2 | | | | | 0.14Multiple Logistische Regression mit PK3 Baseline Wahrscheinlichkeit
plogis(0.13)[1] 0.5324543Multiple Logistische Regression mit PK3 Average counterfactual comparison M+1
m2 |>
avg_comparisons(
variables = list(
Online_news_sites_c = 1, Twitter_c = 1, Facebook_c = 1, Gender = "reference", Age10_c = 1
), # Welcher Vergleich? Hier: Mittelwert und M + 1 und Faktor
type = "response" # In welcher Metrik? Hier: Wahrscheinlichkeit
)
Term Contrast Estimate Std. Error z Pr(>|z|) S
Age10_c +1 0.04062 0.01128 3.60 < 0.001 11.6
Facebook_c +1 -0.02352 0.00814 -2.89 0.00388 8.0
Gender female - male -0.26638 0.03064 -8.69 < 0.001 58.0
Online_news_sites_c +1 0.04060 0.00877 4.63 < 0.001 18.0
Twitter_c +1 -0.00669 0.01592 -0.42 0.67452 0.6
2.5 % 97.5 %
0.0185 0.06273
-0.0395 -0.00755
-0.3264 -0.20633
0.0234 0.05780
-0.0379 0.02452
Type: responseMultiple Logistische Regression mit PK3 Odds Ratios
m2 |>
report_table(metrics = "R2", include_effectsize = FALSE, exponentiate = TRUE)Parameter | Coefficient | 95% CI | z | p | Fit
------------------------------------------------------------------------
(Intercept) | 1.14 | [0.95, 1.36] | 1.42 | 0.157 |
Online news sites c | 1.21 | [1.12, 1.32] | 4.50 | < .001 |
Twitter c | 0.97 | [0.83, 1.12] | -0.42 | 0.675 |
Facebook c | 0.89 | [0.82, 0.97] | -2.82 | 0.005 |
Gender [female] | 0.30 | [0.23, 0.40] | -8.35 | < .001 |
Age10 c | 1.21 | [1.09, 1.35] | 3.55 | < .001 |
| | | | |
Tjur's R2 | | | | | 0.14Multiple Logistische Regression mit PK3 R2 Tabelle
m2 |>
report_table(metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | p | Fit
--------------------------------------------------------------------------
(Intercept) | 0.13 | [-0.05, 0.31] | 1.42 | 0.157 |
Online news sites c | 0.19 | [ 0.11, 0.28] | 4.50 | < .001 |
Twitter c | -0.03 | [-0.18, 0.12] | -0.42 | 0.675 |
Facebook c | -0.11 | [-0.19, -0.04] | -2.82 | 0.005 |
Gender [female] | -1.20 | [-1.48, -0.92] | -8.35 | < .001 |
Age10 c | 0.19 | [ 0.09, 0.30] | 3.55 | < .001 |
| | | | |
Tjur's R2 | | | | | 0.14Multiple Logistische Regression mit PK3 R2 Text
m2 |>
report_performance() |>
cat()The model’s explanatory power is moderate (Tjur’s R2 = 0.14)
Varianten der Modellspezifikation für Test mit 5 Fragen
m5q_v1 <- glm(
cbind(Political_knowledge, 5 - Political_knowledge) ~ Age10_c, # Zahl der richtig und falsch beantworteten Fragen in cbind()
family = binomial(link = "logit"), # eine Frage ist binär, logit-Link = Logistische Regression
data = d
)
m5q_v2 <- glm(
Political_knowledge / 5 ~ Age10_c, # Anteil der richtig beantworteten Fragen als aV
family = binomial(link = "logit"), # eine Frage ist binär, logit-Link = Logistische Regression
weights = rep(5, nrow(d)), # Zahl der Fragen als Gewicht
data = d
)Ausgabe Logistische Regression mit 5 Fragen
m5q_v1 |>
report_table(metrics = "R2", include_effectsize = FALSE)Warning: Can't calculate accurate R2 for binomial models that are not Bernoulli
models.Warning: Models of class `glm` are not yet supported.Parameter | Coefficient | 95% CI | z | p | Fit
---------------------------------------------------------------
(Intercept) | 0.45 | [0.40, 0.51] | 15.34 | < .001 |
Age10 c | 0.25 | [0.21, 0.29] | 11.78 | < .001 |
| | | | | Intercept Logistische Regression mit 5 Fragen
round(plogis(coef(m5q_v1))[1], 2)(Intercept)
0.61 Plot Lineare Regression mit Political_knowledge
d |>
ggplot(aes(Age10, Political_knowledge)) +
geom_point(position = position_jitter(width = 0, height = 0.1, seed = 1), shape = 1) +
labs(x = "Alter in Jahrzehnten (Age / 10)", y = "Richtig beantwortete Fragen") +
geom_smooth(method = "lm", se = F, linewidth = 2)`geom_smooth()` using formula = 'y ~ x'
Plot Logistische Regression mit Political_knowledge
d |>
ggplot(aes(Age10, Political_knowledge / 5, weight = 5)) +
geom_point(position = position_jitter(width = 0, height = 0.02, seed = 1), shape = 1) +
labs(x = "Alter in Jahrzehnten (Age / 10)", y = "Richtig beantwortete Fragen") +
geom_smooth(
method = "glm", method.args = list(family = "binomial"), se = F,
linewidth = 2
)`geom_smooth()` using formula = 'y ~ x'
Mittelwert und Varianz der Beispiel-Variable
d2_RU <- d2 |>
mutate(
ym = floor_date(ymd_hms(created_time), "month"),
ym_int = as.integer(as.factor(ym)) - min(as.integer(as.factor(ym)), na.rm = TRUE)
) |>
filter(uni == "UC San Francisco") |>
na.omit()Warning: There was 1 warning in `mutate()`.
ℹ In argument: `ym = floor_date(ymd_hms(created_time), "month")`.
Caused by warning:
! 33 failed to parse.d2_RU |>
summarise(
M = mean(comments_count),
Var = var(comments_count),
n = n()
)# A tibble: 1 × 3
M Var n
<dbl> <dbl> <int>
1 3.55 28.8 128Plot Lineare Regression mit comments_count - Simuliertes Beispiel
d_fake <- d2_RU |>
mutate(
lambda = exp(2 + 0.08 * ym_int + rnorm(n = n(), mean = 0, sd = 0.3)),
comments_count_fake = rpois(n = n(), lambda = lambda)
)
plt_count_fake <- d_fake |>
ggplot(aes(ym, comments_count_fake)) +
geom_point(position = position_jitter(width = 0, height = 0.1, seed = 1), shape = 1) +
labs(x = "Zeit (in Monatsschritten)", y = "Zahl der Kommentare")
plt_count_fake +
geom_smooth(method = "lm", linewidth = 2, se = F)`geom_smooth()` using formula = 'y ~ x'
Plot Poisson- und negativ-binomiale Regression mit comments_count - Simuliertes Beispiel
plt_count_fake +
geom_smooth(
method = "glm", method.args = list(family = "poisson"),
se = F, linewidth = 1, color = "green", linetype = 2
) +
geom_smooth(
method = MASS::glm.nb,
se = F, linewidth = 1
)`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
Plot Lineare Regression mit comments_count
plt_count <- d2_RU |>
ggplot(aes(ym, comments_count)) +
geom_point(position = position_jitter(width = 0, height = 0.1, seed = 1), shape = 1) +
labs(x = "Zeit (in Monatsschritten)", y = "Zahl der Kommentare")
plt_count +
geom_smooth(method = "lm", linewidth = 2, se = F)`geom_smooth()` using formula = 'y ~ x'
Plot Negativ-binomiale Regression mit comments_count
plt_count +
geom_smooth(
method = MASS::glm.nb,
se = F, linewidth = 2
)`geom_smooth()` using formula = 'y ~ x'
Null-Modell Negativ-binomiale Regression mit comments_count
nbm0 <- glm.nb( # Negativ-binomiale Regression
comments_count ~ 1, # Null-Modell nur mit Intercept (Konstante)
data = d2_RU
)Null-Modell Negativ-binomiale Regression mit comments_count Tabelle
nbm0 |>
report_table(metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | df | p | Fit
----------------------------------------------------------------------------
(Intercept) | 1.27 | [1.02, 1.52] | 10.06 | Inf | < .001 |
| | | | | |
R2_Nagelkerke | | | | | | 3.38e-16exp Funktion
exp(coef(nbm0))(Intercept)
3.554687 Bivariate negativ-binomiale Regression mit comments_count und Zeit
d2_RU <- d2_RU |>
mutate(ym_int_half = ym_int - (max(ym_int) + 1) / 2)
nbm1 <- glm.nb( # Negativ-binomiale Regression
comments_count ~ ym_int_half, # Regressionsgleichung
data = d2_RU
)Bivariate negativ-binomiale Regression mit comments_count und Zeit Tabelle
nbm1 |>
report_table(exponentiate = TRUE, metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | df | p | Fit
-----------------------------------------------------------------------
(Intercept) | 3.22 | [2.57, 4.07] | 9.97 | Inf | < .001 |
ym int half | 1.05 | [1.03, 1.07] | 4.20 | Inf | < .001 |
| | | | | |
R2_Nagelkerke | | | | | | 0.21Bivariate negativ-binomiale Regression comments_count und Zeit Average counterfactual comparison M+1
nbm1 |>
avg_comparisons(
variables = list(ym_int_half = c(0, 1)), # Welcher Vergleich? Hier: Mitte der Zeit und M + 1 Monat
type = "response" # In welcher Metrik? Hier: Erwartete Zahl der Kommentare
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.154 0.0415 3.72 <0.001 12.3 0.073 0.236
Term: ym_int_half
Type: response
Comparison: 1 - 0Bivariate negativ-binomiale Regression mit comments_count und Zeit Average counterfactual comparison M-1
nbm1 |>
avg_comparisons(
variables = list(ym_int_half = c(0, -1)), # Welcher Vergleich? Hier: Mitte der Zeit und M - 1 Monat
type = "response" # In welcher Metrik? Hier: Erwartete Zahl der Kommentare
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.147 0.0381 3.86 <0.001 13.1 0.0726 0.222
Term: ym_int_half
Type: response
Comparison: 0 - -1Bivariate negativ-binomiale Regression mit comments_count und Zeit Average counterfactual comparison Junge
nbm1 |>
avg_comparisons(
variables = list(ym_int_half = c(min(d2_RU$ym_int_half), min(d2_RU$ym_int_half) + 1)), # Welcher Vergleich? Hier: minimum + 1 Monat
type = "response" # In welcher Metrik? Hier: Erwartete Zahl der Kommentare
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.0664 0.00824 8.05 <0.001 50.1 0.0502 0.0826
Term: ym_int_half
Type: response
Comparison: -17 - -18Bivariate negativ-binomiale Regression mit comments_count und Zeit Average counterfactual comparison Alte
nbm1 |>
avg_comparisons(
variables = list(ym_int_half = c(max(d2_RU$ym_int_half) - 1, max(d2_RU$ym_int_half))), # Welcher Vergleich? Hier: maximum - 1 Monat
type = "response" # In welcher Metrik? Hier: Erwartete Zahl der Kommentare
)
Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
0.327 0.143 2.29 0.022 5.5 0.0472 0.606
Term: ym_int_half
Type: response
Comparison: 17 - 16Multiple negativ-binomiale Regression mit comments_count
d2_RU <- d2_RU |>
mutate(
ym_int_half = ym_int - (max(ym_int) + 1) / 2,
word_count_100_c = (word_count - mean(word_count, rm.na = TRUE)) / 100
)
nbm2 <- glm.nb(
comments_count ~ ym_int_half + topic_research + topic_teaching + word_count_100_c, # Regressionsgleichung
data = d2_RU
)Multiple negativ-binomiale Regression mit comments_count Tabelle
nbm2 |>
report_table(exponentiate = TRUE, metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | df | p | Fit
------------------------------------------------------------------------------
(Intercept) | 2.54 | [1.82, 3.61] | 5.31 | Inf | < .001 |
ym int half | 1.05 | [1.03, 1.08] | 4.40 | Inf | < .001 |
topic research [yes] | 1.23 | [0.77, 1.96] | 0.89 | Inf | 0.376 |
topic teaching [yes] | 1.81 | [1.00, 3.45] | 1.92 | Inf | 0.055 |
word count 100 c | 1.12 | [0.31, 4.29] | 0.19 | Inf | 0.848 |
| | | | | |
R2_Nagelkerke | | | | | | 0.26Multiple negativ-binomiale Regression mit comments_count Average counterfactual comparison
nbm2 |>
avg_comparisons(
variables = list(
ym_int_half = c(0, 1), topic_research = "reference",
topic_teaching = "reference",
word_count_100_c = 1
), # Welche Vergleiche?
type = "response" # In welcher Metrik? Hier: Erwartete Zahl der Likes
)
Term Contrast Estimate Std. Error z Pr(>|z|) S 2.5 %
topic_research yes - no 0.683 0.7762 0.880 0.379 1.4 -0.8385
topic_teaching yes - no 2.516 1.6437 1.531 0.126 3.0 -0.7057
word_count_100_c +1 0.388 2.1409 0.181 0.856 0.2 -3.8086
ym_int_half 1 - 0 0.157 0.0411 3.819 <0.001 12.9 0.0763
97.5 %
2.204
5.737
4.584
0.237
Type: responseMultiple negativ-binomiale Regression mit comments_count R2 Tabelle
nbm2 |>
report_table(exponentiate = TRUE, metrics = "R2", include_effectsize = FALSE)Parameter | Coefficient | 95% CI | z | df | p | Fit
------------------------------------------------------------------------------
(Intercept) | 2.54 | [1.82, 3.61] | 5.31 | Inf | < .001 |
ym int half | 1.05 | [1.03, 1.08] | 4.40 | Inf | < .001 |
topic research [yes] | 1.23 | [0.77, 1.96] | 0.89 | Inf | 0.376 |
topic teaching [yes] | 1.81 | [1.00, 3.45] | 1.92 | Inf | 0.055 |
word count 100 c | 1.12 | [0.31, 4.29] | 0.19 | Inf | 0.848 |
| | | | | |
R2_Nagelkerke | | | | | | 0.26Multiple negativ-binomiale Regression mit comments_count R2 Text
nbm2 |>
report_performance() |>
cat()The model’s explanatory power is substantial (Nagelkerke’s R2 = 0.26)
Poisson-Regression
pm2 <- glm(
comments_count ~ ym_int_half + topic_research + topic_teaching + word_count_100_c, # Regressionsgleichung
family = poisson(link = "log"), # aV Zählvariable, log-Link = Poisson-Regression
data = d2_RU
)Modellvergleich Koeffizienten
compare_models(Poisson = pm2, "Negativ-binomial" = nbm2, exponentiate = TRUE)Parameter | Poisson | Negativ-binomial
------------------------------------------------------------
(Intercept) | 2.46 (2.10, 2.88) | 2.54 (1.80, 3.59)
ym int half | 1.06 (1.05, 1.07) | 1.05 (1.03, 1.08)
topic research [yes] | 1.31 (1.08, 1.58) | 1.23 (0.78, 1.96)
topic teaching [yes] | 1.72 (1.37, 2.17) | 1.81 (0.99, 3.33)
word count 100 c | 1.10 (0.66, 1.83) | 1.12 (0.36, 3.48)
------------------------------------------------------------
Observations | 128 | 128Test der einzelnen Modelle: Poisson-Regression
check_overdispersion(pm2)# Overdispersion test
dispersion ratio = 5.822
Pearson's Chi-Squared = 716.094
p-value = < 0.001Overdispersion detected.Test der einzelnen Modelle: Negativ-binomiale Regression
check_overdispersion(nbm2)# Overdispersion test
dispersion ratio = 0.933
p-value = 0.976No overdispersion detected.Modellvergleich Informationkriterien
compare_performance(pm2, nbm2, metrics = c("AIC", "BIC"))# Comparison of Model Performance Indices
Name | Model | AIC (weights) | BIC (weights)
---------------------------------------------
pm2 | glm | 876.9 (<.001) | 891.1 (<.001)
nbm2 | negbin | 589.2 (>.999) | 606.3 (>.999)Modellvergleich Test
library(lmtest)Loading required package: zoo
Attaching package: 'zoo'The following objects are masked from 'package:base':
as.Date, as.Date.numericlrtest(pm2, nbm2)Likelihood ratio test
Model 1: comments_count ~ ym_int_half + topic_research + topic_teaching +
word_count_100_c
Model 2: comments_count ~ ym_int_half + topic_research + topic_teaching +
word_count_100_c
#Df LogLik Df Chisq Pr(>Chisq)
1 5 -433.44
2 6 -288.59 1 289.7 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 121.4 Hausaufgabe (Freiwillig)
- Reproduzieren Sie die verallgemeinerten linearen Modelle aus der Vorlesung. Interpretieren Sie die Ergebnisse.
Lösung
21.5 Transkript
WarnungHinweise zum automatisiert erstellten Transkript
Das folgende Transkript wurde auf Basis der Aufzeichnung der Vorlesung erstellt. Die vollständigen Aufzeichnungen inklusive der Bildschirminhalte sind in Blackboard🔒 verfügbar. Die Tonspur wurde zuerst mit Hilfe der Werkzeuge des Oral-History.Digital Projekts wörtlich transkribiert. Die wörtliche Transkription wurde in Kombination mit den Vorlesungsfolien mithilfe von Sprachmodellen (v. a. Claude Sonnet 4.5 und GPT 5.2) zu einem übersichtlichen Transkript zusammengefasst. Im Anschluss wurde das Transkript von einer studentischen Hilfskraft überprüft, geglättet und ggf. angepasst. In diesem Prozess kann es an verschiedenen Stellen zu Fehlern kommen. Im Zweifel gilt das gesprochene Wort, und auch beim Vortrag mache ich Fehler.
Ich stelle das Transkript hier als experimentelles, ergänzendes Material zur Dokumentation der Vorlesung zur Verfügung. Noch bin ich mir unsicher, ob es eine sinnvolle Ergänzung ist und behalte mir vor, es weiter zu bearbeiten oder zu löschen.
Verallgemeinerte lineare Modelle
Warum verallgemeinerte lineare Modelle?
Bislang wurde in der Vorlesung immer mit dem linearen Modell gearbeitet: Ob bivariate Regression, multiple Regression, Moderation, Mediation, Strukturgleichungsmodelle oder Mehrebenenmodelle – in all diesen Fällen wurde angenommen, dass die abhängige Variable (aV), also die Variable, die erklärt oder vorhergesagt werden soll, metrisch und kontinuierlich ist. Das bedeutet, dass sie theoretisch jeden beliebigen Wert annehmen kann und dass Abstände zwischen Werten gleich interpretierbar sind (z. B. ist der Wert 4 doppelt so groß wie der Wert 2). Diese Annahme ist bei vielen Variablen, etwa bei Likert-Skalen, näherungsweise akzeptabel, weil man von einer dahinterliegenden, latenten kontinuierlichen Skala ausgeht. Es gibt jedoch Variablentypen, bei denen diese Annahme eindeutig nicht erfüllt ist.
Zwei solche Fälle stehen im Zentrum dieser Einheit:
- Binäre Variablen: Sie können nur zwei Ausprägungen annehmen, typischerweise codiert als 0 und 1 (z. B. ‘kommt in einem Social-Media-Post vor’ vs. ‘kommt nicht vor’, oder ‘Wissensfrage richtig beantwortet’ vs. ‘falsch beantwortet’). Diese Variablen sind diskret; ein Fall kann nicht ‘ein bisschen richtig’ sein.
- Zählvariablen (Count-Variablen): Das sind nicht-negative ganze Zahlen, die entstehen, wenn etwas gezählt wird, z. B. die Anzahl der Likes oder Kommentare zu einem Social-Media-Post. Auch hier sind nur diskrete, nicht-negative Werte möglich; negative oder gebrochene Zählwerte ergeben keinen Sinn.
Wenn man solche Variablen mit einem normalen linearen Modell modelliert, entstehen zwei Probleme: Erstens ist die Beziehung zwischen den Prädiktoren und der aV in der Regel nicht linear, weil die aV nur bestimmte diskrete Werte oder einen begrenzten Wertebereich annehmen kann. Zweitens sagt das lineare Modell mitunter Werte vorher, die es in der Realität gar nicht geben kann – zum Beispiel negative Wahrscheinlichkeiten oder negative Zählwerte. Genau deshalb braucht man verallgemeinerte lineare Modelle (Generalized Linear Models, GLM): Sie behalten die Grundlogik der linearen Regression bei, verwenden aber eine passende Verteilungsannahme für die aV sowie eine sogenannte Link-Funktion, die die lineare Vorhersage in den jeweils gültigen Wertebereich überträgt.
Die Einheit gliedert sich in zwei inhaltliche Blöcke: zunächst logistische Modelle für binäre abhängige Variablen, danach Poisson- und negativ-binomiale Modelle für Zählvariablen als aV. Beide Modellarten teilen viele Prinzipien mit der linearen Regression (z. B. die Logik von Nullmodell und Prädiktoren, die statistische Inferenz mit Standardfehlern, Teststatistiken und Konfidenzintervallen), unterscheiden sich aber in der Verteilungsannahme und vor allem in der Interpretation der Koeffizienten.
Binäre Variablen als abhängige Variable: Logistische Modelle
Ausgangsbeispiel: Wissensfragen bei Van Erkel & Van Aelst (2021)
Als Anwendungsbeispiel dient erneut die aus der Vorlesung bekannte Studie von Van Erkel und Van Aelst (2021) zum politischen Surveillance-Wissen. Befragten wurden sechs politische Wissensfragen gestellt (von denen eine Frage, PK5, nicht in den späteren Wissensindex eingeht). Für jede Frage wurde erfasst, ob sie richtig oder falsch beantwortet wurde. Der Anteil richtiger Antworten schwankte deutlich zwischen den Fragen: Bei der leichtesten Frage (PK4) antworteten 87.9 % richtig, bei der schwierigsten (PK6) nur 24.1 %. Diese einzelnen Fragen sind binäre Variablen (richtig = 1, falsch = 0). Bildet man aus den fünf in den Index eingehenden Fragen eine Summe, erhält man eine additive Wissens-Zählvariable mit Werten von 0 bis 5, die in der Verteilung eine Glockenform mit Schwerpunkt bei 3 zeigt. Damit liefert dieselbe Studie sowohl ein Beispiel für eine binäre aV (einzelne Frage) als auch für eine Zählvariable (Summenindex), auf die im späteren Verlauf der Einheit zurückgekommen wird.
Warum die lineare Regression bei binärer aV versagt
Um das Problem zu veranschaulichen, wird zunächst ein simuliertes Beispiel mit einem besonders starken (künstlichen) Zusammenhang gezeigt: das Alter einer Person (in Jahrzehnten gemessen, also durch 10 geteilt) als Prädiktor für die richtige Beantwortung einer Wissensfrage (0 = falsch, 1 = richtig). Legt man hier – wie bei der linearen Regression üblich – einfach eine gerade Linie durch die Punktwolke, ergeben sich zwei Probleme:
- Die vorhergesagten Werte können außerhalb des sinnvollen Bereichs von 0 und 1 liegen, das heißt es werden negative oder über 1 liegende ‘Wahrscheinlichkeiten’ vorhergesagt. Das ist mathematisch und inhaltlich nicht sinnvoll, da eine Wahrscheinlichkeit niemals kleiner als 0 oder größer als 1 sein kann.
- Ein linearer Zusammenhang ist bei einer binären aV konzeptionell nicht plausibel, da es keine ‘richtigeren’ oder ‘falscheren’ Antworten gibt, sondern nur die beiden diskreten Zustände richtig/falsch.
Dieses Vorgehen – eine 0/1-aV direkt linear zu modellieren – wird auch als Linear Probability Model (LPM) bezeichnet. Es ist streng genommen fehlerhaft, wird in der Praxis und in der Forschungsliteratur aber gelegentlich dennoch eingesetzt, weil es einfacher zu interpretieren ist und bei gemäßigten Zusammenhängen oft nur wenig von den Ergebnissen der korrekteren logistischen Regression abweicht.
Die logistische Regression und die Logit-Funktion
Die Lösung besteht darin, statt einer geraden Linie eine S-förmige Kurve durch die Datenpunkte zu legen. Diese S-Kurve sorgt dafür, dass die vorhergesagten Werte immer zwischen 0 und 1 bleiben – egal wie extrem der Wert des Prädiktors ist. Diese Transformation wird über die sogenannte Logit-Funktion realisiert:
- Die Regressionsgleichung selbst bleibt wie gewohnt aufgebaut (Intercept plus Regressionsgewicht mal Prädiktor), sie wird aber nicht mehr direkt als vorhergesagter Wert interpretiert, sondern zunächst als sogenannter Logit.
- Über die Logit-Funktion \(P(Y_i) = 1 / (1 + e^{-(b_0 + b_1 \cdot X_i)})\) wird dieser Logit in eine Wahrscheinlichkeit zwischen 0 und 1 zurückgerechnet.
- In der Mitte der Verteilung verläuft die Kurve nahezu linear, an den Rändern (nahe 0 und nahe 1) flacht sie zunehmend ab, sodass extreme Werte des Prädiktors nie zu Wahrscheinlichkeiten unter 0 oder über 1 führen.
Man muss die Formel der Logit-Funktion nicht auswendig lernen; entscheidend ist das Verständnis, dass die logistische Regression genau dieses Problem der linearen Regression löst, indem sie die lineare Vorhersage durch eine Transformation in den gültigen Wahrscheinlichkeitsbereich überträgt.
Gemeinsamkeiten und Unterschiede: lineare vs. logistische Regression
Logistische Regression teilt viele Grundprinzipien mit der linearen Regression, unterscheidet sich aber in wichtigen Punkten:
- Gleiche Annahmen bezüglich Unabhängigkeit der Fälle und (fehlender) Multikollinearität der Prädiktoren.
- Gleiche Logik der Modellspezifikation: Man beginnt mit einem Nullmodell (nur Intercept) und erweitert es schrittweise um Prädiktorvariablen.
- Gleiche Logik der statistischen Inferenz: Es werden Standardfehler, Teststatistiken (hier: \(z\)-Werte statt \(t\)-Werte) und Konfidenzintervalle berechnet und interpretiert.
- Statistisch unterschiedliche, aber konzeptionell ähnliche Gütemaße (Pseudo-\(R^2\)) kommen zum Einsatz.
- Die Koeffizienten werden grundlegend anders interpretiert als in der linearen Regression, da sie sich nicht mehr auf der Wahrscheinlichkeitsskala, sondern auf der Logit-Skala befinden.
Das Nullmodell und der Intercept
Wie bei der linearen Regression beginnt man mit einem Nullmodell, das nur aus dem Intercept (der Regressionskonstante) besteht. Am Beispiel der dritten Wissensfrage (PK3) aus dem Datensatz von Van Erkel und Van Aelst ergibt sich ein Intercept von \(-0.39\) auf der Logit-Skala. Dieser Wert lässt sich nicht direkt als Wahrscheinlichkeit lesen, sondern muss über die inverse Logit-Funktion (auch plogis-Funktion genannt) zurücktransformiert werden. Rechnet man \(1 / (1 + e^{-(-0.39)})\) beziehungsweise nutzt die Funktion plogis(-0.39), erhält man rund \(0.41\) bzw. 40.7 %. Das entspricht exakt dem deskriptiven Anteil richtiger Antworten bei dieser Frage. Das Nullmodell einer logistischen Regression bildet also – genau wie bei der linearen Regression – schlicht den Mittelwert der aV ab, hier eben in Form der beobachteten Erfolgsquote.
Als Gütemaß für das Nullmodell wird Tjur’s \(R^2\) ausgewiesen, das hier folgerichtig den Wert 0 annimmt, da ohne Prädiktoren noch keine Varianz erklärt werden kann.
Koeffizienten in der logistischen Regression interpretieren
Die Koeffizienten einer logistischen Regression liegen auf der Logit-Skala und sind daher nicht direkt als Veränderung von Wahrscheinlichkeiten interpretierbar. Das liegt daran, dass die Veränderung der Wahrscheinlichkeit bei einer Erhöhung des Prädiktors um eine Einheit nicht überall in der Verteilung gleich groß ist (siehe die S-Form der Kurve): In der Mitte der Verteilung ist die Steigung am größten, an den Rändern flacht sie ab. Zudem hängt die Interpretation zusätzlich vom Wert des Intercepts und der übrigen Prädiktoren ab. Es gibt drei gängige Möglichkeiten, logistische Koeffizienten dennoch verständlich zu interpretieren:
Geteilt-durch-4-Regel
Die einfachste Faustregel teilt den Koeffizienten durch 4. Dieser Wert (\(B/4\)) stellt eine Obergrenze für die tatsächliche Veränderung der Wahrscheinlichkeit dar, wenn sich der Prädiktor um eine Einheit ändert, ausgedrückt in Prozentpunkten. Die Abweichung von dieser Obergrenze ist am Rand der Verteilung des Prädiktors größer als in deren Mitte, da dort die S-Kurve stärker abflacht.
Am Beispiel des Modells mit dem (mittelwertzentrierten) Alter als Prädiktor für PK3 ergibt sich ein Koeffizient von \(0.28\). Geteilt durch 4 ergibt das etwa \(0.07\), also 7 Prozentpunkte. Die Interpretation lautet: Personen, die sich im Alter um 10 Jahre (eine Einheit des in Jahrzehnten codierten Prädiktors) unterscheiden, unterscheiden sich in etwa um 7 Prozentpunkte in der Wahrscheinlichkeit, die Frage richtig zu beantworten – die ältere Person hat die höhere Wahrscheinlichkeit. Der Intercept von \(-0.39\) entspricht dabei, wie oben gezeigt, einer Baseline-Wahrscheinlichkeit von rund 40 % für eine Person mit Durchschnittsalter.
Average Counterfactual Comparison
Eine genauere, aber rechenaufwendigere Alternative ist der Vergleich vorhergesagter Wahrscheinlichkeiten für hypothetische (kontrafaktische) Personen, wie er bereits aus der Einheit zur multiplen linearen Regression bekannt ist. Dabei wird zum Beispiel berechnet, wie stark sich die vorhergesagte Wahrscheinlichkeit ändert, wenn man eine Person um 10 Jahre älter macht. Für das Altersbeispiel zeigt sich: Im Durchschnitt beantwortet eine um 10 Jahre ältere Person die Frage mit einer um etwa 6 bis 7 Prozentpunkte höheren Wahrscheinlichkeit richtig, wenn man den Vergleich um den Mittelwert herum oder bei jüngeren Personen betrachtet. Am Rand der Verteilung, also bei bereits sehr alten Personen, fällt dieser Unterschied kleiner aus (rund 4 bis 5 Prozentpunkte), weil die S-Kurve dort bereits abflacht. Diese Methode ist somit präziser als die Geteilt-durch-4-Regel, weil sie explizit berücksichtigt, dass der Effekt je nach Position auf der Kurve unterschiedlich groß ausfällt.
Odds Ratios (exponenzierte Koeffizienten)
Die dritte, in Publikationen sehr verbreitete Möglichkeit ist die Angabe von Odds Ratios, also von exponenzierten Koeffizienten (\(Exp(B)\)). Um Odds Ratios zu verstehen, muss man zunächst den Begriff der Odds (auf Deutsch: Chance oder Risiko) von dem Begriff der Wahrscheinlichkeit abgrenzen:
- \(Odds = P(Y=1) / (1 - P(Y=1))\). Odds sind also das Verhältnis der Wahrscheinlichkeit, dass ein Ereignis eintritt, zur Wahrscheinlichkeit, dass es nicht eintritt.
- Man spricht von ‘Chance’, wenn das Ereignis positiv konnotiert ist (z. B. die Frage richtig zu beantworten), und von ‘Risiko’, wenn es negativ konnotiert ist (z. B. an einer Krankheit zu erkranken).
- Sind beide möglichen Ausgänge gleich wahrscheinlich (jeweils 50 %), betragen die Odds genau 1.
Eine Odds Ratio ist der Quotient aus zwei solchen Odds, also die Odds einer Person geteilt durch die Odds einer anderen Person (bzw. einer Person mit einem um eine Einheit höheren Prädiktorwert). Odds Ratios lassen sich folgendermaßen interpretieren:
- \(Exp(B) = 1\) bedeutet: kein Unterschied zwischen den verglichenen Personen.
- \(Exp(B) > 1\) bedeutet: eine erhöhte Chance bzw. ein erhöhtes Risiko für die Person mit dem höheren Prädiktorwert.
- \(Exp(B) < 1\) bedeutet: eine niedrigere Chance bzw. ein niedrigeres Risiko.
- \(Exp(B) = 0.5\) bedeutet: Mit jeder Einheit mehr auf dem Prädiktor besteht nur noch ein halb so großes Risiko bzw. eine halb so große Chance.
- \(Exp(B) = 2\) bedeutet: Mit jeder Einheit mehr auf dem Prädiktor verdoppelt sich die Chance bzw. das Risiko.
Für das Altersbeispiel ergibt sich ein Odds Ratio von \(1.32\) für den Prädiktor Alter (in Jahrzehnten). Das bedeutet: Eine Person, die 10 Jahre älter ist, hat eine etwa 1.3-mal so große Chance, die Frage richtig zu beantworten, wie eine 10 Jahre jüngere Person. Wichtig zu beachten ist, dass der exponenzierte Intercept in der Regel keine sinnvoll interpretierbare Bedeutung hat und deshalb meist nicht berichtet wird.
Modell mit mehreren Prädiktoren
Wie in der multiplen linearen Regression lässt sich auch die logistische Regression um mehrere Prädiktoren erweitern. Im Beispiel wird die dritte Wissensfrage (PK3) durch die Nutzung von Online-Nachrichtenseiten, Twitter, Facebook, das (selbstberichtete) Geschlecht sowie das Alter vorhergesagt. Alle metrischen Prädiktoren, für die der Wert 0 keine sinnvolle Bedeutung hätte, werden vorab mittelwertzentriert (also der Mittelwert wird von jedem Wert abgezogen), damit der Intercept sinnvoll interpretierbar bleibt. Die binäre Geschlechtsvariable (0 = männlich, 1 = weiblich) wird dagegen nicht zentriert, da sie bereits einen sinnvollen Nullpunkt (männlich) besitzt.
In diesem Modell zeigt sich, dass die Twitter-Nutzung statistisch nicht signifikant mit der richtigen Beantwortung der Frage zusammenhängt – hier kann kein verlässlicher Effekt in der Grundgesamtheit angenommen werden. Für alle übrigen Prädiktoren (Online-Nachrichten-Nutzung, Facebook-Nutzung, Geschlecht, Alter) zeigen sich hingegen statistisch signifikante Zusammenhänge, die inhaltlich interpretiert werden:
- Der Intercept liegt bei \(0.13\) Logits, was einer Wahrscheinlichkeit von rund 53 % entspricht (berechnet über
plogis(0.13)). Das ist die vorhergesagte Wahrscheinlichkeit für Männer mit Durchschnittsalter und durchschnittlicher Nutzung der drei Mediumsvariablen, die richtige Antwort zu geben – also für Personen, die auf allen Prädiktoren den Wert 0 aufweisen. - Da diese Baseline-Wahrscheinlichkeit von 53 % nahe der Mitte der S-Kurve liegt, funktioniert die Geteilt-durch-4-Regel hier besonders gut als Annäherung.
- Für die Nutzung von Online-Nachrichtenseiten (gemessen auf einer 7er-Skala von 1 = keine Nutzung bis 7 = sehr viel Nutzung) ergibt sich ein Koeffizient von \(0.19\) mit einem Konfidenzintervall von etwa \(0.11\) bis \(0.28\) und einem statistisch signifikanten Effekt (\(z = 4.5\); \(p < 0.001\)). Zwei Personen, die sich nur in dieser Nutzung um einen Skalenpunkt unterscheiden, sich aber auf allen anderen Prädiktoren gleichen, unterscheiden sich um etwa 5 Prozentpunkte (\(0.19 / 4\)) in der Wahrscheinlichkeit, die Frage richtig zu beantworten – die Person mit stärkerer Nutzung hat die höhere Wahrscheinlichkeit.
- Für das Geschlecht (Koeffizient \(-1.20\)) ergibt sich nach der Geteilt-durch-4-Regel ein Unterschied von etwa 30 Prozentpunkten: Frauen beantworten die Frage – bei gleichem Alter und gleicher Mediennutzung wie vergleichbare Männer – mit einer um rund 30 Prozentpunkte geringeren Wahrscheinlichkeit richtig.
Das grundlegende Prinzip der Interpretation ist dabei identisch zur multiplen linearen Regression: Man vergleicht stets zwei Personen (oder Fälle), die sich auf genau einem Prädiktor um eine Einheit unterscheiden, auf allen übrigen im Modell enthaltenen Prädiktoren aber exakt gleiche Werte aufweisen. Der Unterschied liegt lediglich darin, dass bei der logistischen Regression diese Vergleiche über die Geteilt-durch-4-Regel, die Average Counterfactual Comparison oder Odds Ratios in Prozentpunkte bzw. Chancenverhältnisse übersetzt werden müssen, statt sie – wie in der linearen Regression – direkt als Punkte auf der aV lesen zu können.
Auch für dasselbe Modell mit mehreren Prädiktoren lassen sich alternativ Odds Ratios berichten: Für die Online-Nachrichten-Nutzung ergibt sich ein Odds Ratio von \(1.21\) (die Chance einer richtigen Antwort steigt also um den Faktor 1.2 pro Skalenpunkt), für das weibliche Geschlecht ein Odds Ratio von \(0.30\) (Frauen haben demnach nur eine etwa 0.3-mal so große Chance wie gleich alte Männer mit identischer Mediennutzung, die Frage richtig zu beantworten).
Modellgüte: Pseudo-R²
Da es in der logistische Regression – anders als in der linearen Regression – kein eindeutiges \(R^2\) gibt, wurden verschiedene sogenannte Pseudo-\(R^2\)-Maße entwickelt, die die Modellgüte durch den Vergleich von Nullmodell und dem zu beurteilenden Modell mit Prädiktoren abschätzen. Verbreitete Varianten sind McFadden’s \(R^2\), Nagelkerke’s \(R^2\), Cox & Snell’s \(R^2\) sowie Tjur’s \(R^2\). Für logistische Regressionsmodelle wird in dieser Vorlesung Tjur’s \(R^2\) empfohlen. Die Interpretation erfolgt grundsätzlich wie bei einem gewöhnlichen \(R^2\) (höhere Werte bedeuten mehr erklärte Varianz), wobei je nach Maß und Modell die theoretischen Grenzwerte 0 und 1 nicht notwendigerweise erreichbar sind. Im Beispielmodell mit mehreren Prädiktoren für PK3 ergibt sich ein Tjur’s \(R^2\) von \(0.14\), was als moderate Erklärungskraft des Modells eingeordnet wird.
Ausblick: die binomiale Zählvariable
Als Brücke zum nächsten Themenblock wird eine besondere Art von Zählvariable eingeführt: die binomiale Zählvariable. Das ist eine ganze, nicht-negative Zahl mit einem bekannten Maximum, die als Anzahl der Erfolge in einer bekannten Anzahl binärer Versuche verstanden werden kann. Ein Beispiel ist der aus den sechs (bzw. fünf indexrelevanten) Wissensfragen gebildete additive Wissensindex mit Werten von 0 bis 5: Er zählt, wie viele von fünf möglichen binären Wissensfragen richtig beantwortet wurden.
Solche binomialen Zählvariablen lassen sich weiterhin mit einer logistischen Regression modellieren. Dabei gibt man entweder die Anzahl der richtig und falsch beantworteten Fragen gemeinsam als aV ein (über die Funktion cbind()) oder man verwendet den Anteil richtig beantworteter Fragen als aV und die Gesamtzahl der Fragen als Gewicht (weights). In beiden Fällen wird geschätzt, mit welcher Wahrscheinlichkeit eine weitere ‘durchschnittliche’ Frage aus dem Test richtig beantwortet wird.
Im Beispiel mit dem (zentrierten) Alter als Prädiktor für den Wissensindex ergibt sich ein Intercept von \(0.45\) Logits, was einer Wahrscheinlichkeit von 61 % entspricht: Eine Person mit Durchschnittsalter beantwortet demnach eine durchschnittliche Frage mit 61-prozentiger Wahrscheinlichkeit richtig. Der Alters-Koeffizient von \(0.25\) führt nach der Geteilt-durch-4-Regel zu einer Interpretation von etwa 6 Prozentpunkten: Zwei Personen, die sich im Alter um 10 Jahre unterscheiden, unterscheiden sich um rund 6 Prozentpunkte in der Wahrscheinlichkeit, eine durchschnittliche Frage richtig zu beantworten. Vergleicht man in diesem Beispiel die lineare mit der logistischen Modellierung des Wissensindex, zeigt sich, dass bei diesem eher gemäßigten Zusammenhang praktisch kein Unterschied zwischen beiden Modellansätzen sichtbar wird – anders als im stark simulierten Ausgangsbeispiel.
Zählvariablen als abhängige Variable: Poisson- und negativ-binomiale Modelle
Warum die lineare Regression bei Zählvariablen versagt
Auch hier wird das Problem zunächst an einem simulierten Beispiel mit bekanntem, datengenerierendem Poisson-Prozess illustriert: die Zahl der Kommentare zu Posts der UC San Francisco im Zeitverlauf (gemessen in Monatsschritten). Legt man eine gewöhnliche lineare Regressionsgerade durch diese Datenpunkte, ergeben sich zwei Probleme:
- Die vorhergesagte Kommentarzahl kann kleiner als 0 werden – eine negative Anzahl von Kommentaren ist aber unmöglich.
- Weder der lineare Zusammenhang noch die Homoskedastizität (also die Annahme gleich großer Streuung über den gesamten Wertebereich) sind bei Zähldaten typischerweise gegeben, da die Streuung mit steigendem Erwartungswert häufig zunimmt.
Damit verstößt die lineare Regression bei Zählvariablen gegen zentrale eigene Annahmen.
Poisson- und negativ-binomiale Regression: Log-Link-Funktion
Die Lösung besteht darin, statt einer identitätsbasierten linearen Vorhersage eine Log-Link-Funktion zu verwenden. Die Regressionsgleichung wird dabei nicht mehr direkt auf die aV, sondern auf deren natürlichen Logarithmus angewendet:
- \(\log(Y_i) = b_0 + b_1 \cdot X_i\)
- Umgekehrt gilt: \(Y_i = e^{b_0 + b_1 \cdot X_i}\)
Die Exponentialfunktion transformiert die Vorhersagewerte in einen Bereich von 0 bis unendlich, sodass niemals negative Zählwerte vorhergesagt werden können. Im Beispiel führt dies zu einer nach oben gekrümmten Kurve, die sehr gut zum tatsächlichen (exponentiell anwachsenden) Verlauf der Kommentarzahlen im simulierten Datensatz passt.
Poisson- und negativ-binomiale Regression teilen viele Eigenschaften mit der linearen Regression, unterscheiden sich aber ebenfalls in bestimmten Punkten:
- Gleiche Annahmen zur Unabhängigkeit der Fälle und zur (fehlenden) Multikollinearität.
- Spezifische Verteilungsannahme: Bei Poisson entspricht die Varianz dem erwarteten Mittelwert, bei negativ-binomial dem Mittelwert multipliziert mit einem zusätzlichen Dispersionsparameter.
- Gleiche Logik der Modellspezifikation mit Nullmodell und Prädiktorvariablen sowie der statistischen Inferenz.
- Statistisch unterschiedliche, aber konzeptionell ähnliche Pseudo-\(R^2\)-Maße.
- Unterschiedliche, multiplikative Interpretation der Koeffizienten nach dem Exponenzieren (\(e^{b_1}\)).
Bei realistischen (weniger extremen) Zusammenhängen – wie sie in den tatsächlichen Daten von Fähnrich et al. (2020) vorliegen – sind die Unterschiede zwischen linearer und negativ-binomialer Regression deutlich weniger auffällig als im stark vereinfachten simulierten Beispiel, auch wenn die negativ-binomiale Regression weiterhin die theoretisch korrektere Wahl bleibt.
Das Nullmodell und der Intercept
Auch hier beginnt man mit einem Nullmodell, das nur den Intercept enthält (\(\log(Y_i) = b_0\)). Für die Kommentarzahl der UC-San-Francisco-Posts ergibt sich dabei ein Intercept von \(1.27\) auf der logarithmierten Skala. Exponenziert man diesen Wert (\(\exp(1.27)\)), erhält man rund \(3.55\) – das entspricht exakt dem beobachteten Mittelwert der Kommentarzahl in den Daten. Wie bei der linearen und der logistischen Regression bildet also auch hier das Nullmodell schlicht den Mittelwert der aV ab, in diesem Fall den erwarteten Mittelwert der Zählvariable. Als Gütemaß des Nullmodells wird ein \(R^2\) von 0 ausgewiesen (hier als Nagelkerke’s \(R^2\) bezeichnet).
Koeffizienten interpretieren
Koeffizienten von Modellen mit Log-Link sind ebenfalls schwer direkt interpretierbar, weil der multiplikative Effekt eines Prädiktors immer vom jeweiligen Ausgangsniveau der aV abhängt. Es werden zwei Interpretationsmöglichkeiten unterschieden:
Exponenzierte Koeffizienten als Wachstumsfaktor bzw. Wachstumsrate
Analog zu den Odds Ratios in der logistischen Regression exponenziert man auch hier die Koeffizienten. Man unterscheidet dabei zwischen dem Wachstumsfaktor (\(e^{b_1}\)) und der Wachstumsrate (\(e^{b_1} - 1\)). Am Beispiel eines Modells, das die Kommentarzahl über den Beobachtungszeitraum hinweg vorhersagt (der Zeitprädiktor wurde so kodiert, dass die Mitte des Beobachtungszeitraums den Wert 0 hat), ergibt sich ein exponenzierter Intercept von \(3.22\) und ein exponenzierter Zeitkoeffizient von \(1.05\):
- Der exponenzierte Intercept (\(3.22\)) gibt an: In der Mitte des Beobachtungszeitraums liegt die erwartete Zahl der Kommentare bei \(3.2\).
- Der Wachstumsfaktor von \(1.05\) bedeutet: Die erwartete Zahl der Kommentare wächst pro Monat auf das 1.05-Fache.
- Ausgedrückt als Wachstumsrate (\(1.05 - 1 = 0.05\)) bedeutet das: Die erwartete Zahl der Kommentare wächst pro Monat um 5 %.
Für das umfangreichere Modell mit mehreren Prädiktoren (Zeit, zwei Themen-Dummy-Variablen ‘Forschung’ und ‘Lehre’ sowie die Wortlänge des Posts in Hundert-Wort-Schritten, ebenfalls zentriert) ergibt sich ein exponenzierter Intercept von \(2.54\): In der Mitte des Beobachtungszeitraums liegt die erwartete Kommentarzahl für einen durchschnittlich langen Post ohne die beiden Themen bei \(2.5\). Der exponenzierte Zeitkoeffizient bleibt bei \(1.05\), das heißt, die erwartete Kommentarzahl wächst – bei sonst gleichen Eigenschaften des Posts – weiterhin um etwa 5 % pro Monat. Die Koeffizienten für die beiden Themen-Variablen und die Wortlänge sind in diesem Modell dagegen statistisch nicht signifikant, sodass für sie kein verlässlicher Effekt angenommen werden kann.
Average Counterfactual Comparison
Auch für Zählvariablen lässt sich – wie in der multiplen linearen Regression – der Effekt eines Prädiktors über den Vergleich hypothetischer (kontrafaktischer) Fälle an unterschiedlichen Stellen der Verteilung genauer bestimmen, statt nur mit einem einzigen multiplikativen Faktor zu arbeiten. Am Zeitbeispiel zeigt sich, dass der geschätzte Zuwachs an Kommentaren pro Monat je nach Position im Beobachtungszeitraum unterschiedlich groß ausfällt: Zu Beginn des Untersuchungszeitraums steigt die erwartete Kommentarzahl um etwa \(0.07\) pro Monat, in der Mitte des Zeitraums um \(0.15\) zusätzliche Kommentare pro Monat, und am Ende des Zeitraums (vom vorletzten zum letzten Monat) sogar um rund \(0.33\) Kommentare. Das macht deutlich, dass ein konstanter prozentualer Wachstumsfaktor in absoluten Kommentarzahlen zu sehr unterschiedlichen Zuwächsen führt, je nachdem, wie hoch das Ausgangsniveau bereits ist – ein typisches Merkmal exponentiellen Wachstums. Für die übrigen Prädiktoren im Mehrprädiktoren-Modell (Themen, Wortlänge) zeigen sich in dieser Betrachtung keine statistisch signifikanten Unterschiede.
Modellgüte: Pseudo-R²
Auch für Poisson- und negativ-binomiale Regression existieren verschiedene Pseudo-\(R^2\)-Maße (McFadden, Nagelkerke, Cox & Snell, Tjur), die die Passung von Nullmodell und geschätztem Modell vergleichen. Für diese beiden Modelltypen wird in der Vorlesung Nagelkerke’s \(R^2\) empfohlen. Auch hier gilt: Die Interpretation erfolgt analog zur linearen Regression, wobei die theoretischen Grenzwerte 0 und 1 je nach Maß und Modell mitunter nicht erreichbar sind. Im Beispielmodell mit mehreren Prädiktoren für die Kommentarzahl ergibt sich ein Nagelkerke’s \(R^2\) von \(0.26\), was als substanzielle Erklärungskraft des Modells eingestuft wird.
Modellvergleich: Poisson- oder negativ-binomiale Regression?
Am Ende des Blocks wird direkt verglichen, welches der beiden Modelle für die Kommentar-Daten besser geeignet ist. Vergleicht man die (exponenzierten) Koeffizienten eines Poisson-Modells mit denen des entsprechenden negativ-binomialen Modells für dieselben Prädiktoren, fallen die geschätzten Koeffizienten selbst recht ähnlich aus. Deutlich unterschiedlich ist jedoch die Inferenzstatistik: Im Poisson-Modell erscheinen die Unterschiede nach Themen (Forschung, Lehre) statistisch signifikant, im negativ-binomialen Modell dagegen nicht.
Um zu entscheiden, welches Modell verlässlicher ist, wird ein Overdispersion-Test durchgeführt, der prüft, ob die Annahme gleicher Varianz und gleichen Mittelwerts (Poisson) erfüllt ist. Für das Poisson-Modell ergibt sich ein Dispersionsverhältnis von \(5.82\) mit einem statistisch signifikanten Testergebnis (\(p < 0.001\)) – die Varianz in den Daten ist also weit größer als vom Poisson-Modell angenommen, es liegt eine deutliche Überdispersion vor. Für das negativ-binomiale Modell liegt das Dispersionsverhältnis dagegen bei \(0.93\) und ist nicht signifikant, das Modell passt also gut zur beobachteten Streuung.
Ergänzend werden Informationskriterien wie AIC (Akaike Information Criterion) und BIC (Bayesian Information Criterion) zum Modellvergleich herangezogen: Beide Kriterien fallen für das negativ-binomiale Modell (\(AIC = 589.2\)) deutlich niedriger aus als für das Poisson-Modell (\(AIC = 876.9\)), was ebenfalls klar für das negativ-binomiale Modell spricht (niedrigere Werte zeigen die bessere Modellpassung an).
Die inhaltliche Schlussfolgerung lautet: Da die Verteilungsannahme der Poisson-Regression (Varianz entspricht dem Erwartungswert) hier eindeutig verletzt ist, sind die Standardfehler des Poisson-Modells zu klein geschätzt, wodurch die Inferenzstatistik verzerrt ist. Die im Poisson-Modell gefundenen ‘statistisch signifikanten’ Unterschiede nach Themen wären daher mit hoher Wahrscheinlichkeit eine fälschliche Ablehnung der Nullhypothese. Die negativ-binomiale Regression passt wesentlich besser zu den vorliegenden Daten, und die dort gefundenen nicht signifikanten Ergebnisse für die Themen-Variablen gelten als die verlässlicheren Befunde. Für Social-Media-Daten sollte man also grundsätzlich beide Modelle schätzen und anhand von Overdispersion-Test sowie AIC/BIC das besser passende Modell auswählen.
Fazit
Die lineare Regression lässt sich verallgemeinern, um unterschiedliche Verteilungen der abhängigen Variable sowie passende Link-Funktionen zwischen den Prädiktoren und der aV zu ermöglichen. Damit können auch binäre Variablen und Zählvariablen als abhängige Variablen korrekt modelliert werden, ohne gegen deren inhärente Grenzen (z. B. den Wertebereich 0 bis 1 bei Wahrscheinlichkeiten oder die Nicht-Negativität bei Zähldaten) zu verstoßen.
Diese größere Flexibilität hat allerdings ihren Preis: Verallgemeinerte lineare Modelle sind komplexer als die lineare Regression, und ihre Koeffizienten lassen sich weniger intuitiv interpretieren – man benötigt zusätzliche Schritte wie die Geteilt-durch-4-Regel, Average Counterfactual Comparisons oder das Exponenzieren von Koeffizienten zu Odds Ratios bzw. Wachstumsfaktoren/-raten. Bei der Wahl des passenden Modells stellt sich daher immer die Abwägungsfrage: Wann ist die einfachere lineare Regression noch gut genug als Annäherung, und wann sollte man auf ein verallgemeinertes lineares Modell zurückgreifen? Bei sehr starken, deutlich nichtlinearen Zusammenhängen (wie im simulierten Beispiel) sind die Unterschiede gravierend; bei moderaten, realistischen Zusammenhängen (wie in den empirischen Beispielen) fallen sie hingegen oft klein aus.
Wichtig ist außerdem, dass sich vieles, was in den vorherigen Einheiten der Vorlesung zur linearen Regression erarbeitet wurde, direkt auf verallgemeinerte lineare Modelle übertragen lässt. Dazu zählen unter anderem kausale Annahmen, die Konzepte von Moderation und Interaktion, Mediation, Pfadmodelle, Strukturgleichungsmodelle sowie Mehrebenenmodelle. Verallgemeinerte lineare Modelle sind also keine völlig neue Modellklasse, sondern eine konsequente Erweiterung der aus der gesamten Vorlesung bekannten linearen Modelllogik auf andere Verteilungen der abhängigen Variable.
comments_count deskriptiv