---
fontsize: 8pt
bibliography: 12_Referenzen.bib
citation_package: natbib
output:
  beamer_presentation:
    keep_tex: true
    includes:
      in_header: 12_header.tex
---

```{r, include = F}
source("12_R_common.R")
```

#  {.plain}
\center
```{r, echo = FALSE, out.width = "20%"}
knitr::include_graphics("12_Abbildungen/mvda_12_otto.png")
```

\vspace{2mm}

\Huge
Multivariate Datenanalyse
\vspace{6mm}

\Large
MSc Psychologie WiSe 2021/22

\vspace{6mm}
\large
Prof. Dr. Dirk Ostwald

# {.plain}
\vfill
\center
\huge
\textcolor{black}{(12) Kanonische Korrelationsanalyse}
\vfill


#
\setstretch{3}
\vfill
\large
**Vorbemerkungen**

Modellformulierung

Modellschätzung

Modellevaluation


\vfill

# Vorbemerkungen

\textcolor{darkblue}{Datenanalyseszenarien}

\small
\renewcommand{\arraystretch}{2}
\begin{tabular}{lll}
UV
& AV
& Datenanalysemethoden
\\\hline
Univariat
& Univariat
& Korrelation, Einfache Regression,T-Tests
\\
Univariat
& Multivariat
& T$^2$-Tests, MANOVA
\\
Multivariat
& Univariat
& Multiple Korrelation, Multiple Regression, ALM, SVMs, NNs
\\
Multivariat
& Multivariat
& Kanonische Korrelation, MALM, NNs
\end{tabular}

# Vorbemerkungen
\textcolor{darkblue}{Datenanalyseszenarien}
\vspace{2mm}

\small
\begin{minipage}{0.5\linewidth}
\center
\begin{tabular}{c|c}
UV         & AV          \\
$x_1$      & $y_1$       \\\hline
$x_{11}$   & $y_{11}$    \\
$x_{12}$   & $y_{12}$    \\
$x_{13}$   & $y_{13}$    \\
$\vdots$   & $\vdots$    \\
$\vdots$   & $\vdots$    \\
$\vdots$   & $\vdots$    \\
$x_{1n}$   & $y_{1n}$    \\
\end{tabular}
\vspace{2mm}

Korrelation

Einfache Regression

T-Tests
\end{minipage}
\begin{minipage}{0.5\linewidth}
\center
\begin{tabular}{ccc|c}
          & UV        &          & AV          \\
$x_1$     & $\cdots$  & $x_m$    & $y_1$       \\\hline
$x_{11}$  & $\cdots$  & $x_{m1}$ & $y_{11}$    \\
$x_{12}$  & $\cdots$  & $x_{m2}$ & $y_{12}$    \\
$x_{13}$  & $\cdots$  & $x_{m3}$ & $y_{13}$    \\
$\vdots$  & $\vdots$  & $\vdots$ & $\vdots$    \\
$\vdots$  & $\vdots$  & $\vdots$ & $\vdots$    \\
$\vdots$  & $\vdots$  & $\vdots$ & $\vdots$    \\
$x_{1n}$  & $\cdots$  & $x_{mn}$ & $y_{1n}$    \\
\end{tabular}
\vspace{2mm}

Multiple Korrelation

Multiple Regression

ALM, SVMs, NNs

\end{minipage}


# Vorbemerkungen
\textcolor{darkblue}{Datenanalyseszenarien}

\small
\begin{minipage}{0.4\linewidth}
\center
\begin{tabular}{c|ccc}
UV        &          & AV        &             \\
$x_1$     & $y_1$    & $\cdots$  & $y_m$       \\\hline
$x_{11}$  & $y_{12}$ & $\cdots$  & $y_{m1}$    \\
$x_{12}$  & $y_{13}$ & $\cdots$  & $y_{m2}$    \\
$x_{13}$  & $y_{14}$ & $\cdots$  & $y_{m3}$    \\
$\vdots$  & $\vdots$ & $\vdots$  & $\vdots$    \\
$\vdots$  & $\vdots$ & $\vdots$  & $\vdots$    \\
$\vdots$  & $\vdots$ & $\vdots$  & $\vdots$    \\
$x_{1n}$  & $y_{1n}$ & $\cdots$  & $y_{mn}$    \\
\end{tabular}
\vspace{2mm}

T$^2$-Tests

MANOVA

\end{minipage}\hspace{2mm}
\begin{minipage}{0.4\linewidth}
\center
\begin{tabular}{ccc|ccc}
          & UV        &             &            & AV        &             \\
$x_1$     & $\cdots$  & $x_{m_x}$   & $y_1$      & $\cdots$  & $y_{m_y}$   \\\hline
$x_{11}$  & $\cdots$  & $x_{m_x1}$  & $y_{11}$   & $\cdots$  & $y_{m_y1}$  \\
$x_{12}$  & $\cdots$  & $x_{m_x2}$  & $y_{12}$   & $\cdots$  & $y_{m_y2}$  \\
$x_{13}$  & $\cdots$  & $x_{m_x3}$  & $y_{13}$   & $\cdots$  & $y_{m_y3}$  \\
$\vdots$  & $\vdots$  & $\vdots$    & $\vdots$   & $\vdots$  & $\vdots$    \\
$\vdots$  & $\vdots$  & $\vdots$    & $\vdots$   & $\vdots$  & $\vdots$    \\
$\vdots$  & $\vdots$  & $\vdots$    & $\vdots$   & $\vdots$  & $\vdots$    \\
$x_{1n}$  & $\cdots$  & $x_{m_xn}$  & $y_{1n}$   & $\cdots$  & $y_{m_yn}$  \\
\end{tabular}
\vspace{2mm}

Kanonische Korrelation

MALM, NNs
\end{minipage}

# Vorbemerkungen
\textcolor{darkblue}{Überblick und Notation}
\footnotesize

Wir folgen in der Darstellung @mardia_1979, Chapter 10.

Die Datenvektoren $x_{1i},...,x_{m_xi}, i = 1,...,n$ werden als u.i.v. Realisierungen eines Zufallsvektors $X$ interpretiert.
Die Datenvektoren $y_{1i},...,y_{m_yi}, i = 1,...,n$ werden als u.i.v. Realisierungen eines Zufallsvektors $Y$ interpretiert.

Die "erste kanonische Korrelation" ist die maximale Korrelation von Linearkombinationen
von $X$ und $Y$; wir bezeichnen die Linearkombinationen von $X$ und $Y$ mit
Vektoren $a \in \mathbb{R}^{m_x}$ und $b \in \mathbb{R}^{m_y}$ mit
\begin{equation}
\xi      = a^TX = a_1X_1 + \cdots + a_{m_x}X_{m_x}
\mbox{ und }
\upsilon = b^TY =  b_1Y_1 + \cdots + b_{m_y}Y_{m_y}
\end{equation}
$\xi$ und $\upsilon$ sind dann als Linearkombinationen von Zufallsvariablen selbst Zufallsvariablen;
Die Korrelation von $\xi$ und $\upsilon$ bezeichnen wir mit $\rho(\xi,\upsilon)$

Wenn die Zufallsvektoren $X$ als unabhängige Variable und $Y$ als abhängige Variable
interpretiert werden, dann kann $\xi= a^TX$ als "bester Prädiktor" und
$\upsilon  = b^TY$ als "am besten prädizierbares Kriterium" interpretiert werden.
Kanonische Korrelationsanalyse fragt damit nach Parametern $a \in \mathbb{R}^{m_x}$
und $b\in \mathbb{R}^{m_y}$ für die $\rho(\xi,\upsilon)$ maximal ist.

Für Skalare $\alpha,\beta\in \mathbb{R}$ sind die Korrelationen $\rho\left(a^TX,b^TY\right)$ und
$\rho\left(\left(\alpha a^T\right)X,\left(\beta b^T\right)Y\right)$ allerdings
identisch (siehe unten). Man sucht deshalb Parameter
$a \in \mathbb{R}^{m_x}$ und $b\in \mathbb{R}^{m_y}$  für die $\rho(\xi,\upsilon)$ maximal ist und für die
$a^TX$ und $b^TY$ jeweils eine Varianz von 1 haben, also $\mathbb{V}(\xi) = \mathbb{V}(\upsilon) = 1$ gilt.


# Vorbemerkungen

\textcolor{darkblue}{Wahrscheinlichkeitstheorie}

\small
\begin{definition}[Korrelation]
Die Korrelation zweier Zufallsvariablen $\xi$ und $\upsilon$ ist definiert als
\begin{equation}
\rho(\xi,\upsilon) := \frac{\mathbb{C}(\xi,\upsilon)}{\mathbb{S}(\xi)\mathbb{S}(\upsilon)},
\end{equation}
wobei
\begin{equation}
\mathbb{C}(\xi,\upsilon) := \mathbb{E}((\xi-\mathbb{E}(\xi))(\upsilon - \mathbb{E}(\upsilon)))
\end{equation}
die Kovarianz von $\xi$ und $\upsilon$ und
\begin{equation}
\mathbb{S}(\xi) := \sqrt{\mathbb{V}(\xi)} = \sqrt{\mathbb{C}(\xi,\xi)}
\mbox{ und }
\mathbb{S}(\upsilon) := \sqrt{\mathbb{V}(\upsilon)} = \sqrt{\mathbb{C}(\upsilon,\upsilon)}
\end{equation}
die Standardabweichungen von $\xi$ und $\upsilon$, respektive, bezeichnen.
\end{definition}

\footnotesize
Bemerkungen

* $\rho \in [-1,1]$.
* Wenn $\xi$ und $\upsilon$ unabhängig sind, dann gilt $\rho(\xi,\upsilon) = 0$.
* $\rho(\xi,\upsilon)$ misst dem Grad des linear-affinen Zusammenhangs $\upsilon = a\xi + b$.

# Vorbemerkungen
\textcolor{darkblue}{Wahrscheinlichkeitstheorie}

\small
\begin{theorem}[Kovarianz und Korrelation bei linear-affinen Transformationen]
\justifying
\normalfont
$\rho(\xi,\upsilon)$ sei die Korrelation von $\xi$ und $\upsilon$ und es seien
$\alpha,\beta,\gamma,\delta \in \mathbb{R}$. Dann gelten
\begin{equation}
\mathbb{C}(\alpha\xi + \beta, \gamma\upsilon + \delta) = \alpha\beta\mathbb{C}(\xi,\upsilon)
\end{equation}
und
\begin{equation}
\rho(\alpha\xi + \beta, \gamma\upsilon + \delta) = \rho(\xi,\upsilon).
\end{equation}
\end{theorem}

\footnotesize
Bemerkung

* Die Varianzen von $a^TX$ und $b^Ty$ und die Varianzen von Linearkombinationen
von $X$ und $Y$ mit beliebigen skalaren Vielfachen von $a$ und $b$ sind im Sinne
der ersten Aussage zur Kovarianz verschieden.

* Insbesondere gilt einen $m_x$-dimensionalen Zufallsvektoren $X$ und einen
$m_y$-dimensionalen Zufallsvektor $Y$, $a \in \mathbb{R}^{m_x},b \in \mathbb{R}^{m_y}$,
$\alpha,\beta \in \mathbb{R}$ und die Linearkombinationen
$\xi := a^TX$ und $\upsilon := b^TY$ also auch
\begin{align}
\begin{split}
\rho(\xi,\upsilon) & = \rho(\alpha\xi,\beta\upsilon)  \\
\Leftrightarrow
\rho(a^TX, b^TY)   & = \rho(\alpha(a^TX),\beta(b^TY)) \\
\Leftrightarrow
\rho(a^TX, b^TY)   & = \rho((\alpha a^T)X),(\beta b^T)Y).
\end{split}
\end{align}
Die Korrelation von $a^TX$ und $b^Ty$ und die Korrelationen von Linearkombinationen
von $X$ und $Y$ mit beliebigen skalaren Vielfachen von $a$ und $b$ sind also gleich.


# Vorbemerkungen
\textcolor{darkblue}{Wahrscheinlichkeitstheorie}

\footnotesize
\underline{Beweis}

Es gilt zunächst
\begin{align}
\begin{split}
\mathbb{C}(\alpha\xi+\beta,c\upsilon+d)
& = \mathbb{E}((\alpha\xi+\beta-\mathbb{E}(\alpha\xi+\beta))(\gamma\upsilon+\delta-\mathbb{E}(\gamma\upsilon+\delta)))    \\
& = \mathbb{E}((\alpha\xi+\beta-\alpha\mathbb{E}(\xi)-\beta)(\gamma\upsilon+\delta-\gamma\mathbb{E}(\gamma\upsilon)-\delta))   \\
& = \mathbb{E}(\alpha(\xi-\mathbb{E}(\xi))(\gamma(\upsilon -\gamma\mathbb{E}(\upsilon)))             \\
& = \mathbb{E}(\alpha\gamma((\xi-\mathbb{E}(\xi))(\upsilon -\gamma\mathbb{E}(\upsilon))))            \\
&  = \alpha\gamma\mathbb{C}(\xi,\upsilon)
\end{split}
\end{align}
Also folgt
\begin{align}
\begin{split}
\rho(\alpha\xi + \beta, \gamma\upsilon + \delta)
& = \frac{\mathbb{C}(\alpha\xi+\beta,\gamma\upsilon+\delta)}{\sqrt{\mathbb{V}(\alpha\xi+\beta)}\sqrt{\mathbb{V}(\gamma\upsilon+\delta)}} \\
& = \frac{\alpha\gamma\mathbb{C}(\xi,\upsilon)}{\sqrt{\alpha^2\mathbb{V}(\xi)}\sqrt{\gamma^2\mathbb{V}(\upsilon)}}     \\
& = \frac{\alpha\gamma\mathbb{C}(\xi,\upsilon)}{\alpha\mathbb{S}(\xi)c\mathbb{S}(\upsilon)}         \\
& = \frac{\mathbb{C}(\xi,\upsilon)}{\mathbb{S}(\xi)\mathbb{S}(\upsilon)}             \\
& = \rho(\xi,\upsilon).
\end{split}
\end{align}


# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}
\footnotesize
\begin{definition}[Symmetrische Quadratwurzel einer Matrix]
\justifying
$A \in \mathbb{R}^{m \times m}$ sei eine invertierbare symmetrische Matrix mit
positiven Eigenwerten. Dann sind für $r \in \mathbb{N}^0$ und $s \in \mathbb{N}$
die rationalen Potenzen von $A$ einer orthonormalen Matrix $Q \in \mathbb{R}^{m \times m}$
der Eigenvektoren von $A$ und einer Diagonalmatrix  $\Lambda = \mbox{diag}(\lambda_i) \in \mathbb{R}^{m \times m}$
der zugehörigen Eigenwerte $\lambda_1,...,\lambda_m$ von $A$ definiert als
\begin{equation}
A^{r/s} = Q \Lambda^{r/s} Q^T \mbox{ mit } \Lambda^{r/s} = \mbox{diag}\left(\lambda_i^{r/s}\right).
\end{equation}
Der Spezialfall $r:= 1, s := 2$ wird als symmetrische Quadratwurzel von $A$ bezeichnet
und hat die Form
\begin{equation}
A^{1/2} = Q\Lambda^{1/2}Q^T \mbox{ mit } \Lambda^{1/2} = \mbox{diag}\left(\lambda_i^{1/2}\right).
\end{equation}
\end{definition}
\setstretch{1}
Bemerkungen

* Offenbar gilt
\begin{equation}
\left(A^{1/2} \right)^2 = Q\Lambda^{1/2}Q^TQ\Lambda^{1/2}Q^T = Q\Lambda^{1/2}\Lambda^{1/2}Q^T =  Q\Lambda Q^T = A.
\end{equation}

* Weiterhin gilt
\begin{equation}
\left(A^{-1/2} \right)^2 = Q\Lambda^{-1/2}Q^TQ\Lambda^{-1/2}Q^T = Q\Lambda^{-1}Q^T = A^{-1}.
\end{equation}
Die vorletzte Gleichung mag überraschen, aber es gilt ja zum Beispiel
\begin{equation}
4^{-1/2}\cdot 4^{-1/2} = \frac{1}{\sqrt{4}}\cdot\frac{1}{\sqrt{4}} = \frac{1}{4} = 4^{-1}.
\end{equation}

# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

\footnotesize
\begin{theorem}[Eigenwerte und Eigenvektoren von Matrixprodukten]
\justifying
\normalfont
Für $A \in \mathbb{R}^{n \times m}$ und $B \in \mathbb{R}^{m \times n}$ sind
die Eigenwerte von $AB \in \mathbb{R}^{n \times n}$ und $BA \in \mathbb{R}^{m \times m}$
gleich. Weiterhin gilt, dass für einen Eigenvektor $v$ zu einem von Null
verschiedenen Eigenwert $\lambda$ von $AB$ $w := Bv$ ein Eigenvektor von $BA$ ist.
\end{theorem}

Bemerkungen

* Für einen Beweis verweisen wir auf @mardia_1979, S. 468.

\vspace{1mm}
\tiny
\setstretch{1.2}
```{r}
A   = matrix(1:6, nrow = 2,  byrow = T)           # Matrix A \in \mathbb{R}^{2 x 3}
B   = matrix(1:6, ncol = 2,  byrow = T)           # Matrix B \in \mathbb{R}^{3 x 2}
EAB = eigen(A %*% B)                              # Eigenanalyse von AB \in \mathbb{R}^{2 \times 2}
EBA = eigen(B %*% A)                              # Eigenanalyse von BA \in \mathbb{R}^{3 \times 3}
w   = B %*% EAB$vectors[,1]                       # Eigenvektor von BA
cat("Eigenwerte von AB :"  , EAB$values[1:2],
    "\nEigenwerte von BA :", EBA$values[1:2],
    "\nBAw mit w = Bv    :", B %*% A %*% w,
    "\nlw mit w = Bv     :", EBA$values[1] * w)

```

# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

\footnotesize
\begin{theorem}[Eigenwert und Eigenvektor eines Matrixvektorprodukts]
\justifying
\normalfont
Für $A \in \mathbb{R}^{n \times m}, B \in \mathbb{R}^{p \times n}, a \in \mathbb{R}^m$
und $b \in \mathbb{R}^p$ gilt, dass der einzige von Null verschiedene Eigenwert von
$Aab^TB \in \mathbb{R}^{n \times n}$ gleich $b^T BAa$ mit zugehörigem Eigenvektor $Aa$ ist.
\end{theorem}

Bemerkungen

* Für einen Beweis verweisen wir auf @mardia_1979, S. 468.

\vspace{1mm}
\tiny
\setstretch{1.2}
\vspace{1mm}
\tiny
\setstretch{1.2}
```{r}
A      = matrix(1:6, nrow = 2,  byrow = T)     # Matrix A \in \mathbb{R}^{2 x 3}
B      = matrix(1:8, ncol = 2,  byrow = T)     # Matrix B \in \mathbb{R}^{4 x 2}
a      = matrix(1:3, nrow = 3,  byrow = T)     # Vektor a \in \mathbb{R}^{3 x 1}
b      = matrix(1:4, nrow = 4,  byrow = T)     # Vektor b \in \mathbb{R}^{4 x 1}
EAabTB = eigen(A %*% a %*% t(b) %*% B)         # Eigenanalyse von Aab^TB \in \mathbb{R}^{4 x 4}
cat("Eigenwerte von AabTB :", EAabTB$values,
    "\nbTBAa                :", t(b) %*% B %*% A %*% a,
    "\nAa                   :", A %*% a,
    "\n(AabTB)Aa            :",(A %*% a %*% t(b) %*% B) %*% A %*% a,            # Mv
    "\n(bTBAa)Aa            :",as.vector((t(b) %*% B %*% A %*% a)) * (A %*% a)) # = \lambda v
```

# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}
\footnotesize
\begin{theorem}[Maximierung quadratischer Formen mit Nebenbedingungen]
\justifying
\normalfont
$A \in \mathbb{R}^{m \times m}, B \in \mathbb{R}^{m \times m} \mbox{ p.d.}$
seien symmetrische Matrizen und $\lambda_1$ sei der größte Eigenwert von $B^{-1}A$
mit assoziertem Eigenvektor $v_1 \in \mathbb{R}^m$. Dann ist $\lambda_1$ eine Lösung des
Optimierungsproblems
\begin{equation}\label{eq:opt_1}
\max_{x} x^TAx \mbox{ unter der Nebenbedingung } x^TBx = 1.
\end{equation}
\end{theorem}

Bemerkungen

* Das Theorem ist direkt durch die kanonische Korrelationsanalyse motiviert.
* Nach Wortlaut des Theorems gilt also
\begin{equation}
v_1 = \argmax_{x} x^TAx \mbox{ unter der Nebenbedingung } x^TBx = 1
\end{equation}
* Nach Wortlaut des Theorems gilt weiterhin
\begin{equation}
\lambda_1 = \max_{x} x^TAx \mbox{ unter der Nebenbedingung } x^TBx = 1.
\end{equation}


# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

\footnotesize
\underline{Beweis}

$B^{1/2}$ sei die symmetrische Quadratwurzel von $B$ und es sei
\begin{equation}
y := B^{1/2}x \Leftrightarrow x = B^{-1/2}y
\end{equation}
Dann kann mit der symmetrischen Matrix
\begin{equation}
K := B^{-1/2}AB^{-1/2} \in \mathbb{R}^{m \times m}
\end{equation}
das Optimierungsproblem \eqref{eq:opt_1} geschrieben werden als
\begin{equation}\label{eq:opt_2}
\max_{y} y^T K y
\mbox{ unter der Nebenbedingung }
y^Ty = 1.
\end{equation}
Dies gilt, weil
\begin{equation}
\max_{x} x^TAx
\Leftrightarrow
\max_{y} \left(B^{-1/2}y\right)^TA\left(B^{-1/2}y\right)
\Leftrightarrow
\max_{y} y^TB^{-1/2}AB^{-1/2}y
\Leftrightarrow
\max_{y} y^TKy
\end{equation}
und
\begin{equation}
x^TBx = 1
\Leftrightarrow y^T B^{-1/2}BB^{-1/2}y = 1
\Leftrightarrow y^Ty = 1.
\end{equation}
Weil $K$ eine symmetrische Matrix ist, existiert die Orthonormalzerlegung (vgl. (2) Matrizen)
\begin{equation}
K = Q\Lambda Q^T,
\end{equation}
wobei die Spalten der orthogonalen Matrix $Q$ die Eigenvektoren von $K$ und die
Diagonalemente von $\Lambda$ die zugehörigen Eigenwerte von $K$ sind.



# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

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

Mit der orthogonalen Matrix $Q$ aus obiger Orthornomalzerlegung sei nun
\begin{equation}
z := Q^Ty \Leftrightarrow y := Qz.
\end{equation}
Dann kann das Optimierungsproblem \eqref{eq:opt_2} geschrieben werden als
\begin{equation}\label{eq:opt_3}
\max_{z} \sum_{i = 1}^m \lambda_i z_i^2 \mbox{ unter der Nebenbedingung } z^Tz = 1,
\end{equation}
weil
\begin{equation}
\max_{y} y^TKy
\Leftrightarrow
\max_{z} (Qz)^TK(Qz)
\Leftrightarrow
\max_{z} z^TQ^TQ\Lambda Q^TQz
\Leftrightarrow
\max_{z} z^T\Lambda z
\Leftrightarrow
\max_{z} \sum_{i=1}^m \lambda_i z_i^2
\end{equation}
und
\begin{equation}
y^Ty = 1
\Leftrightarrow
(Qz)^T Qz = 1
\Leftrightarrow
z^T Q^TQz = 1
\Leftrightarrow
z^T z = 1.
\end{equation}


# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

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

Die Eigenwerte von $K$ seien nun absteigend sortiert, also $\lambda_1 \ge \cdots \ge \lambda_m$.
Dann gilt für das Optimierungsproblem \eqref{eq:opt_3}, dass
\begin{equation}
\max_{z} \sum_{i = 1}^m \lambda_i z_i^2 \le \lambda_1,
\end{equation}
weil
\begin{equation}
\max_{z} \sum_{i = 1}^m \lambda_i z_i^2
\le
\max_{z} \sum_{i = 1}^m \lambda_1 z_i^2
=
\lambda_1 \max_{z} \sum_{i = 1}^m z_i^2
=
\lambda_1
\end{equation}
wobei sich die letzte Gleichung aus der Nebenbedingung $z^Tz=1$ ergibt. Schließlich
gilt
\begin{equation}
\max_{z} \sum_{i = 1}^m \lambda_i z_i^2  = \lambda_1,
\end{equation}
für $z := e_1 = (1,0,...,0)^T$. Zusammenfassend heißt das, dass $z = e_1$ eine Lösung
des Optimierungsproblem \eqref{eq:opt_3} ist und das $\lambda_1$ das entsprechende
Maximum ist.

# Vorbemerkungen
\textcolor{darkblue}{Lineare Algebra}

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

Damit ergibt sich aber sofort, dass dann
\begin{equation}
y = Qz = Qe_1 =  q_1 \mbox{ und } x = B^{-1/2}q_1
\end{equation}
Lösungen der äquivalenten Optimierungsprobleme \eqref{eq:opt_2} und \eqref{eq:opt_1},
respektive, sind. Nach Konstruktion ist $q_1$ ein Eigenvektor von $B^{-1/2}AB^{-1/2}$ und nach
obigem Theorem zu Eigenwerten und Eigenvektoren von Matrixprodukten damit auch
ein Eigenvektor von
\begin{equation}
B^{-1/2}B^{-1/2}A = B^{-1}A
\end{equation}
und die zugehörigen Eigenwerte sind gleich. Damit aber folgt, dass der
größte Eigenwert von $B^{-1}A$ und sein assoziierter Eigenvektor eine Lösung von
\begin{equation}
\max_{x} x^TAx \mbox{ unter der Nebenbedingung } x^TBx = 1.
\end{equation}
ist.
$\hfill\Box$



#
\setstretch{3}
\vfill
\large
Vorbemerkungen

**Modellformulierung**

Modellschätzung

Modellevaluation


\vfill

# Modellformulierung
\textcolor{darkblue}{Überblick}
\setstretch{1.6}

\small
Zur Entwicklung der kanonischen Korrelationsanalyse werden $X$ und $Y$ als
\begin{equation}
Z := \begin{pmatrix} X \\ Y \end{pmatrix}
\end{equation}
zusammengefasst.

Wir nehmen durchgängig an, dass $\mathbb{E}(Z) = 0_{m}$ mit $m = m_x + m_y$.

Der mathematische Fokus ist auf der Kovarianzmatrix $\mathbb{C}(Z)$.

\footnotesize
* Kovarianzen von Linearkombinationen von $X$ und $Y$ ergeben sich aus Matrixprodukten von $\mathbb{C}(Z)$.

* Die Matrixtheoreme aus den Vorbemerkungen können auf diese Matrixprodukte angewendet werden.

\small
Generell wird im folgenden ein restringierter Optimierungsansatz mithilfe der 
Lagrangefunktion zugunsten der Eigenanalyse von Matrixprodukten supprimiert. 
Für den Lagrangeansatz, siehe zum Beispiel @anderson_2003, Kapitel 12.


# Modellformulierung
\footnotesize
\begin{theorem}[Kovarianzmatrizen von Zufallsvektoren]
\normalfont
\justifying
$Z$ sei ein $m$-dimensionaler Zufallsvektor mit Erwartungswert $\mathbb{E}(Z) = 0_m$ und
es sei
\begin{equation}
Z = \begin{pmatrix} X \\ Y \end{pmatrix}
\mbox{ mit }
\mathbb{E}(Z)  := 0_m
\end{equation}
ein $m_x + m_y$-dimensionaler Zufallsvektor und sein Erwartungswertvektor,
respektive. Dann kann die $m \times m$ Kovarianzmatrix  $Z$ geschrieben werden als
\begin{equation}
\mathbb{C}(Z) =
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\in \mathbb{R}^{m \times m}
\end{equation}
wobei
\begin{align}
\begin{split}
\Sigma_{xx} & := \mathbb{E}\left(XX^T\right) \in \mathbb{R}^{m_x \times m_x}\\
\Sigma_{xy} & := \mathbb{E}\left(XY^T\right) \in \mathbb{R}^{m_x \times m_y}\\
\Sigma_{yx} & := \mathbb{E}\left(YX^T\right) \in \mathbb{R}^{m_y \times m_x}\\
\Sigma_{yy} & := \mathbb{E}\left(YY^T\right) \in \mathbb{R}^{m_x \times m_y}
\end{split}
\end{align}
\end{theorem}

# Modellformulierung
\footnotesize
\underline{Beweis}

\renewcommand{\arraystretch}{1.2}
Nach Definition der Kovarianzmatrix eines Zufallsvektors (vgl. (3) Wahrscheinlichkeitshtheorie) gilt
\begin{align}
\begin{split}
\mathbb{C}(Z)
& = \mathbb{E}\left((Z - \mathbb{E}(Z))(Z - \mathbb{E}(Z))^T \right) \\
& = \mathbb{E}\left((Z - 0_m)(Z - 0_m)^T \right) \\
& = \mathbb{E}\left(ZZ^T\right)\\
& = \mathbb{E}\left(\begin{pmatrix} X \\ Y \end{pmatrix} \begin{pmatrix} X^T & Y^T \end{pmatrix} \right) \\
& = \mathbb{E}\left(\begin{pmatrix} XX^T & XY^T \\ YX^T & YY^T \end{pmatrix}\right)
\\
& =
\begin{pmatrix}
\mathbb{E}\left(XX^T\right) & \mathbb{E}\left(XY^T\right) \\
\mathbb{E}\left(YX^T\right) & \mathbb{E}\left(YY^T\right)
\end{pmatrix}
\\
& =
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\end{split}
\end{align}
$\hfill\Box$

# Modellformulierung
\footnotesize
\begin{theorem}[Linearkombinationen von Zufallsvektorpartitionen]
\justifying
\normalfont
Es sei
\begin{equation}
Z = \begin{pmatrix} X \\ Y \end{pmatrix}
\mbox{ mit }
\mathbb{E}(X) = 0_m
\mbox{ und }
\mathbb{C}(Z) =
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\end{equation}
ein $m$-dimensionaler partitionierter Zufallsvektor sowie sein
Erwartungswertvektor und seine Kovarianzmatrix, respektive. Weiterhin seien für
$a \in \mathbb{R}^{m_x}$ und $b\in\mathbb{R}^{m_y}$
die Zufallsvariablen
\begin{equation}
\xi := a^T X \mbox{ und } \upsilon := b^T Y
\end{equation}
als Linearkombinationen der Komponenten von $X$ und $Y$ definiert.
Dann gelten
\begin{itemize}
\item[(1)] $\mathbb{V}(\xi) = a^T\Sigma_{xx}a$
\item[(2)] $\mathbb{V}(\upsilon) = b^T\Sigma_{yy}b$
\item[(2)] $\rho(\xi,\upsilon) = a^T \Sigma_{xy}b$, wenn  $\mathbb{V}(\xi) = 1$ und $\mathbb{V}(\upsilon) = 1$.
\end{itemize}
\end{theorem}

Bemerkungen

* Die Varianz der Zufallsvariable $a^TX$ ergibt sich als "doppelte Linearkombination" von $\Sigma_{xx}$.
* Die Varianz der Zufallsvariable $b^TY$ ergibt sich als "doppelte Linearkombination" von $\Sigma_{yy}$.
* Die Korrelation der Zufallsvariablen  $a^TX$ und $b^TY$ ergibt sich  "doppellte Linearkombination" von $\Sigma_{xy}$.

# Modellformulierung
\footnotesize
\underline{Beweis von (1) und (2)}

Wir betrachten zunächst die Varianz von $\xi$. Mit dem Varianzverschiebungssatz gilt
\begin{align}
\begin{split}
\mathbb{V}(\xi)
& = \mathbb{E}\left(\xi \xi \right) - \mathbb{E}(\xi)\mathbb{E}(\xi) \\
& = \mathbb{E}\left((a^TX) (a^TX)\right) - \mathbb{E}\left(a^TX\right)\mathbb{E}\left(a^TX\right) \\
& = \mathbb{E}\left((a^TX) (a^TX)^T\right) - \mathbb{E}\left(a^TX\right)\mathbb{E}\left(a^TX\right) \\
& = \mathbb{E}\left(a^TXX^Ta\right) - \mathbb{E}\left(a^TX\right)\mathbb{E}\left(a^TX\right) \\
& = a^T\mathbb{E}\left(XX^T\right)a - a^T\mathbb{E}(X)a^T\mathbb{E}(X) \\
& = a^T\mathbb{E}\left(XX^T\right)a - a^T0_{m_x}a^T0_{m_x} \\
& = a^T\Sigma_{xx}a. \\
\end{split}
\end{align}
Der Beweis zur Varianz von $\upsilon$ folgt dann analog.

# Modellformulierung
\footnotesize
\underline{Beweis von (3)}

Mit der Definition der Korrelation von Zufallsvariablen und mit $\mathbb{V}(\xi) = \mathbb{V}(\upsilon) = 1$
und dem Kovarianzverschiebungssatz gilt
\begin{align}
\begin{split}
\rho(\xi,\upsilon)
& = \frac{\mathbb{C}(\xi,\upsilon)}{\sqrt{\mathbb{V}(\xi)}\sqrt{\mathbb{V}(\upsilon)}} \\
& = \frac{\mathbb{C}(\xi,\upsilon)}{\sqrt{1}\sqrt{1}} \\
& = \mathbb{C}(\xi,\upsilon) \\
& = \mathbb{E}(\xi\upsilon) - \mathbb{E}(\xi)\mathbb{E}(\upsilon) \\
& = \mathbb{E}\left((a^T X)(b^TY)\right) - \mathbb{E}(a^TX)\mathbb{E}(b^TY) \\
& = \mathbb{E}\left((a^T X)(b^TY)^T\right) - \mathbb{E}(a^TX)\mathbb{E}(b^TY) \\
& = \mathbb{E}\left(a^T XY^Tb \right) - \mathbb{E}(a^TX)\mathbb{E}(b^TY) \\
& = a^T\mathbb{E}\left(XY^T \right)b - a^T\mathbb{E}(X)b^T\mathbb{E}(Y) \\
& = a^T\mathbb{E}\left(XY^T \right)b - a^T0_{m_x}b^T0_{m_y} \\
& = a^T\Sigma_{xy}b. \\
\end{split}
\end{align}
$\hfill\Box$


# Modellformulierung
\setstretch{1.2}
\vspace{1mm}
\footnotesize
\begin{definition}[Kanonische Koeffizientenvektoren, Variate, Korrelationen]
\justifying
Es seien
\begin{equation}
Z = \begin{pmatrix} X \\ Y \end{pmatrix}
\mbox{ mit }
\mathbb{E}(Z) := 0_m
\mbox{ und }
\mathbb{C}(Z) :=
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\in \mathbb{R}^{m \times m}
\end{equation}
ein $m$-dimensionaler partitionierter Zufallsvektor sowie sein Erwartungswert
und seine Kovarianzmatrix, respektive. Weiterhin sei
\begin{equation}
K := \Sigma_{xx}^{-1/2}\Sigma_{xy}\Sigma_{yy}^{-1/2} \in \mathbb{R}^{m_x \times m_y}
\end{equation}
mit der Singulärwertzerlegung
\begin{equation}
K = A \Lambda B^T,
\end{equation}
wobei
\begin{equation}
A       := \begin{pmatrix} \alpha_1 & \cdots & \alpha_k \end{pmatrix} \in \mathbb{R}^{m_x \times m_y}
\mbox{ und }
B       := \begin{pmatrix} \beta_1  & \cdots &  \beta_k \end{pmatrix} \in \mathbb{R}^{m_y \times m_y}
\end{equation}
die orthogonalen Matrix der Eigenvektoren von $KK^T$  und die orthogonale Matrix
der Eigenvektoren von $K^TK$, respektive, bezeichnen und
\begin{equation}
\Lambda := \mbox{diag}\left(\lambda^{1/2}_1,...,\lambda_k^{1/2}\right) \in \mathbb{R}^{m_y \times m_y},
\end{equation}
die Diagonalmatrix der Quadratwurzeln der zugehörigen Eigenvektoren bezeichnet.
Schließlich seien für $i = 1,...,k$
\begin{equation}
a_i := \Sigma_{xx}^{-1/2}\alpha_i  \in \mathbb{R}^{m_x} \mbox{ und } b_i := \Sigma_{yy}^{-1/2}\beta_i \in \mathbb{R}^{m_y}.
\end{equation}
Dann heißen für $i = 1,...,k$
\begin{itemize}
\itemsep0mm
\item[(1)] $a_i \in \mathbb{R}^{m_x}$ und $b_i \in \mathbb{R}^{m_y}$ die \textit{$i$ten kanonischen Koeffizientenvektoren},
\item[(2)] die Zufallsvektoren $\xi_i := a_i^TX$ und $\upsilon_i := b_i^TY$ die $i$ten \textit{$i$ten kanonischen Variaten} und
\item[(3)] $\rho_i := \lambda_i^{1/2}$ die \textit{$i$te kanonische Korrelation}.
\end{itemize}
\end{definition}

# Modellformulierung
\footnotesize
\begin{theorem}[Eigenschaften kanonischer Korrelationen und Variaten]
\normalfont
\justifying
Es seien
\begin{equation}
Z = \begin{pmatrix} X \\ Y \end{pmatrix}
\mbox{ mit }
\mathbb{E}(Z) := 0_m
\mbox{ und }
\mathbb{C}(Z) :=
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\in \mathbb{R}^{m \times m}
\end{equation}
ein $m$-dimensionaler partitionierter Zufallsvektor sowie sein Erwartungswert
und seine Kovarianzmatrix, respektive. Weiterhin seien für $i = 1,...,k$ die kanonischen Koeffizientenvektoren $a_i, b_i$,
die kanonischen Variaten $\xi,\upsilon_i$ und die kanonischen Korrelationen $\rho_i$
definiert wie oben. Dann gilt, dass für $1 \le r \le k$ das Maximum des $r$ten
restringierten Optimierungsproblems
\begin{equation}
\phi_r = \max_{a,b} a^T\Sigma_{xy}b
\end{equation}
unter den Nebenbedingungen
\begin{equation}
a^T\Sigma_{xx}a   = 1,
\quad
b^T\Sigma_{yy}b   = 1,
\quad
a_i^T\Sigma_{xx}a = 0 \mbox{ für } i = 1,...,r-1
\end{equation}
(1) den Wert $\phi_r = \rho_r$ hat und (2) bei $a = a_r$ und $b = b_r$ angenommen wird.
\end{theorem}
\vspace{-1mm}
\setstretch{1.2}
Bemerkungen
\vspace{-2mm}
\begin{itemize}
\item $\phi_1$ ist die größtmögliche Korrelation von $\xi = a^TX$ und $\upsilon = b^TY$ unter den Nebenbedingungen
\begin{itemize}
\footnotesize
\item[$\circ$] $\mathbb{V}(\xi) =  a^T\Sigma_{xx}a = 1$ und $\mathbb{V}(\upsilon) = b^T\Sigma_{yy}b = 1$
\end{itemize}
\item $\phi_r$ ist die größtmögliche Korrelation von $\xi = a^TX$ und $\upsilon = b^TY$  unter den Nebenbedingungen
\begin{itemize}
\footnotesize
\item[$\circ$] $\mathbb{V}(\xi) = a^T\Sigma_{xx}a = 1$ und $\mathbb{V}(\upsilon) = b^T\Sigma_{yy}b = 1$
\item[$\circ$]  $\mathbb{C}(\xi_i,\xi) = a_i^T\Sigma_{xx}a = 0$  für die ersten $i = 1,...,r-1$ kanonischen Variaten $\xi_i$
\end{itemize}
\end{itemize}


# Modellformulierung
\footnotesize
\underline{Beweis}

Wir betrachten das restringierte Optimierungsproblem
\tiny
\begin{equation}
\phi_r^2 = \max_{a,b} \left(a^T\Sigma_{xy}b\right)^2\,
\mbox{ u.d.N. }
a^T\Sigma_{xx}a   = 1,\,
b^T\Sigma_{yy}b   = 1,\,
a_i^T\Sigma_{xx}a = 0, i = 1,...,r-1
\end{equation}
\footnotesize
Wir folgen @mardia_1979, S. 284 und gehen schrittweise vor, d.h. wir lösen
das restringierte Optimierungsproblem
\tiny
\begin{equation}
\phi_r^2 = \max_{a} \left(\max_{b} \left(a^T\Sigma_{xy}b\right)^2 \mbox{ u.d.N.} b^T\Sigma_{yy}b   = 1\right)
\mbox{ u.d.N. } a^T\Sigma_{xx}a   = 1,\, a_i^T\Sigma_{xy}a = 0,  i = 1,...,r-1
\end{equation}
\footnotesize
von innen nach außen.

Schritt (1)

Wir wählen wir zunächst ein festes $a \in \mathbb{R}^m$ und betrachten das
restringierte Optimierungsproblem
\begin{equation}
\max_{b} \left(a^T\Sigma_{xy}b\right)^2
\mbox{ u.d.N. }
b^T\Sigma_{yy}b   = 1
\end{equation}
Dieses Optimierungsproblem kann geschrieben werden als
\begin{equation}\label{eq:kka_opt_1}
\max_{b} b^T\Sigma_{yx}aa^T\Sigma_{xy}b
\mbox{ u.d.N. }
b^T\Sigma_{yy}b = 1,
\end{equation}
weil gilt, dass
\begin{equation}
\left(a^T\Sigma_{xy}b\right)^2
= \left(a^T\Sigma_{xy}b\right)\left(a^T\Sigma_{xy}b\right)
= \left(a^T\Sigma_{xy}b\right)^T a^T\Sigma_{xy}b
= b^T\Sigma_{yx}aa^T\Sigma_{xy}b.
\end{equation}


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

Das Optimierungsproblem \eqref{eq:kka_opt_1} kann nun mithilfe des Theorems
zur Maximierung quadratischer Formen mit Nebenbedingen gelöst werden. Im Sinne
dieses Theorems setzen wir dazu
\begin{equation}
A := \Sigma_{yx}aa^T\Sigma_{xy} \mbox{ und } B := \Sigma_{yy}.
\end{equation}
Dann hat \eqref{eq:kka_opt_1} die Form
\begin{equation}\label{eq:kka_opt_2}
\max_{b} b^TAb
\mbox{ unter der Nebenbedingung }
b^TBb = 1,
\end{equation}
Das Maximum von \eqref{eq:kka_opt_2} entspricht nach dem Theorem  zur Maximierung
quadratischer Formen mit Nebenbedingungen dem größten Eigenwert von
\begin{equation}
B^{-1}A = \Sigma_{yy}^{-1}\Sigma_{yx}aa^T\Sigma_{xy}
\end{equation}
Der größte Eigenwert von $\Sigma_{yy}^{-1}\Sigma_{yx}aa^T\Sigma_{xy}$ wiederum
kann mithilfe des Theorems zum Eigenwert und Eigenvektor eines Matrixvektorprodukts
bestimmt werden. Im Sinne dieses Theorems setzen wir dazu
\begin{equation}
A := \Sigma_{yy}^{-1}\Sigma_{yx},\quad b := a,\quad B := \Sigma_{xy}
\end{equation}
und erhalten für den betreffenden Eigenwert
\begin{equation}
\lambda_a = b^TBAa = a^T\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}a.
\end{equation}
als Lösung (Maximum) des restringierten Optimierungsproblems
\begin{equation}
\max_{b} \left(a^T\Sigma_{xy}b\right)^2 \mbox{ u.d.N. } b^T\Sigma_{yy}b = 1
\end{equation}

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

Schritt (2)

Basierend auf Schritt (1) verbleibt die Lösung des restringierten Optimierungsproblem
\begin{equation}\label{eq:kka_opt_3}
\phi_r^2 = \max_{a} a^T\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}a
\mbox{ u.d.N. }
a^T\Sigma_{xx}a   = 1,\,
a_i^T\Sigma_{xx}a = 0,  i = 1,...,r-1
\end{equation}
Dazu halten wir zunächst fest, dass \eqref{eq:kka_opt_3} mit den Definitionen
von $\alpha_i$ und $K$ in der Definition der Kanonischen  Koeffizientenvektoren,
Variaten, und Korrelationen geschrieben werden kann als
\begin{equation}\label{eq:kka_opt_4}
\phi_r^2 = \max_{\alpha} \alpha^T KK^T \alpha
\mbox{ u.d.N. }
\alpha^T \alpha   = 1,\,
\alpha_i^T\alpha  = 0,  i = 1,...,r-1,
\end{equation}
denn
\tiny
\begin{align}
\begin{split}
\phi_r^2 & = \max_{a} a^T\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}a
\mbox{ u.d.N. }
a^T\Sigma_{xx}a   = 1,
a_i^T\Sigma_{xx}a = 0 \Leftrightarrow
\\
\phi_r^2 & = \max_{a} a^T\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}a
\mbox{ u.d.N. }
\alpha^T\Sigma_{xx}^{-1/2}\Sigma_{xx}\Sigma_{xx}^{-1/2}\alpha = 1,
\alpha^T_i\Sigma_{xx}^{-1/2}\Sigma_{xx}\Sigma_{xx}^{-1/2}\alpha = 0
\\
\phi_r^2 & = \max_{a} \alpha^T\Sigma_{xx}^{-1/2}\Sigma_{xy}\Sigma_{yy}^{-1}\Sigma_{yx}\Sigma_{xx}^{-1/2}\alpha
\mbox{ u.d.N. }
\alpha^T\alpha = 1,
\alpha^T_i\alpha = 0
\\
\phi_r^2 & = \max_{a} \alpha^T\Sigma_{xx}^{-1/2}\Sigma_{xy}\Sigma_{yy}^{-1/2}\Sigma_{yy}^{-1/2}\Sigma_{yx}\Sigma_{xx}^{-1/2}\alpha
\mbox{ u.d.N. }
\alpha^T\alpha = 1,
\alpha^T_i\alpha = 0
\\
\phi_r^2 & = \max_{a} \alpha^TKK^T\alpha
\mbox{ u.d.N. }
\alpha^T\alpha = 1,
\alpha^T_i\alpha = 0
\end{split}
\end{align}

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

Dabei sind nach der betreffenden Definition die $\alpha_i$ die Eigenvektoren
von $KK^T$ mit den $i = 1,...,r-1$ größten Eigenwerten. Nach dem Theorem zur
Maximierung quadratischer Formen mit Nebenbedingungen ist die Lösung von
\eqref{eq:kka_opt_4} der größte Eigenwert von $KK^T$ mit seinem assoziierten
Eigenvektor. Die Nebenbedingung $\alpha_i^T\alpha = 0$ schränkt diese Wahl auf
den $r$t-größten Eigenwert und seinen assoziierten Eigenvektor $\alpha_r$ ein.
Mit der Definition von Eigenwerten und Eigenvektoren gilt also
\begin{equation}
\phi_r^2 = \alpha_r^T KK^T \alpha_r = \alpha_r^T \lambda_r \alpha_r = \lambda_r \alpha_r^T\alpha_r = \lambda_r.
\end{equation}
Wir haben also gezeigt, dass das restringierte Optimierungsproblem des Theorems
den Maximumwert $\phi_r = \lambda_r^{1/2}$ hat. Es bleibt zu zeigen, dass dieser
Maximumwert für $a_r$ und $b_r$ angenommen wird.

Schritt (3)

Einsetzen von $a_r$ und $b_r$ in $a^T\Sigma_{xy}b$ ergibt mit
\begin{equation}
K = A\Lambda B^T
\Leftrightarrow KB = A\Lambda B^TB
\Leftrightarrow KB = A\Lambda
\Leftrightarrow K\beta_r = \alpha_r\lambda_r^{1/2}
\end{equation}
dass
\begin{equation}
a^T_r\Sigma_{xy}b_r
= \alpha_r^T\Sigma_{xx}^{-1/2}\Sigma_{xy}\Sigma_{yy}^{-1/2}\beta_r
= \alpha_r^TK\beta_r
= \alpha_r^T\alpha_r\lambda_r^{1/2}
= \rho_r
\end{equation}
Also nimmt $a^T\Sigma_{xy}b$ bei $a_r$ und $b_r$ seinen restringierten Maximalwert $\lambda_r$ an.

$\hfill\Box$

# Modellformulierung
\textcolor{darkblue}{Simulationsbeispiel}

\footnotesize
Wir betrachten das Beispiel (vgl. @uurtio_2018a)
\begin{equation}
p(X) = N(x;0_4,I_4) \mbox{ und } p(Y|X) = N(y; LX, G)
\end{equation}
mit
\begin{equation}
L := \begin{pmatrix*}[r] 0.0 & 0.0 & 1.0 & 0.0 \\ 1.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & -1.0 \end{pmatrix*}
\mbox{ und }
G := \begin{pmatrix*}[r] 0.2 & 0.0 & 0.0 \\ 0.0 & 0.4 & 0.0 \\ 0.0 & 0.0 & 0.3 \end{pmatrix*}
\end{equation}
Hier gilt offenbar $m_x = 4, m_y = 3, m = 7$ und
\begin{align}
\begin{split}
Y_1 & = \,\,\,\, X_3 + \varepsilon_1 \\
Y_2 & = \,\,\,\, X_1 + \varepsilon_2 \\
Y_3 & =       -X_4   + \varepsilon_3 \\
\end{split}
\end{align}

mit
\begin{equation}
X_1 \sim N(0,1), X_3 \sim N(0,1), X_4 \sim N(0,1)
\end{equation}
und
\begin{equation}
\varepsilon_1 \sim N(0,0.2), \varepsilon_2 \sim N(0,0.4), \varepsilon_3 \sim N(0,0.3)
\end{equation}

# Modellformulierung
\textcolor{darkblue}{Simulationsbeispiel}

\footnotesize
Mit dem Theorem zu gemeinsamen Normalverteilungen (vgl. Einheit (3) Matrizen)
ergibt sich, dass
\begin{equation}
\begin{pmatrix} X \\ Y \end{pmatrix}
\sim N(0_7,\Sigma)
\end{equation}
mit
\begin{equation}
\Sigma
=
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy}
\end{pmatrix},
\end{equation}
wobei
\begin{equation}
\Sigma_{xx} = I_4, \quad
\Sigma_{xy} = L^T, \quad
\Sigma_{yx} = L \mbox{ und }
\Sigma_{yy} = G + LL^T.
\end{equation}
Explizit ergibt sich also
\begin{equation}
\Sigma
 =
\begin{pmatrix}
I_4 &  L^T \\
L   &  G + LL^T
\end{pmatrix}
=
\begin{pmatrix*}[r]
1.0 & 0.0 & 0.0 &  0.0 & 0.0 & 1.0 &  0.0 \\
0.0 & 1.0 & 0.0 &  0.0 & 0.0 & 0.0 &  0.0 \\
0.0 & 0.0 & 1.0 &  0.0 & 1.0 & 0.0 &  0.0 \\
0.0 & 0.0 & 0.0 &  1.0 & 0.0 & 0.0 & -1.0 \\
0.0 & 0.0 & 1.0 &  0.0 & 1.2 & 0.0 &  0.0 \\
1.0 & 0.0 & 0.0 &  0.0 & 0.0 & 1.4 &  0.0 \\
0.0 & 0.0 & 0.0 & -1.0 & 0.0 & 0.0 &  1.3 \\
\end{pmatrix*}
\end{equation}

# Modellformulierung
\textcolor{darkblue}{Simulationsbeispiel}
\vspace{2mm}

\tiny
\setstretch{1.1}
```{r, warning = F, message=F}
# R Pakete für Matrizenrechnung
library(matlib)
library(expm)

# Modellparameter
L = matrix(c(0,0,1, 0,
             1,0,0, 0,
             0,0,0,-1),
           nrow  = 3,
           byrow = T)
G = diag(c(0.2,0.4,0.3))

# Kovarianzmatrixpartition
Sigma_xx = diag(4)
Sigma_xy = t(L)
Sigma_yx = L
Sigma_yy = G + L %*% t(L)
Sigma    = rbind(cbind(Sigma_xx, Sigma_xy), cbind(Sigma_yx, Sigma_yy))
print(Sigma)
```
# Modellformulierung
\textcolor{darkblue}{Simulationsbeispiel}
\vspace{2mm}
\tiny
\setstretch{1.2}

```{r}
# Evaluation der iten kanonischen Koeffizientenvektoren und Korrelationen
K      = sqrtm(inv(Sigma_xx)) %*% Sigma_xy %*% sqrtm(inv(Sigma_yy)) # K
ALB    = svd(K)                                                     # K = A\LambdaV
A      = ALB$u                                                      # A
Lambda = ALB$d                                                      # Lambda
B      = ALB$v                                                      # B
rho    = Lambda                                                     # \rho_i = \lambda_i^{1/2}
a      = sqrtm(inv(Sigma_xx)) %*% A                                 # a_i = \Sigma_{xx}^{-1/2}\alpha_i
b      = sqrtm(inv(Sigma_yy)) %*% B                                 # b_i = \Sigma_{yy}^{-1/2}\beta_i
```

\footnotesize
Die kanonische Korrelationen und kanonischen Koeffizientenvektoren  ergeben sich  zu

```{r, echo = F}
# Ausgabe
cat("rho_1 = ", rho[1],", a_1^T = (", a[,1], "), b_1^T = (", b[,1], ")",
    "\nrho_2 = ", rho[2],", a_2^T = (", a[,2],") , b_2^T = (", b[,2], ")",
    "\nrho_3 = ", rho[3],", a_3^T = (", a[,3], "), b_3^T = (", b[,3], ")")

```


# Modellschätzung
\setstretch{3}
\vfill
\large
Vorbemerkungen

Modellformulierung

**Modellschätzung**

Modellevaluation
\vfill


# Modellschätzung

\footnotesize
\begin{definition}[Schätzer kanonischer Korrelationen und Koeffizientenvektoren]
Für $i = 1,...,n$ seien
\begin{equation}
Z_i = \begin{pmatrix} X_i \\ Y_i \end{pmatrix}
\mbox{ mit }
\mathbb{E}(Z_i) := 0_m
\mbox{ und }
\mathbb{C}(Z_i) :=
\begin{pmatrix}
\Sigma_{xx} & \Sigma_{xy} \\
\Sigma_{yx} & \Sigma_{yy} \\
\end{pmatrix}
\in \mathbb{R}^{m \times m}
\end{equation}
unabhängig und identisch verteilte  $m$-dimensionale partitionierte
Zufallsvektoren sowie ihr Erwartungswert und ihre Kovarianzmatrix, respektive, und
\begin{equation}
C :=
\begin{pmatrix}
C_{xx} & C_{xy} \\
C_{yx} & C_{yy} \\
\end{pmatrix}
\in \mathbb{R}^{m \times m}
\end{equation}
sei ihre Stichprobenkovarianzmatrix. Dann sind für $i = 1,...,k := \min \{m_x,m_y\}$
\begin{equation}
\hat{a}_i := C_{xx}^{-1/2}\hat{\alpha}_i \in \mathbb{R}^{m_x}, \quad
\hat{b}_i := C_{yy}^{-1/2}\hat{\beta}_i \in \mathbb{R}^{m_y} \mbox{ und }
\hat{\rho}_i := \hat{\lambda}_i^{1/2}
\end{equation}
Schätzer der $i$ten kanonischen Koeffizientenvektoren und kanonischen Korrelationen,
respektive. Dabei sind mit
\begin{equation}
\hat{K} := C_{xx}^{-1/2}C_{xy}C_{yy}^{-1/2} \in \mathbb{R}^{m_x \times m_y}
\end{equation}
$\hat{\alpha}_i$ und $\hat{\lambda}_i$ der $i$te Eigenvektor und sein zugehöriger
Eigenwert von $\hat{K}\hat{K}^T$ und $\hat{\beta}_i$ der entsprechende Eigenkvektor
von $\hat{K}^T\hat{K}$.
\end{definition}

Bemerkungen

* Zur Modellschätzung wird  $\mathbb{C}(Z)$ also durch $C$ ersetzt.

# Modellschätzung
\textcolor{darkblue}{Simulationsbeispiel}
\vspace{2mm}

\tiny
\setstretch{1.2}
```{r}
# R Pakete
library(MASS)
library(matlib)
library(expm)

# Modellparameter
m_x      = 4
m_y      = 3
k        = min(m_x,m_y)
L        = matrix(c(0,0,1,0,1,0,0,0,0,0,0,-1), nrow = 3,byrow = 3)
G        = diag(c(0.2,0.4,0.3))
Sigma_xx = diag(4)
Sigma_xy = t(L)
Sigma_yx = L
Sigma_yy = G + L %*% t(L)
Sigma    = rbind(cbind(Sigma_xx, Sigma_xy), cbind(Sigma_yx, Sigma_yy))
K        = sqrtm(inv(Sigma_xx)) %*% Sigma_xy %*% sqrtm(inv(Sigma_yy))
ALB      = svd(K)
A        = ALB$u
Lambda   = ALB$d
B        = ALB$v
rho      = Lambda
a        = sqrtm(inv(Sigma_xx)) %*% A
b        = sqrtm(inv(Sigma_yy)) %*% B
```

# Modellschätzung
\textcolor{darkblue}{Simulationsbeispiel}

\vspace{2mm}
\tiny
\setstretch{1.2}

```{r}
# Simulationen
n       = 1e1:1e3
rho_hat = matrix(rep(NaN, length(n)*k)  , nrow = k)
a_1_hat = matrix(rep(NaN, length(n)*m_x), nrow = m_x)
for(i in 1:length(n)){

    # Datengeneration
    Y          = t(mvrnorm(n[i],rep(0, m_x+m_y),Sigma))
    I_n        = diag(n[i])
    J_n        = matrix(rep(1,n[i]^2), nrow = n[i])

    # Stichprobenkovarianzmatrixpartition
    C          = (1/(n[i]-1))*(Y %*% (I_n-(1/n[i])*J_n) %*% t(Y))
    C_xx       = C[1:m_x,1:m_x]
    C_xy       = C[1:m_x,(m_x+1):(m_x+m_y)]
    C_yx       = C[(m_x+1):(m_x+m_y),1:m_x]
    C_yy       = C[(m_x+1):(m_x+m_y),(m_x+1):(m_x+m_y)]

    # Kanonische Korrelationsanalyse
    K_hat        = sqrtm(inv(C_xx)) %*% C_xy %*% sqrtm(inv(C_yy))
    ALB_hat      = svd(K_hat)
    A_hat        = ALB_hat$u
    Lambda_hat   = ALB_hat$d
    B_hat        = ALB_hat$v
    a_hat        = sqrtm(inv(C_xx)) %*% A_hat
    b_hat        = sqrtm(inv(C_yy)) %*% B_hat
    rho_hat[,i]  = as.matrix(Lambda_hat)
    a_1_hat[,i]  = a_hat[,1]
}
```


# Modellschätzung

\textcolor{darkblue}{Simulationsbeispiel}
\vspace{3mm}

```{r, echo = F, eval = F}
graphics.off()
fdir        =  file.path(getwd(), "12_Abbildungen")
dev.new()
par(
family      = "sans",
mfcol       = c(1,1),
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)
cols        = c("gray20", "gray50", "gray80")
matplot(
n,
t(rho_hat),
type = "l",
col  = cols,
lty  = 1,
ylim = c(0.7,1),
ylab = "",
main = "Schätzung kanonischer Korrelationen")
for(i in 1:k){
  abline(h = rho[i], col = cols[i])
}
legend(
"bottomright",
c(TeX("$\\rho_1,\\hat{\\rho}_1$"),
  TeX("$\\rho_2,\\hat{\\rho}_2$"),
  TeX("$\\rho_3,\\hat{\\rho}_3$")),
lty = 1,
col = cols,
bty = "n")
dev.copy2pdf(
file        = file.path(fdir, "mvda_12_modellschätzung_1.pdf"),
width       = 7,
height      = 5)
```

```{r, echo = FALSE, out.width = "85%"}
knitr::include_graphics("12_Abbildungen/mvda_12_modellschätzung_1.pdf")
```


# Modellschätzung

\textcolor{darkblue}{Simulationsbeispiel}
\vspace{3mm}

```{r, echo = F, eval = F}
graphics.off()
fdir        =  file.path(getwd(), "12_Abbildungen")
dev.new()
par(
family      = "sans",
mfcol       = c(1,1),
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)
cols        = c("gray20", "gray50", "darkblue","gray80")
matplot(
n,
t(abs(a_1_hat)),
type = "l",
col  = cols,
lty  = 1,
ylim = c(-.1,1.2),
ylab = "",
main = "Schätzung der Absolutwerte des 1. kanonischen Koeffizientenvektors")
for(i in 1:m_x){
  abline(h = abs(a[i,1]), col = cols[i])
}
legend(
"right",
c(TeX("$|a_{11}|, |\\hat{a}_{11}|$"),
  TeX("$|a_{12}|, |\\hat{a}_{12}|$"),
  TeX("$|a_{13}|, |\\hat{a}_{13}|$"),
  TeX("$|a_{14}|, |\\hat{a}_{14}|$")),
lty = 1,
col = cols,
bty = "n")
dev.copy2pdf(
file        = file.path(fdir, "mvda_12_modellschätzung_2.pdf"),
width       = 7,
height      = 5)
```

```{r, echo = FALSE, out.width = "85%"}
knitr::include_graphics("12_Abbildungen/mvda_12_modellschätzung_2.pdf")
```

# Modellschätzung
\textcolor{darkblue}{Anwendungsbeispiel}
\small

Wir betrachten erneut den Datensatz nach @rudolf_2020 Kapitel 4

Wir betrachten die psychodiagnostischen Daten der $n = 45$ Studierenden

* $X_1$ Intelligenztestscore
* $X_2$ Mathematiktestscore
* $Y_1$ Gewissenhaftigkeitscore
* $Y_2$ Verträglichkeitscore

als 45 unabhängige Realisierungen eines $4$-dimensionaler Zufallsvektors $Z$.
Wir sind also an den kanonischen Korrelationen der Kognitionstestwerte
(Intelligenz, Mathematik) mit den Charaktertestwerten (Gewissenhaftigkeit, Veträglichkeit) interessiert.


# Modellschätzung
\textcolor{darkblue}{Anwendungsbeispiel}
\small

Daten der ersten 12 Proband:innen

\setstretch{1.2}
```{r, echo = F, message = F}
library("foreign")
D           = read.spss(file.path(getwd(), "12_Daten", "studienerfolg.sav"), to.data.frame = T)
Z           = cbind(D$X1,D$X2,D$X3,D$X4)
colnames(Z) = c("x1","x2", "y1", "y2")
knitr::kable(Z[1:12,], "pipe", align = "c")
```

# Modellschätzung
\textcolor{darkblue}{Anwendungsbeispiel}
\small

Kanonische Korrelationsanalyse
\vspace{1mm}

\tiny
\setstretch{.9}
```{r, message = F}
# Datenpräprozessierung
D          = read.spss(file.path(getwd(), "12_Daten", "studienerfolg.sav"), to.data.frame = T)
x          = as.matrix(cbind(D$X1, D$X2))
y          = as.matrix(cbind(D$X3, D$X4))
n          = nrow(x)
m_x        = ncol(x)
m_y        = ncol(y)
Y          = rbind(t(x),t(y))

# Stichprobenkovarianzmatrixpartition
I_n        = diag(n)
J_n        = matrix(rep(1,n^2), nrow = n)
C          = (1/(n-1))*(Y %*% (I_n-(1/n)*J_n) %*% t(Y))
C_xx       = C[1:m_x,1:m_x]
C_xy       = C[1:m_x,(m_x+1):(m_x+m_y)]
C_yx       = C[(m_x+1):(m_x+m_y),1:m_x]
C_yy       = C[(m_x+1):(m_x+m_y),(m_x+1):(m_x+m_y)]

# Kanonische Korrelationsanalyse
K_hat      = sqrtm(inv(C_xx)) %*% C_xy %*% sqrtm(inv(C_yy))
ALB_hat    = svd(K_hat)
A_hat      = ALB_hat$u
Lambda_hat = ALB_hat$d
B_hat      = ALB_hat$v
a_hat      = sqrtm(inv(C_xx)) %*% A_hat
b_hat      = sqrtm(inv(C_yy)) %*% B_hat
rho_hat    = as.matrix(Lambda_hat)
```

```{r, echo = F}
# Ausgabe
cat("rho_hat_1 : ", rho_hat[1],
    "\na_hat_1   : ", a_hat[,1],
    "\nb_hat_1   : ", b_hat[,1],
    "\nrho_hat_2 : ", rho_hat[2],
    "\na_hat_2   : ", a_hat[,2],
    "\nb_hat_2   : ", b_hat[,2])
```

# Modellschätzung
\textcolor{darkblue}{Anwendungsbeispiel}
\small

Kanonische Korrelationsanalyse mir R's `cancor()` Funktion

\vspace{1mm}

\tiny
\setstretch{1}

```{r, message = F}
# Datenpräprozessierung
D          = read.spss(file.path(getwd(), "12_Daten", "studienerfolg.sav"), to.data.frame = T)
x          = as.matrix(cbind(D$X1, D$X2))
y          = as.matrix(cbind(D$X3, D$X4))
cca        = cancor(x,y)

```

```{r, echo = F}
print(cca)
```
# Modellschätzung
\textcolor{darkblue}{Anwendungsbeispiel}

\vspace{2mm}

\small
Die geschätzte maximale Korrelation von Linearkombinationen von $(x_1,x_2)$ und $(y_1,y_2)$ ist 0.25.

* $(x_1,x_2)$ und $(y_1,y_2)$ sind multivariat also "gering bis mäßig" korreliert.

Basierend auf der simulationsvalidierten Schätzung ergibt sich

* $\xi       = 0.02X_1 + 0.05X_2$ als "bester Prädiktor"
* $\upsilon  = 0.09Y_1 - 0.01Y_2$ als "am besten prädizierbares Kriterium"

"Mathematikfähigkeiten" scheinen zur Prädiktion der betrachteten "Charaktereigenschaften"
etwas wichtiger als "Intelligenz"; bei den betrachteten "Charaktereigenschaften"
trägt "Gewissenhaftigkeit" mehr zum Kriterium bei als "Verträglichkeit".

# Modellschätzung
\setstretch{3}
\vfill
\large
Vorbemerkungen

Modellformulierung

Modellschätzung

**Modellevaluation**


# Modellevaluation
\textcolor{darkblue}{Überblick}

\footnotesize
Ein Ziel der Modellevaluation bei der Kanonischen Korrelationsanalyse kann das Testen von 
\begin{equation}
H_0 : \Sigma_{xy} = 0_{m_xm_y}
\end{equation}
sein. Diese Nullhypothese besagt, dass zwischen keiner der Variablen $X_1,...,X_{m_x}$
und $Y_1,...,Y_{m_y}$ eine lineare Abhängigkeit besteht. Dies impliziert, dass alle
kanonischen Korrelationen gleich 0 sind, denn es gilt
\begin{equation}
\Sigma_{xy} = 0_{m_xm_y} \Rightarrow K = 0_{m_xm_y} \Rightarrow \Lambda = 0_{m_ym_y}.
\end{equation}
Die Alternativhypothese zu dieser Nullhypothese lautet 
\begin{equation}
H_1 : \Sigma_{{xy}_{i,j}} \neq 0 
\mbox{ für mindestens ein Paar } (i,j)
\mbox{ mit } 1 \le i \le m_x \mbox{ und } 1 \le j \le m_y.
\end{equation}
Die Alternativhypothese besagt also, dass mindestens eine der Variablen $X_1,...,X_{m_x}$
und eine der Variablen $Y_1,...,Y_{m_y}$ linear abhängig sind und damit nicht alle
kanonischen Korrelationen gleich null sind.

Wie bei der einfaktoriellen Varianzanalyse (und generell im Kontext multivariater
allgemeiner lineare Modelle) können kritische Wert-basierte Tests der obigen Nullhypothese
mit verschiedenen Teststatistiken (\textit{Wilks' $\Lambda$, Pillai Statisik}) konstruiert werden,
der Verteilungen wiederum nur in Speziallfällen analytisch beschrieben sind und 
ansonsten mit $f$-Verteilungen approximiert werden.

Wir betrachten hier lediglich das Testen der obigen Nullhypothese mit der 
Wilk's $\Lambda$ Statistik.

# Modellevaluation 
\footnotesize
\begin{theorem}[Wilks' $\Lambda$ für die kanonische Korrelationsanalyse]
\justifying
\normalfont
Es seien das Modell der kanonischen Korrelationsanalyse, die Partition
der Stichprobenkovarianzmatrix, und die Schätzer der kanonischen Korrelationen
definiert wie oben. Dann hat die die Wilks' $\Lambda$ Teststatistik die Form
\begin{equation}
\Lambda = \frac{|C|}{|C_{xx}||C_{yy}|} = \prod_{i=1}^k \left(1 - \hat{\rho}_i^2\right),
\end{equation}
wobei $|\cdot|$ die Determinante bezeichnet. Weiterhin ist für
\begin{equation}
H_0 : \Sigma_{xy} = 0_{m_xm_y}
\end{equation}
die Statistik
\begin{equation}
\tau := \frac{1 - \Lambda^{1/t}}{\Lambda^{1/t}} \frac{\nu_2}{\nu_1}
\end{equation}
mit
\begin{equation}
\nu_1 := m_xm_y
\mbox{ und }
\nu_2 := wt-\frac{1}{2}m_xm_y + 1
\end{equation}
sowie
\begin{equation}
w := n - \frac{1}{2}(m_x + m_y + 3)
\mbox{ und }
t   := \sqrt{\frac{m_x^2m_y^2 - 4}{m_x^2 + m_y^2 - 5}}
\end{equation}
approximativ $f$-verteilt mit den Freiheitsgradparametern $\nu_1$ und $\nu_2$.
\end{theorem}

# Modellevaluation 
Bemerkungen
\setstretch{2.4}
\small

* Wir verzichten auf Beweise aller Aussagen dieses Theorems.
* $\Lambda$ wird klein, wenn die absoluten Werte der Einträge in $C_{xy}$ groß werden. 
* Man denke in diesem Zusammenhang an das Berechnen  der Determinante einer $2 \times 2$ Matrix.
* $\Lambda$ wird klein, wenn zumindest ein $\hat{\rho}_i$ groß ist
* Für $\hat{\rho}_i = 0$ für alle $i = 1,...,k$ gilt $\Lambda = 1$; für $\hat{\rho}_i = 1$ für alle $i = 1,...,k$ gilt $\Lambda = 0$
* Kleine Werte von $\Lambda$ und damit große Werte von $\tau$ sprechen also gegen $H_0$. 

# Modellevaluation 
\textcolor{darkblue}{Simulationsbeispiel}

\vspace{2mm}
\tiny
\setstretch{.9}

```{r, eval = F}
# Modellparameter
library(MASS)
m_x      = 2
m_y      = 2
Sigma_xx = diag(m_x)
Sigma_xy = matrix(rep(0,m_x*m_y), nrow = 2)              # H_0 : \Sigma_{xy} = 0_{m_xm_y}                                      
Sigma_yx = Sigma_xy
Sigma_yy = diag(m_y)
Sigma    = rbind(cbind(Sigma_xx, Sigma_xy), 
                 cbind(Sigma_yx, Sigma_yy))
# Testszenarioparameter
n        = 45                                            # Anzahl Datenpunkte
w        = n-((1/2)*(m_x+m_y+3))                         
t        = sqrt((m_x^2*m_y^2 - 4)/((m_x^2+m_y^2)-5))   
nu_1     = m_x*m_y                                 
nu_2     = w*t-((1/2)*m_x*m_y) + 1                  
alpha_0  = 0.05                                          # Signifikanzlevel
kW       = qf(1-alpha_0,nu_1,nu_2)                       # kritischer Wert 

# Simulationen
sim      = 1e5                                           # Anzahl an Simulationen
tau      = rep(NaN, sim)                                 # Teststatistikarray
phi      = rep(0,   sim)                                 # Testarray
for(s in 1:sim){
    Y      = t(mvrnorm(n,rep(0, m_x+m_y),Sigma))         # Datengeneration
    I_n    = diag(n)
    J_n    = matrix(rep(1,n^2), nrow = n)
    C      = (1/(n-1))*(Y %*% (I_n-(1/n)*J_n) %*% t(Y))  # Stichprobenkovarianzmatrix
    C_xx   = C[1:m_x,1:m_x]
    C_xy   = C[1:m_x,(m_x+1):(m_x+m_y)]
    C_yy   = C[(m_x+1):(m_x+m_y),(m_x+1):(m_x+m_y)]
    Lambda = det(C)/(det(C_xx)*det(C_yy))                # Wilk's Lambda           
    tau[s] = ((1-Lambda^(1/t))/Lambda^(1/t))*(nu_2/nu_1) # \tau
    phi[s] = tau[s] > kW                                 # Test
}
```

# Modellevaluation 
\textcolor{darkblue}{Simulationsbeispiel}
\vspace{2mm}
\center
\small
$n = 45, m_x = 2, m_y = 2, \nu_1 = 4, \nu_2 = 82, k = 2.48, \hat{\alpha} = 0.0498$

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

# Statistik Histogramm und WDF
f_min = 0
f_max = 3
f_res = 1e3
f     = seq(f_min,f_max,len =f_res)
p     = df(f,nu_1,nu_2)
hist(
tau,
breaks      = 100,
col         = "gray90",
prob        = TRUE,
xlim        = c(f_min, f_max),
ylim        = c(0,1),
xlab        = TeX("$\\tau$"),
ylab        = "",
main        = "")
lines(
f,
p,
lwd         = 2,
col         = "darkorange")

fdir        =  file.path(getwd(), "12_Abbildungen")
dev.copy2pdf(
file        = file.path(fdir, "mvda_12_approximative_lambda_verteilungen.pdf"),
width       = 6,
height      = 4)
```


```{r, echo = FALSE, out.width = "70%"}
knitr::include_graphics("12_Abbildungen/mvda_12_approximative_lambda_verteilungen.pdf")
```

# Modellevaluation
\textcolor{darkblue}{Simulationsbeispiel}
\vspace{1mm}

\tiny
\setstretch{.8}
```{r, message = F}
# Datenpräprozessierung
library(foreign)
D          = read.spss(file.path(getwd(), "12_Daten", "studienerfolg.sav"), to.data.frame = T)
x          = as.matrix(cbind(D$X1, D$X2))
y          = as.matrix(cbind(D$X3, D$X4))
n          = nrow(x)
m_x        = ncol(x)
m_y        = ncol(y)
Y          = rbind(t(x),t(y))

# Testszenarioparameter
w        = n-((1/2)*(m_x+m_y+3))                         
t        = sqrt((m_x^2*m_y^2 - 4)/((m_x^2+m_y^2)-5))   
nu_1     = m_x*m_y                                 
nu_2     = w*t-((1/2)*m_x*m_y) + 1                  
alpha_0  = 0.05                                           
k        = qf(1-alpha_0,nu_1,nu_2)                        

# Stichprobenkovarianzmatrixpartition
I_n        = diag(n)
J_n        = matrix(rep(1,n^2), nrow = n)
C          = (1/(n-1))*(Y %*% (I_n-(1/n)*J_n) %*% t(Y))
C_xx       = C[1:m_x,1:m_x]
C_xy       = C[1:m_x,(m_x+1):(m_x+m_y)]
C_yx       = C[(m_x+1):(m_x+m_y),1:m_x]
C_yy       = C[(m_x+1):(m_x+m_y),(m_x+1):(m_x+m_y)]

# Nullhypothesentest
Lambda     = det(C)/(det(C_xx)*det(C_yy))                        
tau        = ((1-Lambda^(1/t))/Lambda^(1/t))*(nu_2/nu_1)  
phi        = tau > k                                      
pval       = 1 - pf(tau,nu_1,nu_2)
```

```{r, echo = F}
# Ausgabe
cat("Lambda :" , Lambda,
    "\ntau    :" , tau,
    "\nk      :" , k,
    "\nphi    :" , as.numeric(phi),
    "\np      :" , pval)
```

\small
\center
$\Rightarrow$ $H_0 : \Sigma_{xy} = 0_{m_xm_y}$ wird nicht verworfen. 

# References
\footnotesize