--- format: pdf: fontsize: 11pt number-sections: true include-in-header: ../MV_Header_Seminar.tex lang: de --- \setcounter{section}{5} # Multivariate Varianzanalyse ## Mathematische Grundlagen Gegeben sei das *Modell der einfaktoriellen multivariaten Varianzanalyse* \begin{equation} \label{eq:MANOVA} y_{ij} = \mu_i + \varepsilon_{ij} \quad \text{mit} \quad \varepsilon_{ij} \sim N(0_m, \Sigma) \quad \text{u.i.v.} \; , \end{equation} wobei $i = 1,\ldots,p$ ein Index über Gruppen, $j = 1,\ldots,n_i$ ein Index über experimentelle Einheiten innerhalb einer Gruppe und $m$ die Dimensionalität eines Datenvektors $y_{ij} \in \mathbb{R}^m$ ist. Wir betrachten die *multivariate Quadratsummenzerlegung* für dieses Modell mit \vspace{-0.5em} \begin{equation} \label{eq:TBW} \begin{split} T &= \sum_{i=1}^p \sum_{j=1}^{n_i} \left(y_{ij}-\bar{y}\right) \left(y_{ij}-\bar{y}\right)^\mathrm{T} \\ B &= \sum_{i=1}^p n_i \left(\bar{y}_{i}-\bar{y}\right) \left(\bar{y}_{i}-\bar{y}\right)^\mathrm{T} \\ W &= \sum_{i=1}^p \sum_{j=1}^{n_i} \left(y_{ij}-\bar{y}_i\right) \left(y_{ij}-\bar{y}_i\right)^\mathrm{T} \end{split} \end{equation} und die darauf basierende *Wilks'-Lambda-Statistik* \begin{equation} \label{eq:Lambda} \Lambda = \frac{|W|}{|B+W|} \; , \end{equation} wobei $B$ die *between sum-of-squares Matrix*, $W$ die *within sum-of-squares Matrix* und $T = B + W$ die *total sum of-squares Matrix* ist. Das Theorem zu speziellen Verteilungen von Wilks'-$\Lambda$-Transformationen gibt Spezialfälle von Datendimension $m$ und Gruppenanzahl $p$ an, in denen die Verteilung einer Transformation von $\Lambda$ unter der Nullhypothese $H_0 : \mu_1 = \ldots = \mu_p$ exakt bekannt ist. Es besagt, dass für den in den ersten beiden Tabellenspalten aufgeführten Spezialfall die in der dritten Tabellenspalte genannte Teststatistik bei Vorliegen der Nullhypothese einer $f$-Verteilung mit den Freiheitsgradparametern in der vierten Tabellenspalte folgt: \begin{center} \renewcommand{\arraystretch}{1.2} \begin{tabular}{cccc} Datendimension $m$ & Gruppenanzahl $p$ & Statistik & $f$-Verteilungsparameter \\ \hline beliebig & 2 & $\frac{1 - \Lambda}{\Lambda} \frac{n-p-m+1}{m}$ & $m, \; n-p-m+1$ \\ beliebig & 3 & $\frac{1 - \sqrt{\Lambda}}{\sqrt{\Lambda}} \frac{n-p-m+1}{m}$ & $2m, \; 2(n-p-m+1)$ \\ 1 & beliebig & $\frac{1 - \Lambda}{\Lambda} \frac{n-p}{p-1}$ & $p-1, \; n-p$ \\ 2 & beliebig & $\frac{1 - \sqrt{\Lambda}}{\sqrt{\Lambda}} \frac{n-p-1}{p-1}$ & $2(p-1), \; 2(n-p-1)$ \end{tabular} \end{center} Damit gilt z.B. für $m>2$ und $p=3$ (zweite Tabellenzeile): Überschreitet $\frac{1 - \sqrt{\Lambda}}{\sqrt{\Lambda}} \frac{n-p-m+1}{m}$ den kritischen Wert $k_{\alpha_0}$, der mit einem Signifikanzniveau $\alpha_0$ aus der inversen kumulativen Verteilungsfunktion der $f$-Verteilung mit Freiheitsgradparametern $2m$ und $2(n-p-m+1)$ errechnet werden kann, wird die Nullhypothese identischer Erwartungswerte über Gruppen abgelehnt. Anderenfalls wird sie nicht abgelehnt. ## Analyse in R Die Datei `FADE_SAME.csv` enthält den in der ersten Seminarsitzung vorgestellten Datensatz. Erklären Sie die Funktion des folgenden R-Codes: \footnotesize ```{r, echo = T, eval = T, results = 'hide'} # Daten einlesen fname = 'FADE_SAME.csv' # Dateiname D = read.csv(fname) # Dataframe # Datenmatrix extrahieren rows = startsWith(D$subject, 'subA') # Personen aus Studie A cols = c('novelty.FADE','novelty.SAME', # Definition Variablen 'memory.FADE', 'memory.SAME') n = sum(rows) # Anzahl Datenpunkte m = length(cols) # Anzahl Variablen Y = t(as.matrix(D[rows,cols])) # m x n Datenmatrix print(dim(Y)) # Überprüfung Datenmatrix ``` \normalsize ## Erste Programmieraufgabe Wir wollen systematische Unterschiede der FADE-Scores und SAME-Scores zwischen den drei Altersgruppen "junge Personen" (bis 35 Jahre), "mittelalte Personen" (50-59 Jahre) und "ältere Personen" (ab 60 Jahre) untersuchen. Wir werden dazu zunächst eine Variable $x$ definieren, die die unabhängige Variable "Altersgruppe" repräsentiert und die Zugehörigkeit der Datenpunkte zu den verschiedenen Studiengruppen kodiert. Gehen Sie dazu wie folgt vor: - Extrahieren Sie die Variable "age" aus dem Dataframe, aber nur für diejenigen Datenpunkte, die durch die Variable \verb|rows| beschrieben werden. - Bilden Sie drei logische Variablen \verb|g1|, \verb|g2|, \verb|g3|, d.h. Vektoren von derselben Länge wie die Datenmatrix, die jeweils \verb|TRUE| sind, wenn der entsprechende Datenpunkt zu Gruppe 1, 2, 3 gehört. Beispielsweise ist \verb|g1 = (x <= 35)|. - Bilden Sie einen Vektor, der für jeden Datenpunkt den entsprechenden Gruppenindex enthält, d.h. einen Vektor, der 1 ist, wenn Altersgruppe \verb|g1| vorliegt, 2 ist, wenn Altersgruppe \verb|g2| vorliegt und 3 ist, wenn Altersgruppe \verb|g1| vorliegt. Sie erreichen dies z.B., indem Sie die logischen Vektoren mit den Gruppenindizes multiplizieren und aufsummieren. Speichern Sie das Ergebnis in die Variable \verb|x|. - Speichern Sie das Maximum von \verb|x|, d.h. die Anzahl der Gruppen, in die Variable \verb|p|. - Geben Sie die ersten 100 Einträge des Vektors \verb|x| aus. Sie sollten folgende Ergebnisse erhalten: \footnotesize ```{r, echo = F, eval = T} ``` \normalsize ## Analyse in R Im zweiten Teil wollen wir nun mit der Datenmatrix $Y$ und dem Gruppenvektor $x$ eine einfaktorielle multivariate Varianzanalyse durchführen. Erklären Sie dazu zunächst den folgenden R-Code und die Funktionen, die er definiert: \footnotesize ```{r, echo = T, eval = T, results = 'hide'} # Parameter schätzen estimate = function(Y,x) { # Datendimensionalität n = ncol(Y) # Anzahl Datenpunkte m = nrow(Y) # Anzahl Variablen p = max(x) # Anzahl Gruppen # Parameterschätzer mu_hat = matrix(rep(0,m*p), nrow = m) # Erwartungswertparameter Sigma_hat = matrix(rep(0,m*m), nrow = m) # Kovarianzmatrixparameter for (i in 1:p) { Y_i = Y[,x==i] # Datenmatrix Gruppe i n_i = ncol(Y_i) # Anzahl Datenpunkte Gruppe i mu_hat[,i] = as.matrix(rowMeans(Y_i)) # Erwartungswerte Gruppe i for (j in 1:n_i) { Sigma_hat = Sigma_hat + (Y_i[,j] - mu_hat[i]) %*% t(Y_i[,j] - mu_hat[i]) } } Sigma_hat = (1/(n-p)) * Sigma_hat # Funktionswerte return(list(mu_hat = mu_hat, Sigma_hat = Sigma_hat)) } # Multivariate Quadratsummenzerlegung sumofsqr = function(Y,x) { # Datendimensionalität n = ncol(Y) # Anzahl Datenpunkte m = nrow(Y) # Anzahl Variablen p = max(x) # Anzahl Gruppen # sum-of-squares Matrizen y_bar = as.matrix(rowMeans(Y)) # Gesamtstichprobenmittel y_i_bar = matrix(rep(0,m*p), nrow = m) # Gruppenstrichprobenmittel B = matrix(rep(0,m*m), nrow = m) # between-group sum-of-squares Matrix W = matrix(rep(0,m*m), nrow = m) # within-group sum-of-squares Matrix for (i in 1:p) { Y_i = Y[,x==i] # Datenmatrix Gruppe i n_i = ncol(Y_i) # Anzahl Datenpunkte Gruppe i y_i_bar[,i] = as.matrix(rowMeans(Y_i)) # Erwartungswerte Gruppe i for (j in 1:n_i) { B = B + (y_i_bar[,i] - y_bar) %*% t(y_i_bar[,i] - y_bar) W = W + (Y_i[,j] - y_i_bar[,i]) %*% t(Y_i[,j] - y_i_bar[,i]) } } T = B + W # total sum-of-squares Matrix # Funktionswerte return(list(T = T, B = B, W = W)) } ``` \normalsize ## Zweite Programmieraufgabe Nutzen Sie die Funktionen der Parameterschätzung und der multivariaten Quadratsummenzerlegung, um die MANOVA-Modellparameter zu schätzen. Gehen Sie dazu wie folgt vor: - Wenden Sie die Funktion \verb|estimate| auf $Y$ und $x$ an und extrahieren Sie die Felder \verb|mu_hat| sowie \verb|Sigma_hat|. - Wenden Sie die Funktion \verb|sumofsqr| auf $Y$ und $x$ an und extrahieren Sie die Felder \verb|T| sowie \verb|B| und \verb|W|. - Geben Sie die Schätzer für die Gruppenerwartungswertparameter und den Kovarianzmatrixparameter aus. Sie sollten folgende Ergebnisse erhalten: \footnotesize ```{r, echo = F, eval = T} ``` \normalsize \pagebreak Führen Sie nun eine einfaktorielle multivariate Varianzanalyse (MANOVA) durch. Gehen Sie dazu wie folgt vor: - Legen Sie das Signifikanzniveau $\alpha_0 = 0.05$ fest. - Berechnen Sie die Wilks'-Lambda-Statistik $\Lambda$ anhand der obigen Formel. Die Determinante ist in R als die Funktion \verb|det()| implementiert. - Suchen Sie aus der Tabelle in Abschnitt 6.1 den für diese Datenanalyse passenden Spezialfall heraus, berechnen Sie aus $\Lambda$ die Teststatistik und speichern Sie das Ergebnis in die Variable \verb|eff|. - Berechnen Sie die zum Spezialfall gehörenden Freiheitsgradparameter und speichern Sie diese in die Variablen \verb|df_1| und \verb|df_2|. - Ermitteln Sie den kritischen Wert $k_{\alpha_0}$ des Tests. Die inverse kumulative Verteilungsfunktion der F-Verteilung ist in R als \verb|qf()| implementiert. Informieren Sie sich über die Eingabeparameter dieser Funktion, um den kritischen Wert korrekt zu bestimmen. - Ermitteln Sie den p-Wert für den Test. Die kumulative Verteilungsfunktion der F-Verteilung ist in R als \verb|pf()| implementiert. Informieren Sie sich über die Eingabeparameter dieser Funktion, um den p-Wert korrekt zu bestimmen. - Ermitteln Sie das Testergebnis. Es ist 1, wenn die Teststatistik den kritischen Wert überschrittet. Andernfalls ist es 0. - Geben Sie die Resultate ihrer Analyse aus. Sie sollten folgende Ergebnisse erhalten: \footnotesize ```{r, echo = F} ``` \normalsize - Überprüfen Sie Ihre multivariate Varianzanalyse, indem Sie folgenden Code mit der in R eingebauten Funktion \verb|Manova()| ausführen: \footnotesize ```{r, eval = F, warning = F, results = 'hide'} # Einfaktorielle multivariate Varianzanalyse (automatisch) library(car) # R-Paket "car" group = as.factor(x) # faktorielle Variable model = lm(t(Y) ~ group) # Modellformulierung result = Manova(model, test.statistic = "Wilks") # Modellevaluation print(result) # Ergebnisausgabe ``` \normalsize \pagebreak ## Lückentext Füllen Sie mit den in der Übung gewonnenen Erkenntnissen den folgenden Lückentext aus und präsentieren Sie die Ergebnisse im Seminar: \vspace{1em} **Lückentext:** Die Anzahl der beobachteten experimentellen Einheiten ist __________, die Anzahl der abhängigen Variablen pro Einheit ist __________ und die Anzahl der Level der unabhängigen Variable ist __________. Die Wilks'-Lambda-Statistik hat den Wert __________, die daraus sich ergebende Teststatistik ist __________. Der kritische Wert für die multivariate ANOVA ergibt sich aus __________ und beträgt __________. Die Nullhypothese identischer Gruppenerwartungswerte wird daher bei einem Signifikanznivea von 5 Prozent __________ (abgelehnt/nicht abgelehnt). Der p-Wert für die multivariate ANOVA ergibt sich aus __________ und beträgt __________. ## Mögliche Klausurfrage Präsentieren Sie im Seminar folgende Klausurfrage und erklären Sie die richtige Antwort: \vspace{1em} **Frage:** Welche Aussage über die Parameterschätzer $\hat{\mu}_i$ und $\hat{\Sigma}$ im Modell der einfaktoriellen multivariaten Varianzanalyse ist **nicht** korrekt? \begin{enumerate}[a)] \item $\hat{\mu}_i$ ist das multivariate Stichprobenmittel der $i$-ten Gruppe. \item $\hat{\mu}_i$ entspricht dem Gesamtstichprobenmittel $\bar{y}$. \item $\hat{\Sigma}$ entspricht einer skalierten Form der der within-group sum-of-squares Matrix $W$. \item $\hat{\Sigma}$ ist ein unverzerrter Schätzer des Kovarianzmatrixparameters $\Sigma$. \end{enumerate} ## Kinderwitz Gestern haben wir zwei Biber beim Essen beobachtet: \begin{turn}{180} Gab Steg. \end{turn}