---
fontsize: 8pt
format:
  beamer:
    include-in-header: ../ALM_Header.tex
bibliography: ../ALM_Referenzen.bib
---


# {.plain}
<!-- Vorlesungstitel -->
\center
```{r, echo = F, out.width = "20%"}
knitr::include_graphics("../OvGU_Logo.png")
```

\vspace{2mm}
\huge
Allgemeines Lineares Modell

\vspace{6mm}
\large
BSc Psychologie, SoSe 2026

\vspace{5mm}
Joram Soch


<!-- Sitzung 11: Einfaktorielle Varianzanalyse -->
# {.plain}

\vfill
\center
\huge
\textcolor{black}{(11) Einfaktorielle Varianzanalyse}
\vfill


<!-- Inhaltsverzeichnis -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

Modellformulierung

Modellschätzung

Modellevaluation

Anwendung/Praxis

Selbstkontrollfragen

\vfill


<!-- Abschnitt 1: Anwendungsszenario -->
#

\large
\setstretch{2.5}
\vfill

**Anwendungsszenario**

Modellformulierung

Modellschätzung

Modellevaluation

Anwendung/Praxis

Selbstkontrollfragen

\vfill


<!-- Anwendungsszenario -->
# Anwendungsszenario

\textbf{\textcolor{darkblue}{Zwei oder mehr Gruppen}} (Stichproben) randomisierter experimenteller Einheiten.

Annahme der unabhängigen und identischen Normalverteilung $N(\mu_i,\sigma^2)$ der Daten.

$\mu_i$ und $\sigma^2$ unbekannt.

Inferentieller Test der Nullhypothese identischer Gruppenerwartungswerte beabsichtigt.

\center
\vspace{1mm}
$\Downarrow$
\vspace{1mm} 

Generalisierung des Zweistichproben-T-Tests bei unabhängigen Stichproben \linebreak
mit einfacher Nullhypothese für mehr als zwei Stichproben
\vfill


<!-- Anwendungsbeispiel -->
# Anwendungsszenario

\textcolor{darkblue}{Anwendungsbeispiel}

```{r, echo = F, out.width = "70%", fig.align = "center"}
knitr::include_graphics("Abbildungen/beispielstudie.pdf")
```

\small
\vspace{2mm}
Pre-Post-BDI-Differenzwertanalyse für drei Gruppen von Patient:innen (F2F, ONL, WLC):
\vspace{-2mm}

* inferentielle Evidenz für Gruppenerwartungswertunterschiede?
* empirische Evidenz für unterschiedliche Therapiewirksamkeit?

```{r, echo = F, eval = F}
# Initialisierung
set.seed(1)                                                                      # Zufallszahlengenerator initialisieren

# Simulationsparameter
n_c    = 3                                                                       # Anzahl Gruppen
n_i    = 40                                                                      # Anzahl Proband:innen pro Gruppe
n      = n_c*n_i                                                                 # Gesamtanzahl Datenpunte

# Simulationsparameter
mu_pre = c(31, 32, 29)                                                           # \mu pre F2F, ONL, WLC
mu_pos = c(27, 26, 30)                                                           # \mu post F2F, ONL, WLC
sigsqr = 10                                                                      # \sigma^2

# Datensimulation
D            = data.frame("ID" = 1:n)                                            # Dataframe-Initialisierung und ID-Variable
D$Setting    = c(rep("F2F",n_i), rep("ONL", n_i), rep("WLC", n_i))               # Bedingung
D$PreBDI     = c(round(rnorm(n_i, mu_pre[1], sqrt(sigsqr))),                     # PreBDI
                 round(rnorm(n_i, mu_pre[2], sqrt(sigsqr))),
                 round(rnorm(n_i, mu_pre[3], sqrt(sigsqr))))
D$PostBDI    = c(round(rnorm(n_i, mu_pos[1], sqrt(sigsqr))),                     # PostBDI
                 round(rnorm(n_i, mu_pos[2], sqrt(sigsqr))),
                 round(rnorm(n_i, mu_pos[3], sqrt(sigsqr))))
D$dBDI       = -(D$PostBDI - D$PreBDI)                                           # -(PostBDI - PreBDI) = PreBDI - PostBDI

# Datenspeicherung
write.csv(D,
    file     = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv")
```


<!-- Abschnitt 2: Modellformulierung -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

**Modellformulierung**

Modellschätzung

Modellevaluation

Anwendung/Praxis

Selbstkontrollfragen

\vfill


<!-- Definition: EVA-Modell in Erwartungswertparameterdarstellung -->
# Modellformulierung

\footnotesize
\begin{definition}[EVA-Modell in Erwartungswertparameterdarstellung]

\justifying
$y_{ij}$ mit $i= 1,...,p$, welche die Gruppen indizieren, und $j = 1,...,n_i$, welche die experimentellen Einheiten innerhalb der Gruppen indizieren, seien Zufallsvariablen, die die Datenpunkte eines Einfaktoriellen-Varianzanalyse-Szenarios (EVA) modellieren. Dann hat das \textit{EVA-Modell in Erwartungswertparameterdarstellung} die strukturelle Form
\begin{equation}
y_{ij} = \mu_i + \varepsilon_{ij}
\quad \mbox{mit} \quad
\varepsilon_{ij} \sim N(0,\sigma^2)
\quad \mbox{u.i.v. für} \quad
i = 1,...,p, \; j = 1,...,n_i
\quad \mbox{mit} \quad
\mu_i \in \mathbb{R}, \; \sigma^2 > 0
\end{equation}

und die Datenverteilungsform
\begin{equation} \label{eq:eva_klassisch}
y_{ij} \sim N(\mu_i,\sigma^2)
\quad \mbox{u.v. für} \quad
i = 1,...,p, \; j = 1,...,n_i
\quad \mbox{mit} \quad
\mu_i \in \mathbb{R}, \; \sigma^2 > 0.
\end{equation}

Wenn $n_i := m$ für alle $i = 1,...,p$, heißt das EVA-Szenario \textit{balanciert}.

\end{definition}

Bemerkungen

* Die Äquivalenz der beiden Modellformulierungen folgt aus den Ergebnissen zu Transformationen der Normalverteilung (vgl. Einheit (8) in *Wahrscheinlichkeitstheorie und Frequentistische Inferenz* und Einheit (4) in *Allgemeines Lineares Modell*).
* Es handelt sich um die Generalisierung des Zweistichproben-T-Test-Modells für unabhängige Stichproben unter Annahme identischer Varianzparameter von $p = 2$ auf ein beliebiges $p \in \mathbb{N}$.
* Bei balancierten Varianzanalyseszenarien besteht jede Datengruppe aus der gleichen Anzahl von Datenpunkten.


<!-- Motivation der Effektdarstellung -->
# Modellformulierung

\textcolor{darkblue}{Motivation der Effektdarstellung}

\footnotesize
Die Erwartungswertparameterdarstellung des EVA Modells ist ein valides ALM, das sich in dieser Form auch in der Literatur findet (z.B. @georgii2009, Kapitel 12.4). Im Sinne der Konsistenz mit mehrfaktoriellen Varianzanalysemodellen bietet sich jedoch eine Reparameterisierung der Erwartungswertparameter an. Kern dieser Reparameterisierung ist es, den Erwartungswertparameter der $i$ten Gruppe als Summe eines *gruppenübergreifenden Erwartungswertparameters* $\mu_0 \in \mathbb{R}$ und eines *gruppenspezifischen Effektparameters* $\alpha_i \in \mathbb{R}$ zu modellieren,
\begin{equation} \label{eq:uber}
\mu_i := \mu_0 + \alpha_i
\quad \mbox{für} \quad
i = 1,...,p.
\end{equation}

\vspace{-2mm}
Dabei modelliert $\alpha_i$ den Unterschied (die Differenz) zwischen dem $i$ten Erwartungswertparameter $\mu_i$ und dem gruppenübergreifenden Erwartungswertparameter $\mu_0$,
\begin{equation}
\alpha_i = \mu_i - \mu_0
\quad \mbox{für} \quad
i = 1,...,p.
\end{equation}

\vspace{-2mm}
Allerdings hat die in dieser Form vorgenommene Reparameterisierung einen entscheidenen Nachteil: es werden $p$ Erwartungswertparameter $\mu_i, i = 1,...,p$ durch die $p + 1$ Parameter $\mu_0$ und $\alpha_i, i = 1,...,p$ dargestellt. Diese Darstellung ist im Allgemeinen nicht eindeutig. Zum Beispiel können die Erwartungswertparameter
\begin{equation}
\mu_1 = 3, \; \mu_2 = 5, \; \mu_3 = 6
\end{equation}
sowohl durch
\begin{equation}
\mu_0 = 0
\quad \mbox{und} \quad
\alpha_1 = 3, \; \alpha_2 = 5, \; \alpha_3 = 6
\end{equation}
als auch durch
\begin{equation}
\mu_0 = 1
\quad \mbox{und} \quad
\alpha_1 = 2, \; \alpha_2 = 4, \; \alpha_3 = 5
\end{equation}
dargestellt werden. Man sagt in diesem Kontext auch, dass das EVA-Modell mit \eqref{eq:uber} *überparameterisiert* ist.

<!-- Zum Beispiel können die Erwartungswertparameter $\mu_1 = 3, \mu_2 = 5, \mu_3 = 6$ sowohl durch den gruppenspezifischen Erwartungswertparameter $\mu_0 = 0$ und die gruppenunspezifischen Effektparameter $\alpha_1 = 3, \alpha_2 = 5, \alpha_3 = 6$ als auch durch den gruppenunspezifischen Erwartungswertparameter $\mu_0 = 1$ und die gruppenspezifischen Effektparameter $\alpha_1 = 2, \alpha_2 = 4, \alpha_3 = 5$ dargestellt werden. Man sagt in diesem Kontext auch, dass das EVA Modell mit \eqref{eq:uber} *überparameterisiert* ist. -->


<!-- Motivation der Effektdarstellung -->
# Modellformulierung
\textcolor{darkblue}{Motivation der Effektdarstellung}

\footnotesize
Datenanalytisch hat die Überparameterisierung eines Varianzanalysemodells den Nachteil, dass aus $p$ geschätzten Erwartungswertparametern $p + 1$ Betaparameterschätzer bestimmt werden müssten, was wie gesehen nicht eindeutig erfolgen kann. Um diese Probleme in der Effektparameterdarstellung des EVA-Modells zu umgehen und sie konsistent auf mehrfaktorielle Varianzanalysemodelle zu übertragen, bietet sich die Einführung einer Nebenbedingung an:
\begin{equation}
\alpha_1 := 0 \; .
\end{equation}

Es wird also ein Effektparameter von vornherein als "gleich Null" angenommen. Für die gruppenspezifischen Erwartungswertparameter ergibt sich damit
\begin{align}
\begin{split}
\mu_1 & := \mu_0 \\
\mu_i & := \mu_0 + \alpha_i
\quad \mbox{für} \quad
i = 2,...,p.
\end{split}
\end{align}

Hierbei wird die erste Gruppe nun als *Referenzgruppe* bezeichnet und die $\alpha_i$ modellieren die Differenz zwischen dem Erwartungswertparameter der $i$ten Gruppe und dem Erwartungswertparameter der ersten Gruppe:
\begin{equation}
\alpha_i = \mu_i - \mu_0 = \mu_i - \mu_1
\quad \mbox{für} \quad
i = 1,...,p.
\end{equation}

$\mu_0$ ist also kein gruppenübergreifender Erwartungswertparameter mehr, sondern identisch mit dem Erwartungswertparameter der ersten Gruppe. Welche tatsächliche experimentelle Gruppe dabei als "erste Gruppe" definiert wird, ist unerheblich. Entscheidend ist, dass die entsprechenden Erwartungswertparameterschätzer $\hat{\mu}_0$,$\hat{\alpha}_2, ..., \hat{\alpha}_p$ korrekt als (1) Erwartungswertparameterschätzer der Referenzgruppe ($\hat{\mu}_0$) und (2) geschätzte Differenz zwischen dem Erwartungswertparameter der Referenzgruppe und dem Erwartungswertparameter der $i$ten Gruppe ($\hat{\alpha}_2, ..., \hat{\alpha}_p$) verstanden werden. Wir formalisieren das oben Gesagte in folgendem Theorem.


<!-- Theorem: EVA-Modell in Effektdarstellung mit Referenzgruppe I -->
# Modellformulierung

\footnotesize
\begin{theorem}[EVA-Modell in Effektdarstellung mit Referenzgruppe I]

\justifying
\normalfont
Gegeben sei das EVA-Modell in Erwartungswertparameterdarstellung. Dann können die Zufallsvariablen, die die Datenpunkte des EVA-Szenarios modellieren, äquivalent in der strukturellen Form
\begin{align} \label{eq:eva_effekt_1}
\begin{split}
y_{1j} & = \mu_0 + \varepsilon_{1j} \hphantom{+ \alpha_i}
\quad \mbox{mit} \quad
\varepsilon_{1j} \sim N(0,\sigma^2)
\quad \mbox{u.i.v. für} \quad
j = 1,...,n_1 \\
y_{ij} & = \mu_0 + \alpha_i + \varepsilon_{ij}
\quad \mbox{mit} \quad
\varepsilon_{ij} \sim N(0,\sigma^2)
\quad \mbox{u.i.v. für} \quad
i = 2,...,p, \; j = 1,...,n_i
\end{split}
\end{align}

mit $\alpha_i := \mu_i - \mu_1$ für $i = 2,...,p$ und in der entsprechenden Datenverteilungsform
\begin{align}\label{eq:eva_effekt_2}
\begin{split}
y_{1j} & \sim N(\mu_0,\sigma^2) \hphantom{+ \alpha_i} 
\quad \mbox{u.i.v. für} \quad
j = 1,...,n_i
\quad \mbox{mit} \quad
\mu_0 \in \mathbb{R}, \; \sigma^2 > 0 \\
y_{ij} & \sim N(\mu_0 + \alpha_i,\sigma^2)
\quad \mbox{u.i.v. für} \quad
i = 2,..., p, \; j = 1,...,n_i
\quad \mbox{mit} \quad
\alpha_i \in \mathbb{R}, \; \sigma^2 > 0 \\
\end{split}
\end{align}

geschrieben werden. Wir nennen \eqref{eq:eva_effekt_1} und \eqref{eq:eva_effekt_2} die strukturelle und die Datenverteilungsform des \textit{EVA-Modells in Effektdarstellung mit Referenzgruppe}.

\end{theorem}

\footnotesize
\underline{Beweis}

Es gilt
\begin{equation}
\mu_i = \mu_0 + \mu_i - \mu_0 \; .
\end{equation}

Die Parameterisierungen mit $\mu_i$ und mit $\mu_0 + \mu_i - \mu_0$ sind also gleich und damit äquivalent. Dann folgt aber auch
\begin{equation}
\mu_i = \mu_0 + (\mu_i - \mu_0) =: \mu_0 + \alpha_i
\quad \mbox{für} \quad
i = 1,...,p \; .
\end{equation}

Mit $\alpha_1 := 0$ gilt dann $\mu_1 = \mu_0$ und  $\mu_i = \mu_0 + \alpha_i$ für $i = 2,...,p$, wie im Theorem behauptet.
$\hfill\Box$


<!-- Theorem: EVA-Modell in Effektdarstellung mit Referenzgruppe II -->
# Modellformulierung

\vspace{2mm}
\footnotesize
\begin{theorem}[EVA-Modell in Effektdarstellung mit Referenzgruppe II]

\justifying
\normalfont
Gegeben sei die strukturelle Form des EVA-Modells in Effektdarstellung mit Referenzgruppe. Dann hat dieses Modell die Designmatrixform
\begin{equation}
y = X\beta + \varepsilon
\quad \mbox{mit} \quad
\varepsilon \sim N(0_n,\sigma^2 I_n),
\end{equation}

wobei $n := \sum_{i=1}^p n_i$ und
\begin{equation*}
\renewcommand{\arraystretch}{0.8}
y :=
\begin{pmatrix*}[c]
y_{11} \\ \vdots \\ y_{1n_1} \\
y_{21} \\ \vdots \\ y_{2n_2} \\
       \\ \vdots \\          \\
y_{p1} \\ \vdots \\ y_{pn_p}
\end{pmatrix*}
\in \mathbb{R}^n, \quad
X :=
\begin{pmatrix}
1      &  0      &        & 0      \\
\vdots &  \vdots & \cdots & \vdots \\
1      &  0      &        & 0      \\
1      &  1      &        & 0      \\
\vdots &  \vdots & \cdots & \vdots \\
1      &  1      &        & 0      \\
       &         &        &        \\
\vdots &  \vdots & \cdots & \vdots \\
       &         &        &        \\
1      &  0      &        & 1      \\
\vdots &  \vdots & \cdots & \vdots \\
1      &  0      &        & 1      \\
\end{pmatrix}
\in \mathbb{R}^{n \times p}, \quad
\beta :=
\begin{pmatrix}
\mu_0    \\
\alpha_2 \\
\vdots   \\
\alpha_p
\end{pmatrix}
\in \mathbb{R}^{p} \quad
\mbox{und} \quad
\sigma^2 > 0.
\end{equation*}

\end{theorem}

\vspace{-2mm}
\underline{Beweis}

Das Theorem ergibt sich direkt mit den Regeln der Matrixmultiplikation.
$\hfill\Box$


<!-- Designmatrix des Modells -->
# Modellformulierung

Designmatrix des EVA-Modells in Effektdarstellung mit Referenzgruppe ($n = 12$, $p = 4$)

```{r, echo = F, eval = F}
# Designmatrixerzeugung
m      = 3                                                # Anzahl Datenpunkte i-te Gruppe
p      = 4                                                # Anzahl Gruppen
n      = p*m                                              # Gesamtanzahl Datenpunkte
Xt     = cbind(matrix(1, nrow = n, ncol = 1),             # n x p Designmatrix
               kronecker(diag(p), matrix(1, nrow = m, ncol = 1)))
X      = Xt[,-2]                                          # eliminiere zweite Spalte
Xp     = X

# Abbildungsparameter
library(plot.matrix)
graphics.off()
dev.new()
par(
    family    = "sans",
    mar       = c(1,1,1,1),
    pty       = "s",
    bty       = "l",
    lwd       = 1,
    las       = 1,
    mgp       = c(2,1,0),
    xaxs      = "i",
    yaxs      = "i",
    font.main = 1,
    cex       = 1,
    cex.main  = 2)

# Designmatrix
plot(Xp,
    col       = gray(seq(0, 1, length.out = 256)),
    border    = "black",
    xlab      = '',
    ylab      = '',
    main      = '',
    breaks    = seq(-2, 1, length.out = 257),
    key       = NULL,
    axis.col  = NULL,
    axis.row  = NULL,
    asp       = 1)

# Speichern
dev.copy2pdf(
    file      = "Abbildungen/eva_designmatrix.pdf",
    width     = 10,
    height    = 10)
```

\vspace{2mm}
```{r, echo = F, out.width = "60%", fig.align = "center"}
knitr::include_graphics("Abbildungen/eva_designmatrix.pdf")
```


<!-- Beispiel -->
# Modellformulierung

\textcolor{darkblue}{Beispiel}
\small

Es seien
\begin{equation}
p = 4
\quad \mbox{und} \quad
n_i := 3
\quad \mbox{für} \quad
i = 1,...,p,
\quad \mbox{also} \quad
n = 12.
\end{equation}

Dann gilt
\begin{equation}
y =  X\beta + \varepsilon
\quad \mbox{mit} \quad
\varepsilon \sim N(0_{12},\sigma^2 I_{12})
\end{equation}

mit
\footnotesize
\begin{equation}
y := 
\begin{pmatrix}
y_{11} \\
y_{12} \\
y_{13} \\
y_{21} \\
y_{22} \\
y_{23} \\
y_{31} \\
y_{32} \\
y_{33} \\
y_{41} \\
y_{42} \\
y_{43}
\end{pmatrix}, \quad
X :=
\begin{pmatrix}
1  &  0  &  0  &  0  \\
1  &  0  &  0  &  0  \\
1  &  0  &  0  &  0  \\
1  &  1  &  0  &  0  \\
1  &  1  &  0  &  0  \\
1  &  1  &  0  &  0  \\
1  &  0  &  1  &  0  \\
1  &  0  &  1  &  0  \\
1  &  0  &  1  &  0  \\
1  &  0  &  0  &  1  \\
1  &  0  &  0  &  1  \\
1  &  0  &  0  &  1  \\
\end{pmatrix}
\in \mathbb{R}^{12 \times 4}, \quad
\beta :=
\begin{pmatrix}
\mu_0    \\
\alpha_2 \\
\alpha_3 \\
\alpha_4 \\
\end{pmatrix}
\in \mathbb{R}^{4} \quad
\mbox{und} \quad
\sigma^2 > 0.
\end{equation}


<!-- Beispiel in R -->
# Modellformulierung

\vspace{2mm}
\textcolor{darkblue}{Beispiel}
\vspace{1mm}

\tiny
\setstretch{1}
```{r, echo = T}
# Modellformulierung
library(MASS)                                             # multivariate Normalverteilung
m      = 3                                                # Anzahl Datenpunkte i-te Gruppe
p      = 4                                                # Anzahl Gruppen
n      = p*m                                              # Gesamtanzahl Datenpunkte
Xt     = cbind(matrix(1, nrow = n, ncol = 1),             # n x p Designmatrix
               kronecker(diag(p), matrix(1, nrow = m, ncol = 1)))
X      = Xt[,-2]                                          # eliminiere zweite Spalte
I_n    = diag(n)                                          # n x n Einheitsmatrix
beta   = matrix(c(2,0,1,-1), nrow = p)                    # \beta = (\mu_0,\alpha_2,\alpha_3,\alpha_4)^T
sigsqr = 5                                                # \sigma^2

# Datenrealisierung
y      = mvrnorm(1, X %*% beta, sigsqr*I_n)               # eine Realisierung des n-dimensionalen ZVs y
print(X)
```


<!-- Abschnitt 3: Modellschätzung -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

Modellformulierung

**Modellschätzung**

Modellevaluation

Anwendung/Praxis

Selbstkontrollfragen

\vfill


<!-- Theorem: Betaparameterschätzer im EVA-Modell -->
# Modellschätzung

\footnotesize
\begin{theorem}[Betaparameterschätzer im EVA-Modell]

\justifying
\normalfont
Gegeben sei die Designmatrixform der EVA in Effektdarstellung mit Referenzgruppe. Dann ergibt sich für den Betaparameterschätzer
\begin{equation}
\hat{\beta} :=
\begin{pmatrix} \hat{\mu}_0 \\ \hat{\alpha}_2        \\ \vdots \\ \hat{\alpha}_p                      \end{pmatrix} =
\begin{pmatrix} \bar{y}_1   \\ \bar{y}_2 - \bar{y}_1 \\ \vdots \\ \bar{y}_p - \bar{y}_1 \end{pmatrix} \; ,
\end{equation}

wobei
\begin{equation}
\bar{y}_i := \frac{1}{n_i}\sum_{j=1}^{n_i}y_{ij}
\end{equation}

das Stichprobenmittel der $i$ten Gruppe bezeichent.

\end{theorem}

Bemerkungen

* \justifying Der Erwartungswertparameter der ersten Gruppe wird durch das Stichprobenmittel der ersten Gruppe geschätzt; der Effektparameter der $i$ten Gruppe wird durch die Differenz des Stichprobenmittels der $i$ten und der ersten Gruppe geschätzt.


<!-- Beweis -->
# Modellschätzung

\vspace{2mm}
\small
\underline{Beweis}

Wir halten zunächst fest, dass
\tiny
\begin{align*}
\begin{split}
X^\mathrm{T} X
& =
\setcounter{MaxMatrixCols}{20}
\begin{pmatrix}
1 & \cdots & 1 & 1 & \cdots & 1 & & \cdots & & 1 & \cdots & 1 \\
0 & \cdots & 0 & 1 & \cdots & 1 & & \cdots & & 0 & \cdots & 0 \\
  & \vdots &   &   & \vdots &   & & \vdots & &   & \vdots &   \\
0 & \cdots & 0 & 0 & \cdots & 0 & & \cdots & & 1 & \cdots & 1 \\
\end{pmatrix}
\begin{pmatrix}
1       &   0        &          &   0       \\
\vdots  &    \vdots  &  \cdots  &   \vdots  \\
1       &   0        &          &   0       \\
1       &   1        &          &   0       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   1        &          &   0       \\
        &            &          &           \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
        &            &          &           \\
1       &   0        &          &   1       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   0        &          &   1       \\
\end{pmatrix} \\
& =
\begin{pmatrix}
n      & n_2    & n_3    & \cdots  & n_p    \\
n_2    & n_2    & 0      & \cdots  & 0      \\
n_3    & 0      & n_3    & \cdots  & 0      \\
\vdots & \vdots & \vdots & \ddots  & \vdots \\
n_p    & 0      & 0      & \cdots  & n_p    \\
\end{pmatrix}.
\end{split}
\end{align*}


<!-- Beweis -->
# Modellschätzung

\small
\underline{Beweis (fortgeführt)}

Die Inverse von $X^\mathrm{T}X$ ist
\begin{equation}
\renewcommand{\arraystretch}{1.5}
(X^\mathrm{T}X)^{-1} =
\begin{pmatrix*}[c]
  \frac{1}{n_1}
& -\frac{1}{n_1}
& \cdots
& -\frac{1}{n_1} \\
  -\frac{1}{n_1}
& \frac{n_1+n_2}{n_1 n_2}
& \cdots
& \frac{1}{n_1}  \\
  \vdots
& \vdots
& \ddots
& \vdots         \\
  -\frac{1}{n_1}
& \frac{1}{n_1}
& \cdots
& \frac{n_1+n_p}{n_1 n_p}
\end{pmatrix*}.
\end{equation}

So gilt zum Beispiel für $p = 3$, dass
\begin{equation}
X^\mathrm{T} X =
\begin{pmatrix}
n   & n_2 & n_3 \\
n_2 & n_2 & 0   \\
n_3 & 0   & n_3 \\
\end{pmatrix}
\quad \mbox{und} \quad
\renewcommand{\arraystretch}{1.5}
(X^\mathrm{T}X)^{-1} =
\begin{pmatrix*}[r]
  \frac{1}{n_1}
& -\frac{1}{n_1}
& -\frac{1}{n_1} \\
  -\frac{1}{n_1}
& \frac{n_1+n_2}{n_1 n_2}
& \frac{1}{n_1}  \\
  -\frac{1}{n_1}
& \frac{1}{n_1}
& \frac{n_1+n_3}{n_1 n_3}
\end{pmatrix*}.
\end{equation}


<!-- Beweis -->
# Modellschätzung

\small
\underline{Beweis (fortgeführt)}

Wir halten weiterhin fest, dass
\tiny
\begin{align}
\begin{split}
\renewcommand{\arraystretch}{1}
X^\mathrm{T}y
& =
\setcounter{MaxMatrixCols}{20}
\begin{pmatrix}
1 & \cdots & 1 & 1 & \cdots & 1 & & \cdots & & 1 & \cdots & 1 \\
0 & \cdots & 0 & 1 & \cdots & 1 & & \cdots & & 0 & \cdots & 0 \\
  & \vdots &   &   & \vdots &   & & \vdots & &   & \vdots &   \\
0 & \cdots & 0 & 0 & \cdots & 0 & & \cdots & & 1 & \cdots & 1 \\
\end{pmatrix}
\begin{pmatrix*}[c]
y_{11}  \\  \vdots  \\  y_{1n_1} \\
y_{21}  \\  \vdots  \\  y_{2n_2} \\  
        \\  \vdots  \\           \\
y_{p1}  \\  \vdots  \\  y_{pn_p} \\
\end{pmatrix*} \\
& =
\begin{pmatrix}
\sum_{i=1}^p \sum_{j=1}^{n_i} y_{ij} \\
\sum_{j=1}^{n_2} y_{2j}              \\
\vdots                               \\
\sum_{j=1}^{n_p} y_{pj}              \\
\end{pmatrix}.
\end{split}
\end{align}


<!-- Beweis -->
# Modellschätzung

\vspace{2mm}
\small
\underline{Beweis (fortgeführt)}

Es ergibt sich also
\tiny
\begin{equation}
\renewcommand{\arraystretch}{2}
\hat{\beta}
= (X^\mathrm{T} X)^{-1}X^\mathrm{T}y
=
\begin{pmatrix*}[c]
   \frac{1}{n_1}
& -\frac{1}{n_1}
& \cdots
& -\frac{1}{n_1} \\
  -\frac{1}{n_1}
& \frac{n_1+n_2}{n_1 n_2}
& \cdots
& \frac{1}{n_1}  \\
  \vdots
& \vdots
& \ddots
& \vdots         \\
 -\frac{1}{n_1}
& \frac{1}{n_1}
& \cdots
& \frac{n_1+n_p}{n_1 n_p}
\end{pmatrix*}
\begin{pmatrix}
\sum_{i=1}^p \sum_{j=1}^{n_i} y_{ij} \\
\sum_{j=1}^{n_2} y_{2j}              \\
\vdots                               \\
\sum_{j=1}^{n_p} y_{pj}              \\
\end{pmatrix}.
\end{equation}
\footnotesize

\small
Für die erste Komponente von $\hat{\beta}$ ergibt sich damit
\tiny
\begin{align}
\begin{split}
\hat{\beta}_1
& =
\frac{1}{n_1}\sum_{i=1}^p \sum_{j=1}^{n_i} y_{ij}
- \frac{1}{n_1}\sum_{j=1}^{n_2} y_{2j}
- \cdots
- \frac{1}{n_1}\sum_{j=1}^{n_p} y_{pj} \\
& =
\frac{1}{n_1}
\left(
\left(
  \sum_{j=1}^{n_1} y_{1j}
+ \sum_{j=1}^{n_2} y_{2j}
+ \cdots
+ \sum_{j=1}^{n_p} y_{pj}
\right)
- \sum_{j=1}^{n_2} y_{2j}
- \cdots
- \sum_{j=1}^{n_p} y_{pj}
\right) \\
& = \frac{1}{n_1}\sum_{j=1}^{n_1} y_{1j} \\
& = \bar{y}_1.
\end{split}
\end{align}


<!-- Beweis -->
# Modellschätzung

\vspace{2mm}
\small
\underline{Beweis (fortgeführt)}

Für die zweite Komponente von $\hat{\beta}$ und analog für alle weiteren ergibt sich
\tiny
\begin{align}
\begin{split}
\hat{\beta}_2
& =
- \frac{1}{n_1}            \sum_{i=1}^p \sum_{j=1}^{n_i} y_{ij}
+ \frac{n_1 + n_2}{n_1n_2} \sum_{j=1}^{n_2} y_{2j}
+ \frac{1}{n_1}            \sum_{j=1}^{n_3} y_{3j}
+ \cdots
+ \frac{1}{n_1}            \sum_{j=1}^{n_p} y_{pj} \\
& =
\frac{n_1 + n_2}{n_1n_2} \sum_{j=1}^{n_2} y_{2j}
-\frac{1}{n_1}
\left(
\left(
  \sum_{j=1}^{n_1} y_{1j}
+ \sum_{j=1}^{n_2} y_{2j}
+ \cdots
+ \sum_{j=1}^{n_p} y_{pj}
\right)
- \sum_{j=1}^{n_3} y_{3j}
- \cdots
- \sum_{j=1}^{n_p} y_{pj}
\right) \\
& =
 \frac{n_1 + n_2}{n_1n_2} \sum_{j=1}^{n_2} y_{2j}
-\frac{1}{n_1}            \sum_{j=1}^{n_1} y_{1j}
-\frac{1}{n_1}            \sum_{j=1}^{n_2} y_{2j} \\
& =
 \frac{n_1+n_2}{n_1n_2} \sum_{j=1}^{n_2} y_{2j}
-\frac{n_2}{n_1n_2}     \sum_{j=1}^{n_2} y_{2j}
-\frac{1}{n_1}          \sum_{j=1}^{n_1} y_{1j} \\
& =
 \frac{n_1}{n_1n_2} \sum_{j=1}^{n_2} y_{2j}
-\frac{1}{n_1}      \sum_{j=1}^{n_1} y_{1j} \\
& =
 \frac{1}{n_2} \sum_{j=1}^{n_2} y_{2j}
-\frac{1}{n_1} \sum_{j=1}^{n_1} y_{1j}
= \bar{y}_2 - \bar{y}_1.
\end{split}
\end{align}
$\hfill\Box$


<!-- Theorem: Varianzparameterschätzer im EVA-Modell -->
# Modellschätzung

\footnotesize
\begin{theorem}[Varianzparameterschätzer im EVA-Modell]

\normalfont
\justifying
Gegeben sei die Designmatrixform der EVA in Effektdarstellung mit Referenzgruppe. Dann ergibt sich für den Varianzparameterschätzer
\begin{equation}
\hat{\sigma}^2
= \frac{\sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij} - \bar{y}_i)^2}{\sum_{i=1}^p n_i - p}
=: s_{1 \ldots p}^2
\end{equation}

wobei
\begin{equation}
\bar{y}_i := \frac{1}{n_i}\sum_{j=1}^{n_i}y_{ij}
\end{equation}

das Stichprobenmittel der $i$ten Gruppe bezeichent.

\end{theorem}

Bemerkungen

* \justifying Der Varianzparameterschätzer ergibt sich wiederum als gebündelte Stichprobenvarianz, generalisiert auf eine beliebige Anzahl von Gruppen $p$. Wir erinnern daran, dass die gebündelte Stichprobenvarianz $s_{1 \ldots p}^2$ im Allgemeinen nicht der Stichprobenvarianz des Datenvektors $s_y^2$ entspricht.


<!-- Beweis -->
# Modellschätzung

\footnotesize
\underline{Beweis}

Der Varianzparameterschätzer ist gegeben durch
\begin{align}
\begin{split}
\hat{\sigma}^{2}
& =
\frac{(y-X\hat{\beta})^\mathrm{T}(y - X\hat{\beta})}{n - p} \\
& =
\frac{1}{n - p}
\left(
\begin{pmatrix*}[c]
y_{11}  \\  \vdots  \\  y_{1n_1} \\
y_{21}  \\  \vdots  \\  y_{2n_2} \\  
        \\  \vdots  \\           \\
y_{p1}  \\  \vdots  \\  y_{pn_p} \\
\end{pmatrix*} -
\begin{pmatrix}
1       &   0        &          &   0       \\
\vdots  &    \vdots  &  \cdots  &   \vdots  \\
1       &   0        &          &   0       \\
1       &   1        &          &   0       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   1        &          &   0       \\
        &            &          &           \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
        &            &          &           \\
1       &   0        &          &   1       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   0        &          &   1       \\
\end{pmatrix}
\begin{pmatrix}
\bar{y}_1             \\
\bar{y}_2 - \bar{y}_1 \\
\vdots                \\
\bar{y}_p - \bar{y}_1
\end{pmatrix}
\right)^\mathrm{T}
\left(
\begin{pmatrix*}[c]
y_{11}  \\  \vdots  \\  y_{1n_1} \\
y_{21}  \\  \vdots  \\  y_{2n_2} \\  
        \\  \vdots  \\           \\
y_{p1}  \\  \vdots  \\  y_{pn_p} \\
\end{pmatrix*} -
\begin{pmatrix}
1       &   0        &          &   0       \\
\vdots  &    \vdots  &  \cdots  &   \vdots  \\
1       &   0        &          &   0       \\
1       &   1        &          &   0       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   1        &          &   0       \\
        &            &          &           \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
        &            &          &           \\
1       &   0        &          &   1       \\
\vdots  &   \vdots   &  \cdots  &   \vdots  \\
1       &   0        &          &   1       \\
\end{pmatrix}
\begin{pmatrix}
\bar{y}_1             \\
\bar{y}_2 - \bar{y}_1 \\
\vdots                \\
\bar{y}_p - \bar{y}_1
\end{pmatrix} \right)
\end{split}
\end{align}


<!-- Beweis -->
# Modellschätzung

\footnotesize
\underline{Beweis (fortgeführt)}

Es ergibt sich also
\begin{align}
\begin{split}
\hat{\sigma}^{2}
& =
\frac{1}{n_1+\ldots+n_p - p}
\begin{pmatrix*}[c]
y_{11}   - \bar{y}_1                           \\
\vdots                                         \\
y_{1n_1} - \bar{y}_1                           \\
y_{21}   - \bar{y}_1 - (\bar{y}_2 - \bar{y}_1) \\
\vdots                                         \\
y_{2n_2} - \bar{y}_1 - (\bar{y}_2 - \bar{y}_1) \\
                                               \\
\vdots                                         \\
                                               \\
y_{p1}   - \bar{y}_1 - (\bar{y}_p - \bar{y}_1) \\
\vdots                                         \\
y_{pn_p} - \bar{y}_1 - (\bar{y}_p - \bar{y}_1) \\
\end{pmatrix*}^\mathrm{T}
\begin{pmatrix*}[c]
y_{11}   - \bar{y}_1                           \\
\vdots                                         \\
y_{1n_1} - \bar{y}_1                           \\
y_{21}   - \bar{y}_1 - (\bar{y}_2 - \bar{y}_1) \\
\vdots                                         \\
y_{2n_2} - \bar{y}_1 - (\bar{y}_2 - \bar{y}_1) \\
                                               \\
\vdots                                         \\
                                               \\
y_{p1}   - \bar{y}_1 - (\bar{y}_p - \bar{y}_1) \\
\vdots                                         \\
y_{pn_p} - \bar{y}_1 - (\bar{y}_p - \bar{y}_1) \\
\end{pmatrix*} \\
& = \frac{\sum_{j=1}^{n_1} (y_{1j} - \bar{y}_1)^2 + \sum_{j=1}^{n_2} (y_{2j} - \bar{y}_2)^2 + \ldots + \sum_{j=1}^{n_p} (y_{pj} - \bar{y}_p)^2}{n_1+\ldots+n_p - p} \\
& = \frac{\sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij} - \bar{y}_i)^2}{\sum_{i=1}^p n_i - p} \\
& =: s_{1 \ldots p}^2.
\end{split}
\end{align}
$\hfill\Box$


<!-- Simulation -->
# Modellschätzung

\vspace{2mm}

\tiny
\setstretch{0.8}
```{r, echo = T}
# Daten einlesen
fname      = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D          = read.table(fname, sep = ",", header = TRUE)

# Datengruppen
y_1        = D$dBDI[D$Setting == "F2F"]                         # BDI-Differenzwerte in der F2F-Gruppe
y_2        = D$dBDI[D$Setting == "ONL"]                         # BDI-Differenzwerte in der ONL-Gruppe
y_3        = D$dBDI[D$Setting == "WLC"]                         # BDI-Differenzwerte in der ONL-Gruppe

# Modellformulierung
p          = 3                                                  # drei Gruppen
m          = length(y_1)                                        # balancierters Design mit n_i = 40
n          = p*m                                                # Datenvektordimension
y          = matrix(c(y_1, y_2, y_3), nrow = n)                 # Datenvektor
Xt         = cbind(matrix(1, nrow = n, ncol = 1),               # Designmatrix
                   kronecker(diag(p), matrix(1, nrow = m, ncol = 1)))
X          = Xt[,-2]                                            # eliminiere zweite Spalte

# Modellschätzung
beta_hat   = solve(t(X) %*% X) %*% t(X) %*% y                   # Betaparameterschätzer
eps_hat    = y - X %*% beta_hat                                 # Residuenvektor
sigsqr_hat = (t(eps_hat) %*% eps_hat) /(n-p)                    # Varianzparameterschätzer
s_sqr_123  = ((m-1)*var(y_1) +                                  # gebündelte Stichprobenvarianz
              (m-1)*var(y_2) +
              (m-1)*var(y_3))/(m+m+m-p)

# Ausgabe
cat(  "hat{beta}                                    : ", beta_hat,
    "\nbar{y}_1,bar{y}_2,bar{y}_3                   : ", c(mean(y_1),mean(y_2),mean(y_3)),
    "\nbar{y}_1,bar{y}_2-bar{y}_1,bar{y}_3-bar{y}_1 : ", c(mean(y_1),mean(y_2)-mean(y_1),mean(y_3)-mean(y_1)),
    "\nhat{sigsqr}                                  : ", sigsqr_hat,
    "\ns_123^2                                      : ", s_sqr_123,
    "\ns_y^2                                        : ", var(y))
```


<!-- Abschnitt 4: Modellevaluation -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

Modellformulierung

Modellschätzung

**Modellevaluation**

Anwendung/Praxis

Selbstkontrollfragen

\vfill


<!-- Überblick -->
# Modellevaluation
<!-- \Leftrightarrow \Theta_0 := \mathbb{R} \times \{0_{p-1}\} -->
<!-- \Leftrightarrow \Theta_1 := \mathbb{R}^p \setminus \mathbb{R} \times \{0_{p-1}\} -->

\small
\textcolor{darkblue}{Vorbemerkungen und Überblick}

\footnotesize
Prinzipiell sind alle Parameterschätzwerte in einem EVA-Modell von Interesse. Mit der T-Teststatistik können einzelne Komponenten oder lineare Kombinationen der Betaparameterschätzwerte im Sinne von Hypothesentests evaluiert werden.

Nichtsdestotrotz steht im EVA-Szenario häufig die Evaluation der Nullhypothese, dass alle Effektparameter gleich Null sind, im Vordergrund. Intuitiv besagt diese Nullhypothese, dass die wahren, aber unbekannten Erwartungswertparameter aller Gruppen identisch sind, formal haben wir
\begin{equation}
H_0 : \alpha_i = 0
\quad \mbox{für alle} \quad
i = 2,...,p
\end{equation}

und
\begin{equation}
H_1 : \alpha_i \neq 0
\quad \mbox{für mindestens ein} \quad
i = 2,...,p \; .
\end{equation}

Zur Überprüfung von $H_0$ wird im Allgemeinen eine $F$-Teststatistik eingesetzt. Dabei wird das volle Modell mit dem reduzierten Modell verglichen, in der die Designmatrix nur die erste Spalte der Designmatrix des vollen Modells enthält.

Im Folgenden entwickeln wir zunächst diese $F$-Teststatistik anhand einer *Quadratsummenzerlegung* der Datenvariabilität in einem EVA-Szenario und betrachten in diesem Zusammenhang auch das dem $R^2$ analoge Effektstärkemaß $\eta^2$. Ausgestattet mit der speziellen Form der $F$-Teststatistik in dem hier betrachteten Szenario diskutieren wir dann den traditionellen EVA-Hypothesentest. Wir verzichten auf eine Analyse der Testgütefunktion und eine Diskussion der Powerfunktion.


<!-- Theorem: Quadratsummenzerlegung bei einfaktorieller Varianzanalyse -->
# Modellevaluation

\vspace{2mm}
\footnotesize
\begin{theorem}[Quadratsummenzerlegung bei einfaktorieller Varianzanalyse]

\justifying
\normalfont
Für $i = 1,...,p$ und $j = 1,...,n_i$ sei $y_{ij}$ die $j$te Datenvariable in der $i$ten Gruppe eines EVA-Szenarios. Weiterhin seien mit $n := \sum_{i=1}^p n_i$
\center
\vspace{1mm}
\begin{tabular}{ll}
$\bar{y} := \frac{1}{n}\sum_{i=1}^p \sum_{j=1}^{n_i} y_{ij}$ & das \textit{Gesamtstichprobenmittel} \\\\
$\bar{y}_i =\frac{1}{n_i}\sum_{j=1}^{n_i} y_{ij}$            & das \textit{$i$te Stichprobenmittel} \\\\
\end{tabular}
\vspace{1mm}

\flushleft
sowie
\center
\vspace{1mm}
\begin{tabular}{ll}
$\mbox{SQT} := \sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij}-\bar{y})^{2}$   & die \textit{total sum of squares}           \\\\
$\mbox{SQB} := \sum_{i=1}^p n_i(\bar{y}_i-\bar{y})^{2}$              & die \textit{``between'' sum of squares}     \\\\
$\mbox{SQW} := \sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^{2}$ & die \textit{``within'' sum of squares} \; . \\\\
\end{tabular}
\vspace{1mm}
\flushleft

Dann gilt
\begin{equation}
\mbox{SQT} = \mbox{SQB} + \mbox{SQW}.
\end{equation}

\end{theorem}

\vspace{-1mm}
Bemerkungen
\vspace{-2mm}

* Im Vergleich mit der Quadratsummenzerlegung mit Ausgleichsgerade (vgl. Einheit (2) in *Allgemeines Lineares Modell*) spielt $\mbox{SQB}$ die Rolle der erklärten Quadratsumme und $\mbox{SQW}$ die Rolle der residuellen Quadratsumme.
* Die Ausdrücke "between" und "within" beziehen sich hierbei darauf, dass $\mbox{SQB}$ die Variabilität zwischen den Gruppen und $\mbox{SQW}$ die Variabilität innerhalb der Gruppen abbildet.


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis}

Es gilt
\footnotesize
\begin{align}
\begin{split}
\mbox{SQT}
& = \sum_{i=1}^p \sum_{j=1}^{n_i}(y_{ij}-\bar{y})^{2} \\
& = \sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i+\bar{y}_i-\bar{y})^2 \\
& = \sum_{i=1}^p \sum_{j=1}^{n_i} \left((y_{ij}-\bar{y}_i)+(\bar{y}_i-\bar{y}) \right)^2 \\
& = \sum_{i=1}^p \sum_{j=1}^{n_i} \left(
      (y_{ij}-\bar{y}_i)^2
      +2(y_{ij}-\bar{y}_i)(\bar{y}_i-\bar{y})
      +(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
       +\sum_{j=1}^{n_i}  2(y_{ij}-\bar{y}_i)(\bar{y}_i-\bar{y})
      +\sum_{j=1}^{n_i} (\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
      +2(\bar{y}_i-\bar{y})\sum_{j=1}^{n_i}(y_{ij}-\bar{y}_i)
      +n_i(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
      +2(\bar{y}_i-\bar{y})\sum_{j=1}^{n_i} \left(y_{ij}-\frac{1}{n_i}\sum_{j=1}^{n_i}{y_{ij}}\right)
      +n_i(\bar{y}_i-\bar{y})^2
    \right)
\end{split}
\end{align}


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis (fortgeführt)}

und weiter
\footnotesize
\begin{align*}
\begin{split}
\mbox{SQT}
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
      +2(\bar{y}_i-\bar{y})\left(\sum_{j=1}^{n_i} y_{ij}-\sum_{j=1}^{n_i}\left(\frac{1}{n_i}\sum_{j=1}^{n_i}{y_{ij}}\right)\right)
      +n_i(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
      +2(\bar{y}_i-\bar{y})\left(\sum_{j=1}^{n_i} y_{ij}-\frac{n_i}{n_i}\sum_{j=1}^{n_i}{y_{ij}}\right)
      +n_i(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2
      +2(\bar{y}_i-\bar{y})\left(\sum_{j=1}^{n_i} y_{ij}-\sum_{j=1}^{n_i}{y_{ij}}\right)
      +n_i(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \left(
      \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2 +n_i(\bar{y}_i-\bar{y})^2
    \right) \\
& = \sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij}-\bar{y}_i)^2 + \sum_{i=1}^p n_i(\bar{y}_i-\bar{y})^2
\end{split}
\end{align*}

\small
und damit
\begin{equation}
\mbox{SQT}  = \mbox{SQB} + \mbox{SQW}.
\end{equation}
$\hfill\Box$


<!-- Definition: Effektstärkemaß eta^2 -->
# Modellevaluation

\footnotesize
\begin{definition}[Effektstärkenmaß $\eta^2$]

\justifying
Für ein EVA-Szenario seien die \textit{between sum of squares} $\mbox{SQB}$ und die \textit{total sum of squares} $\mbox{SQT}$ definiert wie oben. Dann ist das \textit{Effektstärkenmaß} $\eta^2$ definiert als
\begin{equation}
\eta^2 := \frac{\mbox{SQB}}{\mbox{SQT}}
\end{equation}

\end{definition}

Bemerkungen

* $\eta^2$ ist analog zum Bestimmtheitsmaß $R^2$ der Regression definiert.
* $\eta^2$ gibt den Anteil der Varianz zwischen den Gruppen an der Gesamtvarianz der Daten an.
* Mit dem Theorem zur Quadratsummenzerlegung bei EVA folgt sofort $0 \le \eta^2 \le 1$, da
\begin{align}
\begin{split}
\mbox{SQB} &= 0 \Rightarrow \mbox{SQT} = \mbox{SQW} \quad \mbox{und} \quad \eta^2 = 0 \\
\mbox{SQW} &= 0 \Rightarrow \mbox{SQT} = \mbox{SQB} \quad \mbox{und} \quad \eta^2 = 1 \; .
\end{split}
\end{align}


<!-- Theorem: F-Teststatistik -->
# Modellevaluation

\footnotesize
\begin{theorem}[F-Teststatistik der einfaktoriellen Varianzanalyse]

\justifying
\normalfont
Es sei
\begin{equation}
y =  X\beta + \varepsilon
\quad \mbox{mit} \quad
\varepsilon \sim N(0_n,\sigma^2I_n)
\end{equation}

die Designmatrixform der Effektdarstellung mit Referenzgruppe des EVA-Modells und im Sinne der Definition der F-Statistik sei dieses Modell partioniert mit $p_0 := 1$ und $p_1 := p - 1$. Weiterhin seien
\center
\vspace{3mm}
\begin{tabular}{ll}
$\mbox{MSB} := \frac{\mbox{SQB}}{p-1}$ & die \textit{mean between sum of squares} und \\\\
$\mbox{MSW} := \frac{\mbox{SQW}}{n-p}$ & die \textit{mean within sum of squares}      \\\\
\end{tabular}

\flushleft
Dann gilt
\begin{equation}
F = \frac{\mbox{MSB}}{\mbox{MSW}}.
\end{equation}

\end{theorem}

Bemerkungen

* $p_0 := 1$ impliziert, dass das reduzierte Modell die Designmatrix $X_0 := 1_n$ hat.
* $p_0 := 1$ impliziert zudem, dass das reduzierte Modell den Betaparameter $\beta_0 := \mu_0$ hat.
* $p_0 := 1$ impliziert damit auch, dass das reduzierte Modell keine Effektparameter $\alpha_i$ für $i = 2,...,p$ hat.
* Die Zahl $p-1$ wird auch als "Zähler-Freiheitsgrade" bezeichnet.
* Die Zahl $n-p$ wird auch als "Nenner-Freiheitsgrade" bezeichnet.


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis}

Wir halten zunächst fest, dass für den Betaparameterschätzer des reduzierten Modells gilt, dass
\begin{align}
\begin{split}
\hat{\beta}_0
= (X_0^\mathrm{T}X_0)^{-1}X_0^\mathrm{T}y
= (1_n^\mathrm{T} 1_n)^{-1}1_n^\mathrm{T}y
= \frac{1}{n} \sum_{i=1}^p \sum_{j = 1}^{n_i} y_{ij}
= \bar{y}.
\end{split}
\end{align}

Weiterhin ergibt sich
\begin{align}
\begin{split}
\hat{\varepsilon}_0^\mathrm{T} \hat{\varepsilon}_0
= (y - X_0\hat{\beta}_0)^\mathrm{T}(y - X_0\hat{\beta}_0)
= (y - 1_n\bar{y})^\mathrm{T}(y - 1_n\bar{y})
= \sum_{i=1}^p \sum_{j=1}^{n_i} (y_{ij} - \bar{y})^2
= \mbox{SQT}.
\end{split}
\end{align}

Der Betaparameterschätzer des vollständigen Modells ergibt sich wie oben gesehen zu
\tiny
\begin{equation}
\renewcommand{\arraystretch}{1.4}
\hat{\beta} =
\begin{pmatrix}
\hat{\mu}_0    \\
\hat{\alpha}_2 \\
\vdots         \\
\hat{\alpha}_m
\end{pmatrix}
=
\begin{pmatrix}
\frac{1}{n_1}\sum_{j=1}^{n_1} y_{1j}                                        \\
\frac{1}{n_2}\sum_{j=1}^{n_2} y_{2j} - \frac{1}{n_1}\sum_{j=1}^{n_1} y_{1j} \\
\vdots                                                                      \\
\frac{1}{n_m}\sum_{j=1}^{n_m} y_{mj} - \frac{1}{n_1}\sum_{j=1}^{n_1} y_{1j} \\
\end{pmatrix}
=
\begin{pmatrix}
\bar{y}_1             \\
\bar{y}_2 - \bar{y}_1 \\
\vdots                \\
\bar{y}_p - \bar{y}_1 \\
\end{pmatrix},
\end{equation}

\small
sodass


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis (fortgeführt)}

\begin{align*}
\renewcommand{\arraystretch}{1}
\begin{split}
\hat{\varepsilon}^\mathrm{T}\hat{\varepsilon}
& = (y - X\hat{\beta})^\mathrm{T}(y - X\hat{\beta}) \\
& =
\left(
\begin{pmatrix}
y_{11}   \\
\vdots   \\
y_{1n_1} \\
y_{21}   \\
\vdots   \\
y_{2n_2} \\
         \\
\vdots   \\
         \\
y_{p1}   \\
\vdots   \\
y_{pn_p} \\
\end{pmatrix}
-
\begin{pmatrix}
1       &  0       &          &  0       \\
\vdots  &  \vdots  &  \cdots  &  \vdots  \\
1       &  0       &          &  0       \\
1       &  1       &          &  0       \\
\vdots  &  \vdots  &  \cdots  &  \vdots  \\
1       &  1       &          &  0       \\
        &          &          &          \\
\vdots  &  \vdots  &  \cdots  &  \vdots  \\
        &          &          &          \\
1       &  0       &          &  1       \\
\vdots  &  \vdots  &  \cdots  &  \vdots  \\
1       &  0       &          &  1       \\
\end{pmatrix}
\begin{pmatrix}
\bar{y}_1             \\
\bar{y}_2 - \bar{y}_1 \\
\vdots                \\
\bar{y}_p - \bar{y}_1 \\
\end{pmatrix}
\right)^\mathrm{T}
(y - X\hat{\beta})
\end{split}
\end{align*}


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis (fortgeführt)}

und weiter
\begin{align*}
\begin{split}
\hat{\varepsilon}^\mathrm{T}\hat{\varepsilon}
& =
\begin{pmatrix}
y_{11}   - \bar{y}_1                         \\
\vdots                                       \\
y_{1n_1} - \bar{y}_1                         \\
y_{21}   - \bar{y}_1 - \bar{y}_2 + \bar{y}_1 \\
\vdots                                       \\
y_{2n_2} - \bar{y}_1 - \bar{y}_2 + \bar{y}_1 \\
                                             \\
\vdots                                       \\
                                             \\
y_{p1}   - \bar{y}_1 - \bar{y}_p + \bar{y}_1 \\
\vdots                                       \\
y_{pn_p} - \bar{y}_1 - \bar{y}_p + \bar{y}_1 \\
\end{pmatrix}^\mathrm{T}
\begin{pmatrix}
y_{11}   - \bar{y}_1                         \\
\vdots                                       \\
y_{1n_1} - \bar{y}_1                         \\
y_{21}   - \bar{y}_1 - \bar{y}_2 + \bar{y}_1 \\
\vdots                                       \\
y_{2n_2} - \bar{y}_1 - \bar{y}_2 + \bar{y}_1 \\
                                             \\
\vdots                                       \\
                                             \\
y_{p1}   - \bar{y}_1 - \bar{y}_p + \bar{y}_1 \\
\vdots                                       \\
y_{pn_p} - \bar{y}_1 - \bar{y}_p + \bar{y}_1 \\
\end{pmatrix}
=
\begin{pmatrix}
y_{11}   - \bar{y}_1  \\
\vdots                \\
y_{1n_1} - \bar{y}_1  \\
y_{21}   - \bar{y}_2  \\
\vdots                \\
y_{2n_2} - \bar{y}_2  \\
                      \\
\vdots                \\
                      \\
y_{p1}   - \bar{y}_p  \\
\vdots                \\
y_{pn_p} - \bar{y}_p  \\
\end{pmatrix}^\mathrm{T}
\begin{pmatrix}
y_{11}   - \bar{y}_1 \\
\vdots               \\
y_{1n_1} - \bar{y}_1 \\
y_{21}   - \bar{y}_2 \\
\vdots               \\
y_{2n_2} - \bar{y}_2 \\
                     \\
\vdots               \\
                     \\
y_{p1}   - \bar{y}_p \\
\vdots               \\
y_{pn_p} - \bar{y}_p \\
\end{pmatrix} \\
& = \sum_{i=1}^p \sum_{j = 1}^{n_i} (y_{ij} - \bar{y}_i)^2 \\
& = \mbox{SQW}.
\end{split}
\end{align*}


<!-- Beweis -->
# Modellevaluation

\small
\underline{Beweis (fortgeführt)}

Mit dem Theorem zur Quadratsummenzerlegung bei einfaktorieller Varianzanalyse
\begin{equation}
\mbox{SQT} = \mbox{SQB} + \mbox{SQW}
\quad \Leftrightarrow \quad
\mbox{SQB} = \mbox{SQT} - \mbox{SQW}
\end{equation}

folgt sofort, dass
\begin{align}
\begin{split}
   \mbox{SQB}
&= \mbox{SQT} - \mbox{SQW} \\
&= \hat{\varepsilon}_0^\mathrm{T}\hat{\varepsilon}_0 - \hat{\varepsilon}^\mathrm{T}\hat{\varepsilon}.
\end{split}
\end{align}
Dann aber folgt auch direkt, dass
\begin{equation}
  \frac{\mbox{MSB}}{\mbox{MSW}}
= \frac{\frac{\mbox{SQB}}{p-1}}{\frac{\mbox{SQW}}{n-p}}
= \frac{\frac{\hat{\varepsilon}_0^\mathrm{T}\hat{\varepsilon}_0 - \hat{\varepsilon}^\mathrm{T}\hat{\varepsilon}}{p-1}}{\frac{\hat{\varepsilon}^\mathrm{T}\hat{\varepsilon}}{n-p}}
= F
\end{equation}


<!-- Theorem: Effektstärkenmaß eta^2 und F-Teststatistik -->
# Modellevaluation

\footnotesize
\begin{theorem}[Effektstärkenmaß $\eta^2$ und F-Teststatistik]

\justifying
\normalfont
Für ein EVA-Szenario mit $p$ Gruppen und Gesamtdatenpunktanzahl $n$ seien das Effektstärkenmaß $\eta^2$ und die $F$-Teststatistik wie oben definiert. Dann gilt
\begin{equation}
\eta^2 = \frac{F(p-1)}{F(p-1) + (n-p)}
\end{equation}
\end{theorem}

Bemerkungen

* Das Verhältnis von $F$ und $\eta^2$ ist analog zum Verhältnis von $T$ und Cohen's $d$.
* Die gleichzeitige Angabe von $F$ und $\eta^2$ ist redundant.


<!-- Beweis -->
# Modellevaluation

\small
\setstretch{1}
\underline{Beweis}

Wir halten zunächst fest, dass
\begin{equation}
F = \frac{\mbox{SQB}}{\mbox{SQW}}\cdot \frac{n-p}{p-1} 
\quad \Leftrightarrow \quad
\mbox{SQB} = \frac{p-1}{n-p}\cdot \mbox{SQW} \cdot F \; .
\end{equation}

Damit folgt dann
\begin{align}
\begin{split}
\eta^2
& = \frac{\mbox{SQB}}{\mbox{SQT}}  
  = \frac{\mbox{SQB}}{\mbox{SQB}+\mbox{SQW}}
  = \frac{\frac{p-1}{n-p}\cdot \mbox{SQW} \cdot F}{\frac{p-1}{n-p}\cdot \mbox{SQW} \cdot F +\mbox{SQW}} \\
& = \frac{\frac{F(p-1)}{n-p}\cdot \mbox{SQW} }{\frac{F(p-1)}{n-p}\cdot \mbox{SQW}  +\mbox{SQW}}
  = \frac{\frac{F(p-1)}{n-p}\cdot \mbox{SQW} }{\left(\frac{F(p-1)}{n-p} + 1 \right)\cdot \mbox{SQW}} \\
& = \frac{\frac{F(p-1)}{n-p}}{\frac{F(p-1)}{n-p} + \frac{n-p}{n-p}}
  = \frac{\frac{F(p-1)}{n-p}}{\frac{F(p-1) + (n-p)}{n-p}}
  = \frac{F(p-1)}{F(p-1)+(n-p)} \; .
\end{split}
\end{align}
$\hfill\Box$


<!-- Modellevaluation -->
# Modellevaluation

\small
\setstretch{2}
\textcolor{darkblue}{Gliederung} (vgl. [\textcolor{darkblue}{Einheit (12), um 01:58:00}](https://youtu.be/GJsS-AZe9zI?si=9kldNuX_I3bgjUFE&t=7080) in \textit{Wahrscheinlichkeitstheorie und Frequentistische Inferenz})

(1) Statistisches Modell $\checkmark$

(2) Testhypothesen $\checkmark$

(3) Teststatistik $\checkmark$

(4) Test

(5) Testumfangkontrolle

(6) p-Wert


<!-- Modellevaluation: Teststatistik -->
# Modellevaluation (3) Teststatistik

\setstretch{3}
\textcolor{darkblue}{Bezeichnungskonventionen für $f$-Zufallsvariablen}

\small
* Für eine (nichtzentrale) $f$-Zufallsvariable $\xi$ schreiben wir $\xi \sim f(n_1, n_2)$ (oder $\xi \sim f(\delta, n_1, n_2)$).
* Die WDF einer (nichtzentralen) $f$-Zufallsvariable ist $f(\cdot; n_1, n_2)$ (oder $f(\cdot;\delta, n_1, n_2)$).
* Die KVF einer (nichtzentralen) $f$-Zufallsvariable ist $\varphi(\cdot; n_1, n_2)$ (oder $\varphi(\cdot; \delta, n_1, n_2)$).
* Die Teststatistik eines F-Test-Designs/Hypothesentests bezeichnen wir mit $F$.
* Die Realisierung der F-Teststatistik für einen Datensatz bezeichnen wir mit $f$.
\vfill


<!-- Definition: F-Test für einfaktorielle Varianzanalyse -->
# Modellevaluation (4) Test
<!-- \Leftrightarrow \Theta_0 := \mathbb{R} \times \{0_{p-1}\} -->
<!-- \Leftrightarrow \Theta_1 := \mathbb{R}^p \setminus \mathbb{R} \times \{0_{p-1}\} -->

\footnotesize
\begin{definition}[F-Test für einfaktorielle Varianzanalyse]

\justifying
Gegeben sie das Modell der einfaktoriellen Varianzanalyse sowie die zusammengesetzten Null- und Alternativhypothesen
\begin{equation}
H_0 : \alpha_i = 0
\quad \mbox{für alle} \quad
i = 2,...,p
\end{equation}

bzw.
\begin{equation}
H_1 : \alpha_i \neq 0
\quad \mbox{für mindestens ein} \quad
i = 2,...,p \; .
\end{equation}

Weiterhin sei die F-Teststatistik definiert durch
\begin{equation}
F = \frac{\mbox{MSB}}{\mbox{MSW}}
\end{equation}

mit der mean between sum of squares $\mbox{MSB}$ und der mean within sum of squares $\mbox{MSW}$. Dann ist der \textit{F-Test für die einfaktoriellen Varianzanalyse (EVA-F-Test)} definiert als der kritische Wert-basierte Test
\begin{equation}
\phi(y) := 1_{\{F \ge k\}} :=
\begin{cases}
1 & F \ge k \\
0 & F <   k
\end{cases}.
\end{equation}

\end{definition}


<!-- Theorem: Testumfangkontrolle -->
# Modellevaluation (5) Testumfangkontrolle

\footnotesize
\begin{theorem}[Testumfangkontrolle]

\justifying
\normalfont
$\phi$ sei der F-Test zur einfaktoriellen Varianzanalyse. Dann ist $\phi$ ein Level-$\alpha_0$-Test mit Testumfang $\alpha_0$, wenn der kritische Wert definiert ist durch
\begin{equation}
k_{\alpha_0} := \varphi^{-1}(1-\alpha_0; p-1, n-p),
\end{equation}
wobei $\varphi^{-1}(\cdot; p-1, n-p)$ die inverse KVF der $f$-Verteilung mit Freiheitsgradparametern $p-1$ und $n-p$ ist.

\end{theorem}
\vfill


<!-- Visualisierung -->
# Modellevaluation (5) Testumfangkontrolle

\small
\center
Wahl von $k_{\alpha_0} := \varphi^{-1}(1 - \alpha_0; p-1, n-p)$ mit $p = 3, m = 12, \alpha_0 := 0.05$ und Ablehnungsbereich

```{r, echo = F, eval = F}
par(
    family      = "sans",
    mfcol       = c(1,2),
    pty         = "m",
    bty         = "l",
    lwd         = 1,
    las         = 1,
    mgp         = c(2,1,0),
    xaxs        = "i",
    yaxs        = "i",
    font.main   = 1,
    cex         = 1,
    cex.main    = 1.2)

# Parameter
p           = 3                                                                  # Anzahl Gruppen
m           = 12                                                                 # Anzahl experimenteller Einheiten pro Gruppe
n           = p*m                                                                # Gesamtstichprobengröße
alpha_0     = 0.05                                                               # Signifikanzlevel
k_alpha_0   = qf(1 - alpha_0, p-1,n-p)                                           # kritischer Wert
eff         = seq(0,6,length=1e4 )                                               # F-Statistikwerte
Peff        = pf(eff,p-1,n-p)                                                    # F-Statistik KVF für H_0
peff        = df(eff,p-1,n-p)                                                    # F-Statistik WDF für H_0

# KVF-Perspektive
plot(eff, Peff,
    type        = "l",
    xlab        = " ",
    ylab        = " ",
    ylim        = c(0,1),
    main        = TeX("$\\varphi(F;2,33)$"))
lines(k_alpha_0, 0,
    type        = "p",
    pch         = 16,
    xpd         = TRUE)
lines(min(eff), 1 - alpha_0,
    type        = "p",
    pch         = 16,
    xpd         = TRUE)
arrows(
    x0          = min(eff),
    y0          = 1 - alpha_0,
    x1          = k_alpha_0,
    y1          = 1 - alpha_0,
    col         = "darkorange",
    angle       = 45,
    length      = .1)
arrows(
    x0          = k_alpha_0,
    y0          = 1-alpha_0,
    x1          = k_alpha_0,
    y1          = 0,
    col         = "darkorange",
    angle       = 45,
    length      = .1)
text(k_alpha_0-0.5, 0.1, TeX("$\\k_{\\alpha_0}$"), xpd = TRUE)
text(1,             1,   TeX("$1 - \\alpha_0$"),   xpd = TRUE)

# WDF-Perspektive
plot(eff, peff,
    type        = "l",
    ylab        = " ",
    xlab        = " ",
    ylim        = c(0,.4),
    main        = TeX("$f(F;2,33)$"))
polygon(c(eff[eff  >= k_alpha_0], max(eff), k_alpha_0),
        c(peff[eff >= k_alpha_0], 0,        0        ),
    col         = "gray90",
    border      = NA)
lines(seq(k_alpha_0, max(eff), len = 1e2), rep(0,1e2),
    type        = "l",
    lwd         = 5,
    col         = "darkorange")
lines(k_alpha_0, 0,
    type        = "p",
    pch         = 16,
    xpd         = TRUE)
text(k_alpha_0-0.5, 0.03, TeX("$\\k_{\\alpha_0}$") ,                  xpd = TRUE)
text(5,             0.05, TeX("$P(F >= k_{\\alpha_0}) = \\alpha_0$"), xpd = TRUE, cex = 1, col = "gray50")

# PDF-Speicherung
dev.copy2pdf(
    file        = "Abbildungen/eva_testumfang.pdf",
    width       = 8,
    height      = 4)
```

```{r, echo = F, out.width = "100%"}
knitr::include_graphics("Abbildungen/eva_testumfang.pdf")
```
\vfill


<!-- Beweis -->
# Modellevaluation (5) Testumfangkontrolle

\footnotesize
\underline{Beweis}

Die Testgütefunktion des betrachteten Tests im vorliegenden Testszenario ist definiert als
\begin{equation}
q : \mathbb{R}^p \to [0,1], \beta \mapsto q_\phi(\beta) := \mathbb{P}_\beta(\phi = 1).
\end{equation}

Wir haben bereits gesehen (siehe Einheit (8) in *Allgemeines Lineares Modell*), dass die F-Teststatistik für $p_1 = p - 1$ nach einer nichtzentralen $f$-Verteilung verteilt ist, 
\begin{equation}
F \sim f(\delta,p-1, n-p).
\end{equation}

Weiterhin ist der Ablehnungsbereich des hier betrachteten Tests gegeben als $[k,\infty[$. Für die funktionale Form der Testgütefunktion ergibt sich also
\begin{align}
\begin{split}
   \mathbb{P}_\beta(\phi = 1)
&= \mathbb{P}_\beta(F \in [k,\infty[) \\
&= \mathbb{P}_\beta(F \ge k) \\
&= 1 - \mathbb{P}_\beta(F \le k) \\
&= 1 - \varphi(k;\delta,p-1,n-p), \\
\end{split}
\end{align}
wobei $\varphi(k;\delta,p-1,n-p)$ den KVF-Wert der nichtzentralen $f$-Verteilung an Stelle $k$ mit Nichtzentralitätsparameter $\delta$ sowie Freiheitsgradparametern $p-1$ und $n-p$ bezeichnet (vgl. Einheit (8) in *Allgemeines Lineares Modell*).

Damit der betrachtete Test ein Level-$\alpha_0$-Test ist, muss bekanntlich gelten, dass
\begin{equation}
q_\phi(\beta) \le \alpha_0
\quad \mbox{für alle} \quad
\beta \in \Theta_0
\quad \mbox{mit} \quad
\Theta_0 = \mathbb{R} \times \{0_{p-1}\}.
\end{equation}


<!-- Beweis -->
# Modellevaluation (5) Testumfangkontrolle

\vspace{2mm}
\footnotesize
\underline{Beweis}

Der Nichtzentralitätsparameter ist gegeben durch (siehe Einheit (8) in *Allgemeines Lineares Modell*)
\begin{equation}
\delta = \frac{1}{\sigma^2} c^\mathrm{T}\beta \left(c^\mathrm{T}(X^\mathrm{T}X)^{-1}c\right)^{-1}c^\mathrm{T} \beta \; .
\end{equation}

Mit dem dazugehörigen $c$ und unter der Nullhypothese, dass $\beta \in \Theta_0$, also
\begin{equation}
c = \begin{pmatrix} 0 \\ 1_{p-1} \end{pmatrix} \in \mathbb{R}^p
\quad \mbox{und} \quad
\beta = \begin{pmatrix} \mu_0 \\ 0_{p-1} \end{pmatrix} \in \mathbb{R}^p \; ,
\end{equation}

folgt dann aber $\delta = 0$ und somit
\begin{equation}
q_\phi(\beta) = 1 - \varphi(k; p-1, n-p)
\quad \mbox{für alle} \quad
\beta \in \Theta_0 \; ,
\end{equation}

wobei $\varphi(k;p-1,n-p)$ den Wert der KVF der $f$-Verteilung an der Stelle $k$ mit Freiheitsgradparametern $p-1$ und $n-p$ bezeichnet (vgl. Einheit (8) in *Allgemeines Lineares Modell*). Der Testumfang des betrachteten Tests ergibt sich per Definition (siehe Einheit (12) in *Wahrscheinlichkeitstheorie und Frequentistische Inferenz*) als
\begin{equation}
\alpha = \max_{{\beta \in \Theta_0}} q_{\phi}(\beta) = 1 - \varphi(k;p-1,n-p),
\end{equation}
da $q_{\phi}(\beta)$ für $\beta \in \Theta_0$ nicht von $\mu_0$ abhängt. Wir müssen also lediglich zeigen, dass die Wahl von $k_{\alpha_0}$ wie im Theorem garantiert, dass $\phi$ den Testumfang $\alpha_0$ hat. Sei also $k := k_{\alpha_0}$. Dann gilt für alle $\beta \in \Theta_0$
\begin{equation}
q_\phi(\beta) 
= 1-\varphi( \varphi^{-1}(1-\alpha_0; p-1, n-p); p-1, n-p)
= 1-(1-\alpha_0)
= \alpha_0
\end{equation}

und damit ist alles gezeigt.
$\hfill\Box$


<!-- Modellevaluation: p-Wert -->
# Modellevaluation (6) p-Wert

\footnotesize
Nach Definition ist der p-Wert das kleinste Signifikanzlevel $\alpha_0$, bei welchem man die Nullhypothese basierend auf einem vorliegenden Wert der Teststatistik ablehnen würde. Wir wollen einen vorliegenden Wert der $F$-Teststatistik hier mit $f$ bezeichnen.

Bei $F = f$ würde $H_0$ für jedes $\alpha_0$ mit $f \ge \psi^{-1}(1-\alpha_0;p-1,n-p)$ abgelehnt werden. Für ein solches $\alpha_0$ gilt aber
\begin{equation}
\alpha_0 \ge \mathbb{P}(F \ge f),
\end{equation}

denn
\begin{align}
\begin{split}
f & \ge \psi^{-1}(1-\alpha_0;p-1,n-p) \\
\Leftrightarrow \psi(f;p-1,n-p) & \ge \psi(\psi^{-1}(1-\alpha_0;p-1,n-p), p-1,n-p) \\
\Leftrightarrow \psi(f;p-1,n-p) & \ge 1-\alpha_0 \\
\Leftrightarrow \mathbb{P}(F \le f) & \ge 1-\alpha_0 \\
\Leftrightarrow \alpha_0  & \ge 1 - \mathbb{P}(F \le f) \\
\Leftrightarrow \alpha_0  & \ge \mathbb{P}(F \ge f).
\end{split}
\end{align}

Das kleinste $\alpha_0 \in [0,1]$ mit $\alpha_0 \ge \mathbb{P}(F \ge f)$ ist 
dann $\alpha_0 = \mathbb{P}(F \ge f)$, also folgt
\begin{equation}
\mbox{p-Wert} = \mathbb{P}(F \ge f) = 1 - \varphi(f; p-1, n-p).
\end{equation}


<!-- Praktisches Vorgehen -->
# Modellevaluation (6) p-Wert

Praktisches Vorgehen
\footnotesize

* \justifying Man nimmt an, dass ein vorliegender Datensatz aus $p$ Gruppendatensätzen besteht, wobei $y_{11}, ...,y_{1n_1}$, $y_{21}, ...,y_{2n_2}$, ..., $y_{p1}, ...,y_{pn_p}$ Realisationen von $y_{1j} \sim N(\mu_0,\sigma^2)$ und $y_{ij} \sim N(\mu_0 + \alpha_i,\sigma^2)$ für $i = 2,...,p$ mit wahren, aber unbekannten Parametern $\mu_0, \alpha_i, i = 2,...,p$ und $\sigma^2>0$ sind.

* Man möchte entscheiden, ob $H_0 : \alpha_i = 0$ für alle $i = 2,...,p$ eher zutrifft oder eher nicht.

* Man wählt ein Signifikanzniveau $\alpha_0$ und bestimmt den zugehörigen Freiheitsgradparameter-abhängigen kritischen Wert $k_{\alpha_0}$. Zum Beispiel gilt bei Wahl von $\alpha_0 := 0.05, p = 3, m = 12, i = 1,2,3$ und somit $n = 36$, dass $k_{\alpha_0} = \varphi^{-1}(1 - 0.05; 3-1, 36-3) \approx 3.28$ ist.

* Anhand der MSB und MSW berechnet man den realisierten Wert der F-Teststatistik, den wir hier mit $f$ bezeichnen.

* Wenn $f$ größer-gleich $k_{\alpha_0}$ ist, lehnt man die Nullhypothese ab, andernfalls nicht.

* Die oben entwickelte Theorie garantiert dann, dass man im Mittel in höchstens $\alpha_0 \cdot 100$ von 100 Fällen die Nullhypothese fälschlicherweise ablehnt.

* Schließlich ergibt sich der assoziierte p-Wert der realisiertern F-Teststatistik $f$ zu
\begin{equation}
\mbox{p-Wert} = \mathbb{P}(F \ge f) = 1 - \varphi(f; p-1, n-p)
\end{equation}


<!-- Abschnitt 5: Anwendung/Praxis -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

Modellformulierung

Modellschätzung

Modellevaluation

**Anwendung/Praxis**

Selbstkontrollfragen

\vfill


<!-- Anwendungsszenario -->
# Anwendung/Praxis

\small
\vspace{2mm}
\textcolor{darkblue}{Daten einlesen:} $j = 1,...,10$ für jede Gruppe

\tiny
\setstretch{0.6}

\vspace{1mm}
```{r, echo = T}
fname       = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D           = read.table(fname, sep = ",", header = TRUE)
```

\vspace{-1mm}
```{r, echo = F}
knitr::kable(D[c(1:10, 41:50, 81:90),2:6],
             align = "cccrr",
             "pipe")
```


<!-- Visualisierung: Boxplot -->
# Anwendung/Praxis

\textcolor{darkblue}{Boxplot}
\vspace{2mm}

\tiny
\setstretch{1}
```{r, echo = T, eval = F}
# Daten einlesen
fname       = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D           = read.table(fname, sep = ",", header = TRUE)

# Abbildungsparameter
par(                                        # für Details siehe ?par
    family      = "sans",                   # Serif-freier Fonttyp
    pty         = "m",                      # maximale Abbildungsregion
    bty         = "l",                      # L-förmige Box
    lwd         = 1,                        # Liniendicke
    las         = 1,                        # horizontale Achsenbeschriftung
    font.main   = 1,                        # Titel nicht fett
    cex         = 1,                        # Textvergrößerungsfaktor
    cex.main    = 1.2)                      # Titeltextvergrößerungsfaktor

# Boxplot
boxplot(dBDI ~ Setting,                     # BDI ~ Setting enkodiert die Datengruppierung
    data        = D)

# PDF-Speicherung
dev.copy2pdf(
    file        = "Abbildungen/eva_boxplot.pdf",
    width       = 7,
    height      = 4)
```


<!-- Visualisierung: Boxplot -->
# Anwendung/Praxis

\textcolor{darkblue}{Boxplot}

\center
\vspace{-2mm}

```{r, echo = F, out.width = "100%"}
knitr::include_graphics("Abbildungen/eva_boxplot.pdf")
```


<!-- Visualisierung: Balkendiagramm -->
# Anwendung/Praxis

\textcolor{darkblue}{Balkendiagramm}
\vspace{2mm}

\tiny
\setstretch{0.8}
```{r, echo = T, eval = F}
# Daten einlesen
fname       = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D           = read.table(fname, sep = ",", header = TRUE)

# Datenselektion
A = data.frame(F2F = D$dBDI[D$Setting == "F2F"],   # BDI-Daten für F2F
               ONL = D$dBDI[D$Setting == "ONL"],   # BDI-Daten für ONL
               WLC = D$dBDI[D$Setting == "WLC"])   # BDI-Daten für WLC

# deskriptive Statistiken
groupmeans  = colMeans(A)                          # Gruppenmittelwerte
groupstds   = apply(A,2,sd)                        # Gruppenstandardabweichungen

# Balkendiagramm
par(                                               # für Details siehe ?par
    family      = "sans",                          # Serif-freier Fonttyp
    pty         = "m",                             # maximale Abbildungsregion
    bty         = "l",                             # L-förmige Box
    lwd         = 1,                               # Liniendicke
    las         = 1)                               # horizontale Achsenbeschriftung
x = barplot(groupmeans,                            # Ausgabe der x-Ordinaten
    ylim        = c(-5,15),
    col         = "gray90",
    ylab        = "dBDI",
    xlab        = "Condition")
arrows(
    x0          = x,                               # arrow start x-ordinate
    y0          = groupmeans - groupstds,          # arrow start y-ordinate
    x1          = x,                               # arrow end   x-ordinate
    y1          = groupmeans + groupstds,          # arrow end   y-ordinate
    code        = 3,                               # Pfeilspitzen beiderseits
    angle       = 90,                              # Pfeilspitzenwinkel: Linie
    length      = 0.05)                            # Linielänge

# PDF-Speicherung
dev.copy2pdf(
    file        = "Abbildungen/eva_barplot.pdf",
    width       = 7,
    height      = 4)
```


<!-- Visualisierung: Balkendiagramm -->
# Anwendung/Praxis

\textcolor{darkblue}{Balkendiagramm}
\vspace{2mm}

\center
Gruppenmittelwerte $\pm$ Gruppenstandardabweichungen
\vspace{-2mm}

```{r, echo = F, out.width = "100%"}
knitr::include_graphics("Abbildungen/eva_barplot.pdf")
```


<!-- Anwendungsbeispiel -->
# Anwendung/Praxis

\textcolor{darkblue}{Anwendungsbeispiel}

\tiny
\setstretch{0.9}
\vspace{2mm}
```{r, echo = T}
# Daten einlesen
fname      = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D          = read.table(fname, sep = ",", header = TRUE)

# Datengruppen
y_1        = D$dBDI[D$Setting == "F2F"]                # BDI-Differenzwerte in der F2F-Gruppe
y_2        = D$dBDI[D$Setting == "ONL"]                # BDI-Differenzwerte in der ONL-Gruppe
y_3        = D$dBDI[D$Setting == "WLC"]                # BDI-Differenzwerte in der WLC-Gruppe

# Modellformulierung
p          = 3                                         # p = 3; drei Gruppen
m          = length(y_1)                               # m = 40; balancierters Design
n          = p*m                                       # n = 120; Datenvektordimension
y          = matrix(c(y_1, y_2, y_3), nrow = n)        # Datenvektor
Xt         = cbind(matrix(1, nrow = n, ncol = 1),      # Designmatrix vollständiges Modell
                   kronecker(diag(p), matrix(1, nrow = m, ncol = 1)))
X          = Xt[,-2]                                   # eliminiere zweite Spalte
X_0        = X[,1]                                     # Designmatrix reduziertes Modell

# Evaluation der F-Teststatistik
beta_hat   = solve(t(X)   %*% X)   %*% t(X)   %*% y    # Betaparameterschätzer vollständiges Modell
beta_hat_0 = solve(t(X_0) %*% X_0) %*% t(X_0) %*% y    # Betaparameterschätzer reduziertes Modell
eps_hat    = y - X   %*% beta_hat                      # Residuenvektor vollständiges Modell
eps_hat_0  = y - X_0 %*% beta_hat_0                    # Residuenvektor reduziertes Modell
SQT        = t(eps_hat_0) %*% eps_hat_0                # total sum of squares
SQW        = t(eps_hat)   %*% eps_hat                  # within sum of squares
SQB        = SQT - SQW                                 # between sum of squares
DFB        = p - 1                                     # between degrees of freedom
DFW        = n - p                                     # within degrees of freedom
MSB        = SQB/DFB                                   # mean between sum of squares
MSW        = SQW/DFW                                   # mean within sum of squares
Eff        = MSB/MSW                                   # F-Teststatistik
pval       = 1 - pf(Eff, p-1, n-p)                     # p-Wert
```


<!-- Anwendungsbeispiel -->
# Anwendung/Praxis

\textcolor{darkblue}{Anwendungsbeispiel}

\small
\vspace{2mm}
```{r, echo = T}
# Ausgabe
cat(  "DFB :", DFB,
    "\nDFW :", DFW,
    "\nSQB :", SQB,
    "\nSQW :", SQW,
    "\nMSB :", MSB,
    "\nMSW :", MSW,
    "\nF   :", Eff,
    "\np   :", paste(pval))
```


<!-- Anwendungsbeispiel -->
# Anwendung/Praxis

\textcolor{darkblue}{Anwendungsbeispiel}

\small
\vspace{2mm}
```{r, echo = T}
# Daten einlesen
fname   = "Daten/Einfaktorielle_Varianzanalyse_Daten.csv"
D       = read.table(fname, sep = ",", header = TRUE)

# Benutzung von R's aov Funktion und Ausgabe
res.aov = aov(dBDI ~ Setting, data = D)
summary(res.aov)
```
\vfill


<!-- Datensatz: Great 10k -->
# Anwendung/Praxis

\small
Datensatz: Zielzeiten beim Great 10k (08.10.2017), getrennt nach Altersklasse ($n = 4883$)

```{r, echo = F, eval = T}
# Daten einlesen
filename = "Daten/great_10k.csv"
D        = read.csv(filename)

# Daten reduzieren
excl     = c("FU14", "FU16", "FU18", "FU20", "FU23", "F70", "F75", "F80", "MU14", "MU16", "MU18", "MU20", "MU23", "M70", "M75", "M80")
incl     = c("FE", "FH", "F30", "F35", "F40", "F45", "F50", "F55", "F60", "F65", "ME", "MH", "M30", "M35", "M40", "M45", "M50", "M55", "M60", "M65")
D        = D[!(D$AgeGroup %in% excl),]
D        = D[!(D$X5kTime  == "NaN:NaN:NaN"),]
D        = D[!(D$X10kTime == "NaN:NaN:NaN"),]
n        = nrow(D)

# Daten umsortieren
D$AgeGroup = factor(D$AgeGroup, levels = incl)

# Daten umrechnen
t5_hms   = D$X5kTime
t5_min   = rep(0, n)
t10_hms  = D$X10kTime
t10_min  = rep(0, n)
for(i in 1:n){
    t5         = t5_hms[i]
    t10        = t10_hms[i]
    t5_min[i]  = as.numeric(substr(t5,1,2))*60  + as.numeric(substr(t5,4,5))  + as.numeric(substr(t5,7,8))/60
    t10_min[i] = as.numeric(substr(t10,1,2))*60 + as.numeric(substr(t10,4,5)) + as.numeric(substr(t10,7,8))/60
}
D$T5k1   = t5_min
D$T5k2   = t10_min-t5_min
D$T10k   = t10_min
D$Tdiff  = D$T5k2 - D$T5k1
```

```{r, echo = F, eval = F}
# Daten visualisieren
library(vioplot)
par(
    family      = "sans",
    mar         = c(3,3,1,1),
    pty         = "m",
    bty         = "l",
    lwd         = 1,
    las         = 1,
    mgp         = c(2,1,0),
    xaxs        = "i",
    yaxs        = "i",
    font.main   = 1.5,
    cex         = 1.5)

# Violinplot
vioplot(D$T10k ~ D$AgeGroup, D,
    col         = "gray80",
    rectCol     = "black",
    lineCol     = "white",
    colMed      = "black",
    border      = "black",
    pchMed      = 16,
    plotCentre  = "lines",
    xlab        = "Altersklasse",
    ylab        = "Zielzeit [min]",
    ylim        = c(25,85),
    drawRect    = FALSE)

# Datenpunkte
stripchart(D$T10k ~ D$AgeGroup, D,
    method      = "jitter",
    xaxt        = "n",
    vertical    = TRUE,
    pch         = 19,
    col         = "black",
    add         = TRUE,
    cex         = 0.4)

# Mittelwerte
dx = 0.3
for (i in 1:length(incl)) {
    segments(
        x0      = i - dx,
        x1      = i + dx,
        y0      = median(D$T10k[D$AgeGroup==incl[i]]),
        y1      = median(D$T10k[D$AgeGroup==incl[i]]),
        col     = "black",
        lwd     = 2)
}

# Speichern
dev.copy2pdf(
    file        = "Abbildungen/great_10k.pdf",
    width       = 16,
    height      = 9)
```

\vspace{2mm}
```{r, echo = F, out.width = "100%"}
knitr::include_graphics("Abbildungen/great_10k.pdf")
```


<!-- Einfaktorielle Varianzanalyse (alle) -->
# Anwendung/Praxis

\small
Einfaktorielle Varianzanalyse: 10-km-Zeit, alle Altersgruppen ($n = 4883$)

\vspace{2mm}
\footnotesize
```{r, echo = T, eval = T}
# Daten extrahieren
D$T10k     = t10_min                                    # 10-km-Zeit [min]

# Einfaktorielle Varianzanalyse (falsch)
grps       = c("FE", "FH", "F30", "F35", "F40", "F45", "F50", "F55", "F60", "F65",
               "ME", "MH", "M30", "M35", "M40", "M45", "M50", "M55", "M60", "M65")
p          = length(grps)                               # Anzahl Gruppen
n_i        = rep(0,p)                                   # Anzahl Datenpunkte pro Gruppe
y          = matrix(rep(1,0))                           # Initialisierung Datenvektor
for (i in 1:p) {
    y_i    = matrix(D$T10k[D$AgeGroup == grps[i]], ncol = 1)
    n_i[i] = length(y_i)                                # Anzahl Datenpunkte Gruppe i
    y      = rbind(y, y_i)                              # Ergänzung Datenvektor
}
n          = sum(n_i)                                   # Anzahl Datenpunkte
```

\scriptsize
```{r, echo = F, eval = T}
# Einfaktorielle Varianzanalyse
X          = matrix(rep(1,n), nrow = n)                 # Initialisierung Designmatrix
for (i in 2:p) {
    n_i1   = sum(n_i[1:(i-1)])
    n_i2   = n_i[i]                                     # Anzahl Datenpunkte Gruppe i
    n_i3   = sum(n_i[(i+1):p])
    if (i == p) { n_i3 = 0 }                            # Ergänzung Designmatrix
    X      = cbind(X, matrix(c(rep(0,n_i1), rep(1,n_i2), rep(0,n_i3)), nrow = n))
}
X_0        = X[,1]                                      # Designmatrix reduziertes Modell

# Evaluation der F-Teststatistik
beta_hat   = solve(t(X)   %*% X)   %*% t(X)   %*% y     # Betaparameterschätzer vollständiges Modell
beta_hat_0 = solve(t(X_0) %*% X_0) %*% t(X_0) %*% y     # Betaparameterschätzer reduziertes Modell
eps_hat    = y - X   %*% beta_hat                       # Residuenvektor vollständiges Modell
eps_hat_0  = y - X_0 %*% beta_hat_0                     # Residuenvektor reduziertes Modell
SQT        = t(eps_hat_0) %*% eps_hat_0                 # total sum of squares
SQW        = t(eps_hat)   %*% eps_hat                   # within sum of squares
SQB        = SQT - SQW                                  # between sum of squares
DFB        = p - 1                                      # between degrees of freedom
DFW        = n - p                                      # within degrees of freedom
MSB        = SQB/DFB                                    # mean between sum of squares
MSW        = SQW/DFW                                    # mean within sum of squares
Eff        = MSB/MSW                                    # F-Teststatistik
pval       = 1 - pf(Eff, p-1, n-p)                      # p-Wert

# Ausgabe der Ergebnisse
cat(  "Datenpunkte pro Gruppe            : ", n_i,
    "\nBetaparameterschätzer             : ", round(beta_hat,   digits = 1),
    "\nBetaparameterschätzer red. Mod.   : ", round(beta_hat_0, digits = 3),
    "\nresiduelle Quadratsumme           : ", round(SQW,        digits = 3),
    "\nresiduelle Quadratsumme red. Mod. : ", round(SQT,        digits = 3),
    "\nF-Teststatistik                   : ", round(Eff,        digits = 3),
    "\np-Wert                            : ", round(pval,       digits = 3), "\n\n")
```


<!-- Einfaktorielle Varianzanalyse (Frauen) -->
# Anwendung/Praxis

\small
Einfaktorielle Varianzanalyse: 10-km-Zeit, Altersgruppen F30/F35/F40 ($n = 936$)

\vspace{2mm}
```{r, echo = T, eval = T}
# Daten extrahieren
D$T10k = t10_min                                        # 10-km-Zeit [min]

# Einfaktorielle Varianzanalyse
y1     = D$T10k[D$AgeGroup == "F30"]                    # Datenpunkte Gruppe 1
y2     = D$T10k[D$AgeGroup == "F35"]                    # Datenpunkte Gruppe 2
y3     = D$T10k[D$AgeGroup == "F40"]                    # Datenpunkte Gruppe 3
n      = length(y1) + length(y2) + length(y3)           # Anzahl Datenpunkte
p      = 3                                              # Anzahl Gruppen
y      = matrix(c(y1, y2, y3), nrow = n)                # Datenvektor
```

```{r, echo = F, eval = T}
# Einfaktorielle Varianzanalyse
n1     = length(y1)                                     # Anzahl Datenpunkte Gr. 1
n2     = length(y2)                                     # Anzahl Datenpunkte Gr. 2
n3     = length(y3)                                     # Anzahl Datenpunkte Gr. 3
X      = matrix(c(rep(1,n),                             # Designmatrix vollständiges Modell
                  c(rep(0,n1), rep(1,n2), rep(0,n3)),
                  c(rep(0,n1), rep(0,n2), rep(1,n3))), nrow = n)
X_0    = X[,1]                                          # Designmatrix reduziertes Modell

# Evaluation der F-Teststatistik
beta_hat   = solve(t(X)   %*% X)   %*% t(X)   %*% y     # Betaparameterschätzer vollständiges Modell
beta_hat_0 = solve(t(X_0) %*% X_0) %*% t(X_0) %*% y     # Betaparameterschätzer reduziertes Modell
eps_hat    = y - X   %*% beta_hat                       # Residuenvektor vollständiges Modell
eps_hat_0  = y - X_0 %*% beta_hat_0                     # Residuenvektor reduziertes Modell
SQT        = t(eps_hat_0) %*% eps_hat_0                 # total sum of squares
SQW        = t(eps_hat)   %*% eps_hat                   # within sum of squares
SQB        = SQT - SQW                                  # between sum of squares
DFB        = p - 1                                      # between degrees of freedom
DFW        = n - p                                      # within degrees of freedom
MSB        = SQB/DFB                                    # mean between sum of squares
MSW        = SQW/DFW                                    # mean within sum of squares
Eff        = MSB/MSW                                    # F-Teststatistik
pval       = 1 - pf(Eff, p-1, n-p)                      # p-Wert

# Ausgabe der Ergebnisse
cat(  "Datenpunkte pro Gruppe            : ", c(n1,n2,n3),
    "\nBetaparameterschätzer             : ", round(beta_hat,   digits = 3),
    "\nBetaparameterschätzer red. Mod.   : ", round(beta_hat_0, digits = 3),
    "\nresiduelle Quadratsumme           : ", round(SQW,        digits = 3),
    "\nresiduelle Quadratsumme red. Mod. : ", round(SQT,        digits = 3),
    "\nF-Teststatistik                   : ", round(Eff,        digits = 3),
    "\np-Wert                            : ", round(pval,       digits = 3), "\n\n")
```


<!-- Einfaktorielle Varianzanalyse (Männer) -->
# Anwendung/Praxis

\small
Einfaktorielle Varianzanalyse: 10-km-Zeit, Altersgruppen M30/M35/M40 ($n = 1323$)

\vspace{2mm}
```{r, echo = T, eval = T}
# Daten extrahieren
D$T10k = t10_min                                        # 10-km-Zeit [min]

# Einfaktorielle Varianzanalyse
y1     = D$T10k[D$AgeGroup == "M30"]                    # Datenpunkte Gruppe 1
y2     = D$T10k[D$AgeGroup == "M35"]                    # Datenpunkte Gruppe 2
y3     = D$T10k[D$AgeGroup == "M40"]                    # Datenpunkte Gruppe 3
n      = length(y1) + length(y2) + length(y3)           # Anzahl Datenpunkte
p      = 3                                              # Anzahl Gruppen
y      = matrix(c(y1, y2, y3), nrow = n)                # Datenvektor
```

```{r, echo = F, eval = T}
# Einfaktorielle Varianzanalyse
n1     = length(y1)                                     # Anzahl Datenpunkte Gr. 1
n2     = length(y2)                                     # Anzahl Datenpunkte Gr. 2
n3     = length(y3)                                     # Anzahl Datenpunkte Gr. 3
X      = matrix(c(rep(1,n),                             # Designmatrix vollständiges Modell
                  c(rep(0,n1), rep(1,n2), rep(0,n3)),
                  c(rep(0,n1), rep(0,n2), rep(1,n3))), nrow = n)
X_0    = X[,1]                                          # Designmatrix reduziertes Modell

# Evaluation der F-Teststatistik
beta_hat   = solve(t(X)   %*% X)   %*% t(X)   %*% y     # Betaparameterschätzer vollständiges Modell
beta_hat_0 = solve(t(X_0) %*% X_0) %*% t(X_0) %*% y     # Betaparameterschätzer reduziertes Modell
eps_hat    = y - X   %*% beta_hat                       # Residuenvektor vollständiges Modell
eps_hat_0  = y - X_0 %*% beta_hat_0                     # Residuenvektor reduziertes Modell
SQT        = t(eps_hat_0) %*% eps_hat_0                 # total sum of squares
SQW        = t(eps_hat)   %*% eps_hat                   # within sum of squares
SQB        = SQT - SQW                                  # between sum of squares
DFB        = p - 1                                      # between degrees of freedom
DFW        = n - p                                      # within degrees of freedom
MSB        = SQB/DFB                                    # mean between sum of squares
MSW        = SQW/DFW                                    # mean within sum of squares
Eff        = MSB/MSW                                    # F-Teststatistik
pval       = 1 - pf(Eff, p-1, n-p)                      # p-Wert

# Ausgabe der Ergebnisse
cat(  "Datenpunkte pro Gruppe            : ", c(n1,n2,n3),
    "\nBetaparameterschätzer             : ", round(beta_hat,   digits = 3),
    "\nBetaparameterschätzer red. Mod.   : ", round(beta_hat_0, digits = 3),
    "\nresiduelle Quadratsumme           : ", round(SQW,        digits = 3),
    "\nresiduelle Quadratsumme red. Mod. : ", round(SQT,        digits = 3),
    "\nF-Teststatistik                   : ", round(Eff,        digits = 3),
    "\np-Wert                            : ", round(pval,       digits = 3), "\n\n")
```


<!-- Stichprobenverzerrung -->
# Anwendung/Praxis

\textcolor{darkblue}{Stichprobenverzerrung (\textit{selection bias})}

\vspace{2mm}
```{r, echo = F, out.width = "50%", fig.align = "center"}
knitr::include_graphics("Abbildungen/selection_bias.png")
```

\vspace{2mm}
\center
\small
Webcomic zur Stichprobenverzerrung

\footnotesize
(Quelle: [*xkcd*: "Selection Bias" (#2618)](https://xkcd.com/2618/); Lizenz: CC-BY-NC 2.5).


<!-- Überlebendenverzerrung -->
# Anwendung/Praxis

\textcolor{darkblue}{Überlebendenverzerrung (\textit{survivorship bias})}

\vspace{2mm}
```{r, echo = F, out.width = "60%", fig.align = "center"}
knitr::include_graphics("Abbildungen/survivorship_bias.png")
```

\vspace{2mm}
\center
\small
Hypothetische Verteilung von Einschusslöchern bei zurückgekehrten Flugzeugen

\footnotesize
(Quelle: [*Wikimedia Commons*: "Survivorship-bias.svg"](https://commons.wikimedia.org/wiki/File:Survivorship-bias.svg); Lizenz: CC-BY-SA 4.0; vgl. @mangel1984)


<!-- Anhang: Selbstkontrollfragen -->
#

\large
\setstretch{2.5}
\vfill

Anwendungsszenario

Modellformulierung

Modellschätzung

Modellevaluation

Anwendung/Praxis

**Selbstkontrollfragen**

\vfill


<!-- Selbstkontrollfragen -->
# Selbstkontrollfragen

\vspace{2mm}
\footnotesize
\setstretch{1.5}
\begin{enumerate}

\item Erläutern Sie das Anwendungsszenario einer einfaktoriellen Varianzanalyse (EVA).
\item Geben Sie die Definition des EVA-Modells in Erwartungswertparameterdarstellung wieder.
\item Geben Sie die strukturelle Form des EVA-Modells in Effektdarstellung wieder.
\item Erläutern Sie die Motivation für die Reparameterisierung des EVA-Modells.
\item Welche Bedeutung haben die Parameter $\mu_0,\alpha_2,...,\alpha_p$ in der Effektparameterdarstellung des EVA-Modells?
\item Warum gibt es bei $p$ Gruppen eines EVA-Szenarios nur die $p-1$ Effektparameter $\alpha_2,...,\alpha_p$?
\item Geben Sie die Designmatrixform des EVA-Modells in Effektdarstellung wieder.
\item Formulieren Sie die Designmatrix eines EVA-Modells mit $n_i = 3$ und $p = 2$.
\item Formulieren Sie die Designmatrix eines EVA-Modells mit $n_i = 2$ und $p = 5$.
\item Geben Sie das Theorem zur Betaparameterschätzung im EVA-Modell wieder.
\item Mit welchen deskriptiven Statistiken werden die Parameter $\mu_0,\alpha_2,...,\alpha_p$ eines EVA-Modells geschätzt?

\end{enumerate}
\vfill


<!-- Selbstkontrollfragen -->
# Selbstkontrollfragen

\vspace{2mm}
\footnotesize
\setstretch{1.5}
\begin{enumerate}
\setcounter{enumi}{11}

\item Geben Sie das Theorem zur Quadratsummenzerlegung bei einfaktorieller Varianzanalyse wieder.
\item Erläutern Sie die Begriffe total sum of squares, between sum of squares und within sum of squares.
\item Geben Sie die Definition des Effektstärkenmaßes $\eta^2$ an.
\item Wann nimmt das Effektstärkenmaß $\eta^2$ der EVA seinen Minimalwert an und wie lautet dieser?
\item Wann nimmt das Effektstärkenmaß $\eta^2$ der EVA seinen Maximalwert an und wie lautet dieser?
\item Geben Sie das Theorem zur F-Teststatistik der EVA wieder.
\item Erläutern Sie die Begriffe mean between sum of squares und mean within sum of squares.
\item Geben Sie das Theorem zum Zusammenhang von Effektstärkenmaß $\eta^2$  und F-Teststatistik der EVA wieder.
\item Geben Sie die Definition des EVA-F-Test wieder.
\item Erläutern sie die Null- und Alternativhypothesen des EVA-F-Tests.
\item Geben Sie das Theorem zur Testumfangkontrolle der EVA wieder.
\item Skizzieren Sie den Beweis zur Testumfangkontrolle der EVA.
\item Geben Sie den p-Wert zum F-Test der EVA wieder.

\end{enumerate}
\vfill


<!-- Literatur -->
# Referenzen
