set.seed(2025)
# # Wir ziehen 50 Xi1- und Xi2-Werte mit Var = 1
Xi1 <- rnorm(50, 0, 1)
Xi2 <- rnorm(50, 0, 1)Was ist das Incidental-Parameters-Problem?
Grundproblem: Pro Person nur wenige Messungen \(\rightarrow\) für jede Person muss ein eigener Mittelwert geschätzt werden \(\rightarrow\) je mehr Personen, desto mehr zusätzliche Parameter \(\rightarrow\) die Varianzschätzung (innerhalb einer Person) bleibt verzerrt, egal wie groß die Stichprobe wird.
Auf das Problem gestoßen bin ich während der Lektüre von Stigler (2016). Beschrieben wurde es von Neyman und Scott (1948). Da ich dieses Phänomen bis dahin nicht kannte, versuchen wir es im Folgenden mithilfe einer Simulation nachzuvollziehen.
Beschreibung des Problems aus praktischer Sicht
Wir erheben Daten von \(N=50\) Personen in einem Test-Retest-Setting. Es gibt also zwei Beobachtungen (Messungen) pro Person, aber die Anzahl der Personen kann theoretisch gegen unendlich gehen. Da wir davon ausgehen, dass der Messfehler bei allen Personen gleich ist, ist auch die Variation der Messwerte \(\sigma^2\) konstant (eine Annahme, die auch der KTT zugrunde liegt). Diese Varianz ist zufällig und wirkt bei allen Personen gleichermaßen. Schätzen wir die Varianz \(\hat{\sigma}_{ML}^2\) nach der Maximum-Likelihood-Methode (ML), erhalten wir als Schätzer einen Wert, der nur halb so groß ist wie das “wahre” \(\sigma^2\).
Wir versuchen nun mathematisch nachzuvollziehen, warum das Problem zustande kommt. Dabei beschränken wir uns auf den Fall von genau zwei Beobachtungen pro Einheit.
Schritt 1: Das Modell
Für jede Einheit \(i = 1,\dots,N\) werden zwei Beobachtungen erfasst:
\[X_{i1}, X_{i2} \sim \mathcal{N}(\mu_i, \sigma^2)\]
Dabei ist \(\sigma^2\) der interessante Parameter. Die \(\mu_i\) sind sogenannte beiläufige (incidental) Parameter. Es gilt: \(N \to \infty\), aber pro Einheit bleibt es bei nur 2 Beobachtungen.
Schritt 2: Maximum-Likelihood-Schätzer
Der ML-Schätzer des Mittelwerts:
\(\hat{\mu}_i = \frac{X_{i1} + X_{i2}}{2}\)
Der ML-Schätzer für die Varianz:
\[\hat{\sigma}^2_{ML} = \frac{1}{2N} \sum_{i=1}^N \sum_{j=1}^2 (X_{ij} - \hat{\mu}_i)^2\]
# wir ziehen 50 Xi1 und Xi2 Werte mit Var = 1
mu_hat_i <- (Xi1 + Xi2) / 2
sigma_hat_i <- sum(
(Xi1 - mu_hat_i)^2 + (Xi2 - mu_hat_i)^2
) / (2 * 50)
sigma_hat_i[1] 0.5680844
Hier sieht man dass der Schätzer \(\hat{\sigma}^2_{ML} = 0.57\) zu klein ist.
Das ist das Problem – es folgt noch die mathematische Herleitung (zugegeben mit GPT 5), aber das ändert nichts an der Mathematik.
Schritt 3: Fehlerdarstellung
\(X_{ij} = \mu_i + \varepsilon_{ij}\) mit \(\varepsilon_{ij} \sim \mathcal{N}(0,\sigma^2)\)
Dann:
\(\hat{\mu}_i = \mu_i + \frac{\varepsilon_{i1} + \varepsilon_{i2}}{2}\)
Schätzfehler von \(\mu_i\):
\(\hat{\mu}_i - \mu_i = \frac{\varepsilon_{i1} + \varepsilon_{i2}}{2}\)
Residuen:
\(X_{i1} - \hat{\mu}_i = \frac{\varepsilon_{i1} - \varepsilon_{i2}}{2}\)
\(X_{i2} - \hat{\mu}_i = -\frac{\varepsilon_{i1} - \varepsilon_{i2}}{2}\)
Schritt 4: Summe quadrierter Residuen
\(\sum_{j=1}^2 (X_{ij} - \hat{\mu}_i)^2 = \frac{(\varepsilon_{i1} - \varepsilon_{i2})^2}{2}\)
Einsetzen in den Varianzschätzer:
\(\hat{\sigma}^2_{ML} = \frac{1}{4N} \sum_{i=1}^N (\varepsilon_{i1} - \varepsilon_{i2})^2\)
Schritt 5: Erwartungswert
\(\mathrm{Var}(\varepsilon_{i1}) = \mathrm{Var}(\varepsilon_{i2}) = \sigma^2\)
\(\mathrm{Var}(\varepsilon_{i1} - \varepsilon_{i2}) = \sigma^2 + \sigma^2 = 2\sigma^2\)
Erwartungswert:
\(\mathbb{E}[(\varepsilon_{i1} - \varepsilon_{i2})^2] = 2\sigma^2\)
Damit:
\(\mathbb{E}[\hat{\sigma}^2_{ML}] = \frac{1}{4N} \cdot N \cdot 2\sigma^2 = \frac{\sigma^2}{2}\)
Der ML-Schätzer unterschätzt also systematisch die wahre Varianz.
Schritt 6: Inkonsistenz
Nach dem Gesetz der großen Zahlen:
\(\hat{\sigma}^2_{ML} \xrightarrow{P} \frac{\sigma^2}{2} \neq \sigma^2\)
Der ML-Schätzer ist nicht konsistent.
Schritt 7: Korrektur
Ein konsistenter Schätzer wäre:
\(\tilde{\sigma}^2 = 2 \cdot \hat{\sigma}^2_{ML}\)
Dann gilt:
\(\tilde{\sigma}^2 \xrightarrow{P} \sigma^2\)
Simulation
Diese Simulation demonstriert das Neyman-Scott-Problem.
set.seed(2025)
n_pairs <- 50 # Anzahl der Paare (n im Text)
sigma2 <- 1 # wahre Varianz σ^2
n_sim <- 10000 # Anzahl der Simulationen
## Funktion für eine Simulation
get_ml_sigma <- function() {
### Ziehe n verschiedene Mittelwerte mu (könnten auch alle 0 sein);
### sd ist egal
mu <- rnorm(n_pairs, mean = 0, sd = 2)
## Schritt 1
### Datenpaare (Xi1, Xi2) mit der Wurzel aus sigma2 als SD
Xi1 <- rnorm(n_pairs, mean = mu, sd = sqrt(sigma2))
Xi2 <- rnorm(n_pairs, mean = mu, sd = sqrt(sigma2))
## Schritt 2
### Maximum-Likelihood-Schätzer der Mittelwerte:
mu_hat <- (Xi1 + Xi2) / 2
### Maximum-Likelihood-Schätzer für σ^2 aus dem Text:
sum(((Xi1 - mu_hat)^2 + (Xi2 -mu_hat)^2))/(2*n_pairs)
### Das ist die Formel aus Schritt 4
sigma2_ML <- mean((Xi1 - Xi2)^2) / 4
# Unverzerrter Schätzer für σ^2 bei 2 Beobachtungen pro Paar:
sigma2_unbiased <- mean((Xi1 - Xi2)^2) / 2
c(sigma2_ML = sigma2_ML, sigma2_unbiased = sigma2_unbiased)
}
## Viele Wiederholungen
res <- replicate(n_sim, get_ml_sigma())
res <- t(res) # Spalten: sigma2_ML, sigma2_unbiased
## Mittlere Schätzwerte über alle Simulationen
colMeans(res) sigma2_ML sigma2_unbiased
0.5000341 1.0000682
Keine Auswirkung auf statistische Tests
Obwohl es so aussieht, als müsste die Within-ANOVA individuelle Personenmittelwerte schätzen, werden diese nicht zur Schätzung der Fehlervarianz für den Zeiteffekt herangezogen. Die Schätzung des Within-Fehlers bleibt unverzerrt und ist nicht vom Incidental-Parameters-Problem betroffen. Ebenso ist es in linearen Mehrebenen-Modellen: Hier werden die Personenunterschiede typischerweise als Realisationen aus einer gemeinsamen Verteilung modelliert, etwa über einen zufälligen Intercept. Dadurch bleibt die Anzahl der Varianzparameter konstant, und das Neyman-Scott-Problem tritt in dieser Form nicht auf.
Literatur
Stigler, S. M. (2016). The seven pillars of statistical wisdom. Harvard University Press.
Neyman, J., & Scott, E. L. (1948). Consistent Estimates Based on Partially Consistent Observations. Econometrica, 16(1), 1–32. https://doi.org/10.2307/1914288