32 Die repeated / mixed ANOVA
Letzte Änderung am 01. April 2025 um 19:16:37
“I never once failed at making a light bulb. I just found out 99 ways not to make one.” — Thomas A. Edison
Stand des Kapitels: Konstruktion (seit 02.2025)
Dieses Kapitel wird in den nächsten Wochen geschrieben. Ich plane zum Ende des SoSe 2025 eine neue Version des Kapitels erstellt zu haben. Während das Kapitel entsteht, funktioniert so manaches dann nicht so wie es soll.
In diesem Kapitel soll es dann um Messwiederholungen gehen. Wir haben also irgendwas wiederholt gemessen. Damit haben wir auch eine Zeitvariable mit in den Daten. Hier müssen wir dann etwas aufpassen, ob wir wirklich eine repeated oder mixed ANOVA rechnen wollen oder nicht doch irgendwie eine Zeitreihe vorliegen haben. Die Frage, die du dir stellen musst ist hierbei eigentlich recht einfach. Wenn du ein Interesse an einem wie auch immer gearteten Gruppenvergleich hast und du diese Gruppen über die Zeit gemessen hast, dann bist du hier in diesem Kapitel richtig. Wenn es um eine ANOVA geht, dann geht es auch immer um einen Gruppenvergleich. Ansonsten könntest du eine Prognose wollen und das wäre dann eher eine klassische Zeitreihe.
Hier eine kurze Warnung. Häufig wird die repeated und die mixed ANOVA in einen Topf geworfen. Die beiden Verfahren unterscheiden sich aber konzeptionell. Darüber hinaus hast du meistens in den Agrarwissenschaften eine mixed ANOVA vorliegen, wenn du mit Messwiederholungen gearbeitet hast und keine repeated ANOVA.
32.1 Allgemeiner Hintergrund
Die gemischte ANOVA (eng. mixed ANOVA) und die ANOVA für Messwiederholungen (eng. repeated ANOVA) sind eng miteinander verwandt. Vermutlich ist dies auch der Grund, warum es bei beiden Verfahren so viele Verwirrungen gibt. Ebenso finden sich auch verschiedene Namen zu den beiden ANOVA Algorithmen, aber wir bleiben einmal hier bei den beiden Namen. Die beiden Verfahren basieren auf Messwiederholungen. Daher messen wir wiederholt an einer Beobachtung einen Messwert. Der Unterschied liegt in dem Versuchsdesign und welche Beobachtung welche Behandlung erhält. Wir haben haben meistens in den Agrarwissenschaften eine mixed ANOVA vorliegen.
- Was haben die mixed ANOVA und die repeated ANOVA gemeinsam?
-
In beiden Fällen haben wir Messwiederholungen in unseren Daten. Das heißt, wir messen wiederholt an dem gleichen Subjekt einen Messwert. Wir messen also weiderholt die Größe einer Pflanze oder die Milchleistung einer Kuh.
Oder um den Sachstand nochmal mit einem Zitat zu untermauern, habe ich aus dem Handbuch der SPSS Statistik einmal folgende Einordnung gefunden. Hier ist nochmal wichtig, was eigentlich mit den einzelnen Beobachtungen (eng. subject) passiert.
“However, the fundamental difference is that in a mixed ANOVA, the subjects that undergo each condition (e.g., a control and treatment) are different, whereas in a two-way repeated measures ANOVA, the subjects undergo both conditions (e.g., they undergo the control and the treatment).” — Two-way repeated measures ANOVA using SPSS Statistics
Das klingt jetzt etwas kompliziert und deshalb habe ich dir den Zusammenhang einmal in der folgenden Abbildung 32.1 visualisiert. Nehmen wir einmal an, du hast einen Behandlungsfaktor \(f_A\) mit drei Leveln. Jedes Level steht für eine Art der Behandlung. Wenn wir uns jetzt eine mixed ANOVA anschauen, dann haben wir drei Individuen vorliegen, die wiederholt in der gleichen Art der Behandlung gemessen werden. Das Individium mit der \(ID_1\) erhält nur die Behandlung \(A.3\) wiederholt verabreicht oder steht unter dieser Behandlung. Dann messen wir wiederholt den Messwert.
Etwas anders sieht es dann bei der repeated ANOVA aus. Hier haben wir ein Individum was alle Behandlungen über alle Zeitpunkte erhält. Damit erhät das Individium \(ID_1\) erst die Behandlung \(A.3\) viermal und danach dann die beiden anderen Behandlungen auch jeweils viermal wiederholt gemessen. Das ist ein großer konzeptioneller Unterschied in dem Versuchsdesign zwischen der mixed und der repeated ANOVA. In beiden ANOVA Typen wird wiederholt gemessen. In der repeated ANOVA erhält jedes Individum aber auch alle Behandlungen über die Zeit.
Bevor wir jetzt weiter in das Modell und die einzelnen Besonderheiten der mixed ANOVA und repeated ANOVA einsteigen nochmal der Hinweis, dass wir hier natürlich auch wieder nur global einen Fakor testen. Das heißt, du weißt natürlich wieder nicht, welcher paarweise Vergleich innerhalb eines Faktors einen signifikanten Unterschied ausmacht, wenn du eine signifikante ANOVA vorliegen hast. Dann musst du wieder weiter machen mit einem Post-hoc Test. Auch hier hilft nochmal die Vingette von {afex} mit ANOVA and Post-Hoc Contrasts. Oder wie es Gueorguieva und Krystal (2004) mit der wissenschaftlichen Veröffentlichung Move Over ANOVA und der Nutzung von linearen gemischten Modellen wie folgt beschreibt.
“Mixed-effects models use all available data, can properly account for correlation between repeated measurements on the same subject, have greater flexibility to model time effects, and can handle missing data more appropriately. Their flexibility makes them the preferred choice for the analysis of repeated measures data” — Gueorguieva und Krystal (2004)
Warum machen wir das nicht immer? Also wir nutzen immer ein lineares gemischtes Modell? Das hat meistens mit der geringen Anzahl an Beobachtungen in den Gruppen und allgemein mit zu wenigen Gruppen zu tun. Häufig reicht eben das Design und die Wiederholungen nicht aus für ein sauberes gemischtes Modell. Deshab ist dann eine Lösung sich die Sachlage einmal in einer mixed oder repeated ANOVA anzuschauen. Es gilt wie immer, komplexe Fragestellungen kann man meistens nicht in einem einfachen Satz beantworten. Hier kommt es dann auch wieder auf die konkrete Fragestellung an.
Welche Pakete gibt es eigentlich?
Wenn um die Anwendung der mixed oder repeated ANOVA in R geht, dann haben wir eine Menge Pakete zur Auswahl. Wie immer macht die Fragestellung und das gewählte Modell den Großteil der Entscheidungsindung aus. Ich zeige dir später in der Anwendung dann auch alle Pakete einmal, gebe dir dann aber auch immer eine Empfehlung mit. Normalerweise brauchen wir einen normalverteilten Messwert \(y\) und Varianzhomogenität in den Faktoren. In den letzten Jahren wurden aber noch weitere Implementierungen der ANOVA entwickelt, so dass hier auch Alternativen vorliegen.
Gehen wir jetzt mal die Pakete durch. Wir immer gibt es einiges an Möglichkeiten und ich zeige dir hier eben die Auswahl. Es gibt hier das ein oder andere noch zu beachten, aber da gehe ich dann bei den jeweiligen Methoden drauf ein. Es macht eben dann doch einen Unterschied ob ich eine einfaktorielle oder komplexere ANOVA rechnen will. Nicht alles geht in allen R Pakten oder gar Excel.
- Der Standard mit der Funktion
aov()aus{base} -
Die Standardfunktion
aov()erlaubt es eine mixed oder repeated ANOVA direkt auf einem Datensatz zu rechnen. Hier brauchen wir nur ein Modell in der in R üblichen Formelschreibweisey ~ f. Du kannst diesen Ansatz als schnelle ANOVA begreifen. Dein Messwert \(y\) muss hier normalverteilt sein. - Mit der Funktion
anova_test()aus{rstatix} -
Das R Paket
{rstatix}bietet bietet die Funktionanova_test()die es erlaubt einfacher auf die Funktionalität der Funktionaov()und R Funktionen aus dem Paket{car}zuzugreifen. Es ist eben ein klassischer Wrapper um besser rechnen zu können. Aber auch hier musst du verstehen wie die Funktion funktioniert. Ich finde aber die Anwendung relativ einfach und deshalb zeige ich die Funktion hier einmal. - Mit den ANOVA Funktionen aus
{afex} -
Wir können die mixed oder repeated ANOVA anwenden, wenn wir Messwiederholungen vorliegen haben. Daher bietet sich hier insbesondere das R Paket
{afex}an. Du musst bei den Funktionen von{afex}immer eine ID mitlaufen lassen, die angibt welche Individuen wiederholt gemessen wurden. Also hat jede Zeile eine Nummer, die beschreibt welche Beobachtung hier vorliegt. Besonders wichtig bei Messungen über die Zeit. Darüber hinaus kann das Paket sehr gut Interaktionen schätzen und Bedarf dort keiner zusätzlichen Optionen. Auch hier hilft nochmal die Vingette von{afex}mit ANOVA and Post-Hoc Contrasts weiter, wenn es um den Post-hoc Test geht. - Mit den ANOVA Funktionen aus
{MANOVA.RM} -
Das R Paket
{MANOVA.RM}von Friedrich u. a. (2019) bietet eine weitere Möglichkeit die reapted und mixed ANOVA zu rechnen, wenn die Annahme an die Normalverteilung sowie der Varianzhomogenität nicht erfüllt ist. Wir haben es hier dann nochmal mit einer neueren Implementierung zu tun. Eine weitere robuste Alternative zu den Standardimplementierungen. - Mit den ANOVA Funktionen aus
{WRS2} -
Eine Annahme an die mixed oder repeated ANOVA ist, dass wir es mit normalverteilten Messwerten \(y\) sowie Varianzhomogenität in den Faktoren \(f\) vorliegen haben. Das R Paket
{WRS2}mit der hervorragenden Vingette Robust Statistical Methods Using WRS2 erlaubt nun aber diese beiden Annahmen zu umgehen und bietet eine robuste ANOVA an. Robust meint hier, dass wir uns nicht um die Normalverteilung und Varianzhomogenität kümmern müssen.
Wir immer ist es hier ein Tutorium um gleich mit der ANOVA starten zu können. Wenn du tiefer in die Materie einsteigen möchtest und noch andere Dinge rechts und links anschauen willst, dann findest du in den folgenden Quellen nochmal mehr Anregungen. Es geht wie immer mehr, aber wir wollen hier dann ja auch fertig werden und haben beschlossen eine repeated oder mixed ANOVA zu rechnen.
Weitere Tutorien für die repeated & mixed ANOVA
Wir oben schon erwähnt, kann dieses Kapitel nicht alle Themen der repeated und mixed ANOVA abarbeiten. Insbesondere der theoretische Hintergrund ist ja bei mir hier sehr gekürzt. Daher präsentiere ich hier eine Liste von Literatur und Links, die mich für dieses Kapitel hier inspiriert haben. Nicht alles habe ich genutzt, aber vielleicht ist für dich was dabei.
- Im Tutorium zur Repeated Measures ANOVA in R und Mixed ANOVA in R auf der Seite Datanovia wird nochmal mehr auf die möglichen Vortest der Annahmen eingegangen. Auch wird nochmal der Posthoc-Test gezeigt, den ich ja ausgelagert habe.
- Das Kapitel zu den Random Effects aus dem Openbook Data Analysis in R gibt nochmal einen Einblick in die Idee der zufälligen Effekte, die wir hier bewusst etwas unter den Tisch fallen lassen.
- Das Openbook ANOVA and Mixed Models A Short Introduction Using R gibt in dem Kapitel zu Random and Mixed Effects Models nochmal eine Übersicht über die verschiedenen Modelle. Insbesondere der Theorieteil ist gut, da ich ja hier in diesem Kapitel aktuell auf eine theoretische Betrachtung verzichtet habe.
- Wenn dich die Formelschreibweise in R verwirrt, dann hilft dir nochmal der Blogeintrag zu ANOVA and other models, mixed and fixed.
32.2 Genutzte R Pakete
Wir wollen folgende R Pakete in diesem Kapitel nutzen.
R Code [zeigen / verbergen]
pacman::p_load(tidyverse, magrittr, broom, WRS2, scales,
readxl, see, car, patchwork, rstatix, emmeans,
interactions, effectsize, afex, report, multcomp,
janitor, MANOVA.RM, conflicted)
conflicts_prefer(dplyr::mutate)
conflicts_prefer(dplyr::summarize)
conflicts_prefer(dplyr::filter)
conflicts_prefer(dplyr::arrange)
conflicts_prefer(dplyr::select)
conflicts_prefer(magrittr::set_names)
conflicts_prefer(effectsize::eta_squared)
cbbPalette <- c("#000000", "#E69F00", "#56B4E9", "#009E73",
"#F0E442", "#0072B2", "#D55E00", "#CC79A7")An der Seite des Kapitels findest du den Link Quellcode anzeigen, über den du Zugang zum gesamten R-Code dieses Kapitels erhältst.
32.3 Daten
Beginnen wir wie immer mit den Daten. Hier ist die Besonderheit, dass wir immer Messweiderholungen an einer Beobachtung vorliegen haben. Ich habe dann die jeweilige Beobachtung in der Spalte .id jeweils genau benannt. Daher kommt eben die Beobachtung mit der id gleich 1 mehrfach in den Daten vor. Genauer, wir haben eben nicht nur eine Zeile in den Daten für eine Beobachtung sondern mehrere Zeilen für eine Beobachtung vorliegen. Das ist am Anfang immer etwas verwirrend, aber man gewöhnt sich daran.
Eine weitere Besonderheit ist, dass wir häufig die Daten in einem Wide-Format erstellen. Das geißt wir haben pro Zeile eine ID und dann mehrere Spalten wo wir wiederholt messen. Wir nennen dann diese Spalten meistens t0 bis t5, wenn wir sechs Zeitpunkte bemessen. In R brauchen wir dann doch wieder das Long-Format und dementsprechend nutzen wir die Funktion pivot_longer() um uns aus den Daten im Wide-Format die entsprechenden Daten im Long-Format zu bauen. Das Openbook R for Data Science hat hier noch das passende Kapitel zum nachlesen. Du kannst auch das Cheatsheet nochmal anschauen, dort wird alles nochmal visuell erklärt.
32.3.1 Mixed ANOVA
Bei der mixed ANOVA haben wir mindestens ein zweifaktorielles Design vorliegen. Warum ist das so? Zum einen brauchen wir einen Behandlungsfaktor und dazu kommt dann noch der Faktor der Messwiederholung. Wir messen ja wiederholt und dieses wiederholte Messen kommt dann ebenfalls in einen Faktor rein. Somit haben wir einen Behandlungsfaktor und einen Zeitfaktor in den Daten. Manchmal erhöhen wir noch die Komplexität und ergänzen einen zweiten Behandlungsfaktor und haben dann einen dreifaktoriellen Datensatz vorliegen.
Zweifaktoriell
Im Folgenden siehst du einmal einen Datensatz für eine zweifaktorielle mixed ANOVA. In diesem Beispiel schauen wir uns nur Katzenflöhe an. Wir haben in dem Datensatz einmal den Faktor feeding als die Fütterungsarten mit einer Kontrolle als Zuckerwasser, Blut und einmal Ketchup sowie als zweiten Faktor die verschiedenen Zeitpunkte der Messung der Sprungweite in [cm]. Wir haben an insgesamt sechs Zeitpunkten wiederholt die Sprungweite gemessen. Wir wollen nun wissen, in wie weit sich die Sprungweite über die Zeit für die drei Floharten ändert.
R Code [zeigen / verbergen]
mixed_fac2_tbl <- read_excel("data/fleas_complex_data.xlsx",
sheet = "mixed-fac2") |>
select(.id, feeding, t0:t5) |>
pivot_longer(cols = t0:t5,
values_to = "jump_length",
names_to = "time_fct") |>
mutate(feeding = as_factor(feeding),
time_fct = as_factor(time_fct),
jump_length = round(jump_length, 2),
.id = as_factor(.id))In der folgenden Tabelle siehst du einmal die rohen Daten eingelesen. Wir haben den Faktor feeding sowie die sechs Zeitpunkte der Messung der Sprungweite in [cm]. Wir müssen uns dann in R die Daten über die Funktion pivot_longer() in ein Long-Format umwandeln.
| .id | feeding | t0 | t1 | t2 | t3 | t4 | t5 |
|---|---|---|---|---|---|---|---|
| 1 | ctrl | 31.56 | 27.61 | 26.12 | 17.66 | 18.33 | 7.43 |
| 2 | ctrl | 32.87 | 26.12 | 28.23 | 19.25 | 12.01 | 8.75 |
| 3 | ctrl | 33.43 | 31.23 | 26.59 | 24.19 | 13.95 | 0.27 |
| 4 | ctrl | 24.88 | 30.47 | 26.82 | 23.77 | 20.71 | 1.34 |
| ... | ... | ... | ... | ... | ... | ... | ... |
| 21 | ketchup | 48.68 | 47.28 | 33.81 | 21.89 | 27.6 | 33.84 |
| 22 | ketchup | 51.56 | 49.85 | 37.17 | 30.66 | 27.67 | 35.03 |
| 23 | ketchup | 53.51 | 47.6 | 33.31 | 30.46 | 25.32 | 28.31 |
| 24 | ketchup | 49.38 | 42.85 | 31.99 | 34.09 | 27.1 | 34.81 |
Jetzt habe ich dir in der folgenden Abbildung auch nochmal das Versuchsdesign aufgezeigt. Wie du siehst, haben wir immer einen Floh, der eine Behandlung mit einer Fütterungsart erhält. Dieser Flog springt dann zum Beispiel unter der Fütterung mit Ketchup sechs Mal und wir notieren jedesmal die Sprungweite zu den Messzeitpunkten.
feeding. Jeder Floh erhält eine Fütterung und anschließend wird sechs mal die Sprungweite gemessen. Hier sind nur vier Zeitpunkte dargestellt. [Zum Vergrößern anklicken]
In der folgenden Abbildung siehst du dann nochmal den Zusammenhang zwischen der Sprungweite in [cm] der Katzenflöhe unter unterschiedlichen Fütterungsarten. Dabei wurden die Flöhe wiederholt an sechs Zeitpunkten gemessen. Wir du siehst, fällt mit der Zeit die Sprungleistung ab. Die Katzenflöhe springen unabhängig von der Fütterungsart weniger weit mit der Zeit. Am letzten Zeitpunkt scheinen die Flöhe unter Ketchupernährung dann wieder etwas weiter zu springen. Jetzt müssen wir aber einmal statistisch Nachweisen, ob wir hier wirklich einen signifikanten Zusammenhang und Unterschiede zwischen den Fütterungsarten vorliegen haben.
Dreifaktoriell
Für unser dreifaktorielles Beispiel erweiter wir einmal den Datensatz um ein Workout für die Katzenflöhe. Neben den drei Fütterungsarten Kontrolle mit Zuckerwasser, Blut sowie Ketchup mussten die Flöhe als zweiten Faktor noch ein Workout bestehen. Dabei haben nicht alle Katzenflöhe jeweils trainiert, sondern wir haben eine untrainierte Gruppe sowie eine trainierte Gruppe. Beide Gruppen des Workout wurden dann mit allen drei Fütterungsarten gefüttert. Dann wurden insgesamt zu sechs Zeitpunkten die Sprungweite in [cm] bestimmt. Wir haben es daher mit einem dreifaktoriellen Datensatz zu tun, da wir zwei Behandlungsfaktoren sowie den Faktor der Messwiederholungen in unseren Daten vorliegen haben.
R Code [zeigen / verbergen]
mixed_fac3_tbl <- read_excel("data/fleas_complex_data.xlsx",
sheet = "mixed-fac3") |>
select(.id, feeding, workout, t0:t5) |>
pivot_longer(cols = t0:t5,
values_to = "jump_length",
names_to = "time_fct") |>
mutate(feeding = as_factor(feeding),
workout = as_factor(workout),
time_fct = as_factor(time_fct),
jump_length = round(jump_length, 2),
.id = as_factor(.id))In der folgenden Tabelle sind die Daten dann nochmal im Wide-Format dargestellt. Wir haben hier dann pro Floh jeweils eine Faktorkombination aus der Fütterungsart feeding sowie dem Traingsstatus workout. Im Anschluss haben wir dann die Sprungweite von jedem Floh sechsmal gemessen. Wir müssen dann in R die Daten nochmal in das Long-Format über die Funktion pivot_longer() umbauen.
| .id | feeding | workout | t0 | t1 | t2 | t3 | t4 | t5 |
|---|---|---|---|---|---|---|---|---|
| 1 | ctrl | yes | 44.5 | 32.4 | 33.02 | 18.14 | 20.88 | 10.78 |
| 2 | ctrl | yes | 36.33 | 37.05 | 31.79 | 25.72 | 28.01 | 8.79 |
| 3 | ctrl | yes | 31.76 | 31.75 | 29.82 | 26.39 | 27.4 | 15.52 |
| 4 | ctrl | yes | 38.28 | 37.08 | 30.53 | 27.8 | 28.37 | 5.21 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 21 | ketchup | no | 59.12 | 51.66 | 35.64 | 39.56 | 27.77 | 32.88 |
| 22 | ketchup | no | 53.51 | 58.28 | 36.65 | 34.62 | 28.09 | 38.01 |
| 23 | ketchup | no | 62.22 | 48.53 | 46.73 | 34.81 | 34.2 | 29.03 |
| 24 | ketchup | no | 53.4 | 52 | 37.78 | 47.77 | 35.5 | 27.63 |
Manchmal ist es schwer die Zusammenhänge dann in dem Datensatz zu sehen, deshalb hier nochmal das Versuchsdesign visualisert. Jeder Katznfloh springt wiederholt, aber das ganz dann nur unter einer Faktorkombination.
feeding. Der zweite Faktor der Status des Workouts für jeden Floh workout. Jeder Floh erhält eine Fütterung sowie ein Workout und anschließend wird sechs mal die Sprungweite gemessen. Hier sind nur vier Zeitpunkte dargestellt. [Zum Vergrößern anklicken]
Abschließend können wir uns dann nochmal die Daten in der folgenden Abbildung anschauen. Auf den ersten Blick sehen wir, dass die Leistung in der Sprungweite über die Zeit abfällt. Diesen Trend können wir über alle Fütterungsarten beobachten. Was erst auf den zweiten Blick auffällt ist, dass teilweise ja nach Fütterungsart ein Training positiv oder negativ ausfällt. So führt ein Training bei den Kontrollflöhen zu weiteren Sprüngen. Die Katzenflöhe unter einer Blutdiät springen jedoch unter Training nicht so weit wie ohne ein Training.
32.3.2 Repeated ANOVA
Die repeated ANOVA ist in den Agrarwissenschaften selten. In den Sozialwissenschaften oder aber allgemeiner den Humanwissenschaften kommt die repeated ANOVA häufiger vor. Das hat vor allem damit zu tun, dass wir jedes Individuum mit allen Faktorkombinationen behandeln und wiederholt messen. Das wir ein Subjekt wiederholt messen kommt auch häufiger in den Agrarwissenschaften vor, aber selten erhält jedes Tier oder Pflanze alle Behandlungen verabreicht. Deshalb hier einmal die Daten für die repated ANOVA als Veranschaulichung was du als Daten und Experiment vorliegen haben musst.
Einfaktoriell
Im Gegensatz zu mxied ANOVA können wir in der repeated ANOVA auch nur den einfaktoriellen Fall vorliegen haben. Da klingt jetzt erstmal etwas seltsam, aber wir haben in diesem Fall nur die Zeitpunkte als Faktor vorliegen. Wir messen also nur Katzenflöhe ohne eine Behandlung viermal. Wir schauen also, ob sich die Sprungweite von Katzenflöhen über die Zeitpunkte verändert. Wir vergleichen faktisch also nur die Zeitpunkte untereinander.
R Code [zeigen / verbergen]
repeated_fac1_tbl <- read_excel("data/fleas_complex_data.xlsx",
sheet = "repeated-fac1") |>
select(.id, animal, t1:t4) |>
pivot_longer(cols = t1:t4,
values_to = "jump_length",
names_to = "time_fct") |>
mutate(animal = as_factor(animal),
time_fct = as_factor(time_fct),
jump_length = round(jump_length, 2),
.id = as_factor(.id))Auch hier einmal die Daten im Wide-Format, damit die Daten etwas übersichtlicher sind. Wir haben nur die Katzenflöhe vorligen und schauen uns zwölf einzelne Katzenflöhe an. Jeder Floh springt dann vier Mal und wir messen dann die Sprungweite in [cm]. Wir müssen uns dann noch das Long-Format über die Funktion pivot_longer() in R bauen.
| .id | animal | t1 | t2 | t3 | t4 |
|---|---|---|---|---|---|
| 1 | cat | 23.72 | 37.57 | 36.67 | 25.89 |
| 2 | cat | 17.84 | 29.36 | 22.63 | 15.44 |
| 3 | cat | 24.11 | 30.43 | 30 | 28.14 |
| ... | ... | ... | ... | ... | ... |
| 10 | cat | 28.52 | 32.39 | 33.35 | 23.91 |
| 11 | cat | 25.46 | 29.02 | 22.68 | 19.88 |
| 12 | cat | 20.92 | 34.39 | 26.81 | 24.76 |
So ergibt sich dann auch ein sehr simples Versuchsdesign. Wir lassen jeden Katzenfloh viermal springen und messen jeweils die Sprungweite für jeden Floh.
Damit ergibt sich am Ende auch eine sehr übersichtliche Abbildung. Wir sehen einmal die einzelnen Katzenflöhe. Jede Sprungweite über die vier Zeitpunkte habe ich für die individuellen Flöhe mit einer gestrichelten Linie verbunden. Den globalen Trend kannst du in der farbigen Linie nachvollziehen. Am Anfang springen die Flöhe nicht soweit, steigern sich dann zum zweiten Messzeitpunkt um dann wieder stetig abzunehmen. Einzelne Flöhe zeigen dabei eine größere Variabilität als der globale Trend.
Zweifaktoriell
Eigentlich wird die Besonderheit des repeated Design einer ANOVA erst bei der zweifaktoriellen repeated ANOVA so richtig ersichtlich. Hier haben wir es nämlich zum einen mit einem Behandlungsfaktor zu tun und einer Faktor mit den jeweiligen Messzeitpunkten. Wir haben hier wieder die Katzenflöhe mit unterschiedlichen Fütterungsarten ernährt und dann an vier Zeitpunkten gemessen, wie weit die Flöhe springen. Die Besonderheit ist jetzt, dass jeder Katzenfloh alle Behandlungen durchläuft. Das macht einen gewaltigen Unterschied im Design und in der Auswertung.
R Code [zeigen / verbergen]
repeated_fac2_tbl <- read_excel("data/fleas_complex_data.xlsx",
sheet = "repeated-fac2") |>
select(.id, animal, feeding, t1:t4) |>
pivot_longer(cols = t1:t4,
values_to = "jump_length",
names_to = "time_fct") |>
mutate(animal = as_factor(animal),
feeding = as_factor(feeding),
time_fct = as_factor(time_fct),
jump_length = round(jump_length, 2),
.id = as_factor(.id))In der folgenden Tabelle der Rohdaten der Sprungweiten in [cm] unserer Katzenflöhe siehst du nochmal die Besonderheit der Messwiederholung. Wir messen die Sprungweite von jedem Floh \(.id\) viermal für jede Fütterungsart. Wir müssen natürlich zwischen den Fütterungsarten Zeit lassen, damit wir wieder die Fütterungseffekte raus haben, aber ansonsten wird jeder Floh mit allen Ernährungsformen gefüttert. Hier ist dann auch die \(.id\) Spalte so wichtig, da wir ja wissen müssen welche Fütterungsart zu welchem Floh gehört. Abschließend nutzen wir die Funktion pivot_longer() um uns die Daten passend in R zusammenzubauen.
| .id | animal | feeding | t1 | t2 | t3 | t4 |
|---|---|---|---|---|---|---|
| 1 | cat | ctrl | 34.85 | 50.21 | 50.4 | 40.92 |
| 1 | cat | blood | 50.37 | 57.87 | 66.84 | 45.45 |
| 1 | cat | ketchup | 51.01 | 44.54 | 42.65 | 38.74 |
| 2 | cat | ctrl | 34.24 | 43.85 | 37.34 | 30.92 |
| ... | ... | ... | ... | ... | ... | ... |
| 11 | cat | ketchup | 42.29 | 49.83 | 46.08 | 39.97 |
| 12 | cat | ctrl | 42.88 | 46.23 | 45.93 | 29.81 |
| 12 | cat | blood | 49.9 | 56.72 | 66.54 | 47.85 |
| 12 | cat | ketchup | 46.47 | 52.61 | 44.92 | 33.58 |
In der folgenden Abbildung siehst du dann das Versuchsdesign für einen Katzenfloh nochmal dargestellt. Der Floh beginnt mit der Ernährungsform ketchup und springt dann viermal. Jedes Mal messen wir dann die Sprungweite. Nach einer Pause erhält der gleiche Katzenfloh dann die Ernährungsform blood und springt wieder zeitversetzt viermal. Am Ende ernähren wir dann den Floh noch mit Zuckerwasser ctrl und lassen ihn abschließend viermal springen. Damit springt ein einzelner Katzenfloh über den Versuch hinweg zwölfmal unter verschiedenen Versuchsbedingungen. Wie du siehst, ist dieses Design vermutlich nicht so häufig in den Agrarwissenschaften anzutreffen.
Dann schauen wir uns nochmal den Zusammenhang zwischen der Ernährungsart und den einzelnen Zeitpunkten einmal in einer Abbildung an. Die einzelnen Flöhe sind als getrichelte Linien dargestellt. Die durchgezogenen Linien stellen die globalen Mittelwerte an den jeweiligen Zeitpunkten für die Fütterungsart dar. Wir sehen, dass wir Unterschiede in den Sprungweiten über die Zeitpunkte je nach Ernährungsform vorliegen haben. Dabei haben die Kontrolle und der Ketchup ähnliche Verläufe.
Dreifaktoriell
Jetzt wird es richtig wild. Wir schauen uns jetzt nämlich nicht nur einen Behandlungsfaktor an, sondern haben zwei Behandlungen sowie die Messwiederholungen als Faktor über die Zeit. Wir schauen hier einmal auf die Fütterungsart mit der Kontrolle als Zuckerwasser, der Gabe von Blut sowie der Fütterung mit Ketchup. Dazu kommt dann noch eine Workoutsession oder eben keine Workoutsession. Wir haben damit zwei Behandlungsfaktoren vorliegen. Da wir ja nun jeden Floh wiederholt messen unter jeder Faktorkombination, kommt hier einiges an Sprungweitenmessung auf jeden einzelnen Floh zu. Daher sind repeated ANOVAs selten. Wir müssen ja nicht nur jeden Floh wiederholt für die Fütterungen durchmessen sondern eben auch nochnmal alles wiederholt für die beiden Workoutgruppen.
R Code [zeigen / verbergen]
repeated_fac3_tbl <- read_excel("data/fleas_complex_data.xlsx",
sheet = "repeated-fac3") |>
select(.id, animal, feeding, workout, t1:t4) |>
pivot_longer(cols = t1:t4,
values_to = "jump_length",
names_to = "time_fct") |>
mutate(animal = as_factor(animal),
feeding = as_factor(feeding),
workout = as_factor(workout),
time_fct = as_factor(time_fct),
jump_length = round(jump_length, 2),
.id = as_factor(.id))In der folgenden Tabelle siehst du einmal die Daten im Wide-Format. In den ersten beiden Zeilen siehst du auch schon die Besonderheit. Wir Füttern eben einmal mit Bkut und messen dann viermal die Sprungweite für die beiden Gruppen des Workout. Daher messen wir dann jeden Floh insgesamt vierundzwanzigmal mit den jeweiligen Pausen zwischen den Behandlungskombinationen. Da kommen schon eine Menge Sprünge zusammen.
| .id | animal | feeding | workout | t1 | t2 | t3 | t4 |
|---|---|---|---|---|---|---|---|
| 1 | cat | blood | yes | 57.62 | 67.87 | 68.82 | 61.19 |
| 1 | cat | blood | no | 42.32 | 49.28 | 58.87 | 50.75 |
| 1 | cat | ctrl | yes | 49.53 | 54.91 | 46.01 | 34.32 |
| 1 | cat | ctrl | no | 31.02 | 42.39 | 26.64 | 11.59 |
| ... | ... | ... | ... | ... | ... | ... | ... |
| 12 | cat | ctrl | yes | 48.82 | 52.84 | 41.37 | 36.02 |
| 12 | cat | ctrl | no | 32.85 | 32.18 | 29.1 | 16.41 |
| 12 | cat | ketchup | yes | 59.16 | 55.02 | 39.41 | 41.03 |
| 12 | cat | ketchup | no | 35.11 | 36.99 | 27.43 | 30.36 |
Eigentlich sieht man das Versuchsdesign erst so richtig in der folgenden Abbildung. Hier siehst du einmal wie ein einzelner Floh alle 24 Sprünge nacheinander unter der jeweiligen Faktorkombination durchführt. In einem idealen Setting folgt natürlich jeder Floh eine zufällige Reihenfolge der Faktorkombinationen. Bitte beachte, dass zwischen den einzelnen Behandlungskombinationen nach dem vierten Sprung natürlich auch eine Zeit vergehen muss, damit sich die Behandlung wieder ausgeschlichen hat.
Am Ende schauen wir uns auch die Daten einmal in einer Abbildung an. Hier haben wir dann die beiden Workoutgruppen sowie die drei Fütterungsarten in einer Abbildung vereint. Wir sehen, dass die nicht trainierten Flöhe im allgemeinen nicht so weit springen wir die trainierten Flöhe. Ebenso wirkt, dass die Ernährung mit Blut auf jeden Fall die bessere Sprungleistung bewirkt. Teilweise sehen wir, dass sich die Linien der Mittelwerte überschneiden, so dass wir hier eventuell von einer Interaktion ausgehen müssen. Das schauen wir uns dann aber in den entsprechenden statistischen Tests einmal näher an.
32.4 Mixed ANOVA
Beginnen wir mit der Auswertung unserer beispieldaten mit der mixed ANOVA. Wir kriegen hier natürlich nur eine globale Aussage über den Faktor der Behandlung oder der Zeit. Wenn dich die paarweisen Vergleiche interessieren, dann müsst du noch einen Posthoc-Test anschließen. Mehr Informationen zu der Modellierung findest du auch in dem Kapitel zu gemischten Modellen. Dort erfährst du dann auch, wie du dir den Posthoc-Test bauen kannst. Hier bleiben wir dann bei der ANOVA und den entsprechenden Ergebnissen.
32.4.1 Zweifaktoriell
Als erstes schauen wir uns die Auswertung der zweifaktoriellen mixed ANOVA an. Nochmal zur Erinnerung, wir haben einen Behandlungsfaktor mit den drei Arten der Fütterung, Kontrolle, Blut und Ketchup, vorliegen. Darüber hinaus messen wir die Sprungweite an sechs Zeitpunkten widerholt für jeden Katzenfloh. Wir wollen jetzt wissen, ob die Fütterung einen Einfluss auf die Sprungweiten der Katzenflöhe über die Zeit hat. Dabei sind wir eigentlich nicht so sehr an den Zeitpunkten wie an der Behandlung der Fütterung interessiert. In einem zweifaktoriellen Versuch können wir uns dann auch immer die Interaktion mit anschauen um zu sehen, die Behandlung unabhängig von den Zeitpunkten ist.
Die Standardimplementierung ist eigentlich recht einfahc durchzuführen. Wir haben wie immer erst unseren Messwert \(y\) und dann die Formelschreibweise für unser beiden Faktoren und deren Interaktion. Ich mache mir es hier etwas einfacher und schreibe dann noch die Faktoren mit dem * um etwas Platz zu sparen. Die Besonderheit hier ist jetzt, dass wir noch spezifizieren müssen welche Variable in unseren Daten die wiederholt gemessenen Individuen darstellt und über welchen Faktor diese Individuen dann wiederholt gemessen wurden. In unserem Fall stehen die Individuen in der Spalte .id und wir haben über die Spalte time_fct wiederholt die einzelnen Flöhe gemessen. Dann müssen wir das ganze noch in die Funktion Error() schreiben, damit die ANOVA auch erkennt, dass wir hier den Teil mit den Messwiederholungen darstellen.
R Code [zeigen / verbergen]
aov(jump_length ~ feeding * time_fct + Error(.id/time_fct),
data = mixed_fac2_tbl) |>
tidy() |>
na.omit() |>
select(term, df, statistic, p.value) |>
mutate(p.value = pvalue(p.value))# A tibble: 3 × 4
term df statistic p.value
<chr> <dbl> <dbl> <chr>
1 feeding 2 214. <0.001
2 time_fct 5 116. <0.001
3 feeding:time_fct 10 8.49 <0.001 Ich filtere hier einiges weg, da wir hier viele Informationen erhalten, die wir meistens nicht brauchen. Wir sehen, dass einen signifikanten Einfluss der Fütterung auf die Sprungweite haben, da der p-Wert kleiner ist als das Signifikanzniveau \(\alpha\) gleich 5%. Dies ist auch der Fall für die Zeitpunkte und die Interaktion. Da wir hier eine signifikante Interaktion vorliegen haben, können wir die Fütterungsarten nicht unabhängig über die Zeit interpretieren. Es gibt damit Zeitpunkte an denen die Fütterungen unterschiedliche Effekte zeigen.
Das R Paket {rstatix} liefert mit der Funktion anova_test() auch eine Möglichkeit eine mixed ANOVA zu rechnen. Hier ist es dann mit den Optionen etwas anders. Ich mag ja die Formelschreibweise sehr gerne, da weiß ich dann wo was steht. Auf der anderen Seite mag es aber auch einfacher sein, wenn man einfach die Variablen in den Daten einer Option in der Funktion zuordnet. Hier haben wir den Fall, dass wir sagen müssen, was ist unser Messwert? Der Messwert kommt in die Option dv. Dann brauchen wir die individuellen IDs in der Option wid. Die Faktoren, die wir untereinander vergleichen wollen, kommen dann in die Option between. In unserem Fall sind das die Fütterungsarten. Wir wollen wissen ob es zwischen den Fütterungsarten einen Unterschied gibt. Die Zeit als Messwiederholung kommt dann in die Option within. Wir messen ja die Zeitpunkte innerhalb (eng. within) der Individuen.
R Code [zeigen / verbergen]
rstatix_aov <- anova_test(data = mixed_fac2_tbl,
dv = jump_length, wid = .id,
between = feeding, within = time_fct)
get_anova_table(rstatix_aov)ANOVA Table (type II tests)
Effect DFn DFd F p p<.05 ges
1 feeding 2 21 214.133 1.08e-14 * 0.750
2 time_fct 5 105 115.543 4.79e-41 * 0.824
3 feeding:time_fct 10 105 8.488 5.19e-10 * 0.408Numerisch sind die Werte die gleichen Werte, die wir auch in der aov() Funktion erhalten. Im Hintergrund wird auch eine aov() gerechnet, nur die Übergabe der Variablen ist hier anders. Neben der Signifikanz kriegen wir in der Spalte ges auch noch ein Eta\(^2\) Wert wiedergegeben, den wir klassisch interpretieren können. Da wir es in einer mixed ANOVA mit einem etwas anderen Eta\(^2\) Wert zu tun haben, können wir hier erstmal sagen, dass je gößer der Wert desto stärker der Effekt. Die p-Werte sind auch hier für alle Vergleiche hoch signifikant.
Das R Paket {afex} elraubt im besonderen die Modellierung von Messwiederholungen von Beobachtungen. Faktisch ist das Paket um diese Art der Analyse gebaut. Dementsprechend gibt es auch zwei Arten des Funktionsaufrufs. Wir konzentrieren uns hier aus Gründen der konsistents zu aov() auf die Funktion aov_car(). Am Ende kommen hier die gleichen Werte raus wie schon bei den anderen Tabs. Wir haben hier aber die Möglichkeit etwas mehr zu machen. Alles Weitere dann in der Vingette zu {afex} und dem Posthoc-Test. Wenn du viel mit Messwiederholungen arbeitest, dann ist eine Einarbeitung in {afex} sicherlich sinnvoll. Leider bist du bei {afex} an die Normalverteilung und die Varianzhomogenität gebunden. MAn kann eben nicht alles haben.
R Code [zeigen / verbergen]
aov_car(jump_length ~ feeding * time_fct + Error(.id/time_fct),
data = mixed_fac2_tbl)Anova Table (Type 3 tests)
Response: jump_length
Effect df MSE F ges p.value
1 feeding 2, 21 13.29 214.13 *** .750 <.001
2 time_fct 4.11, 86.33 18.70 115.54 *** .824 <.001
3 feeding:time_fct 8.22, 86.33 18.70 8.49 *** .408 <.001
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Auch hier erhalten wir die gleichen numwerischen Werte wie schon in den vorherigen Tabs zu {base} und {rstatix}. Die Signifikanz ist für alle Faktoren sowie der Interaktion gegeben. Wir erhalten auch hier in der Spalte ges dann das Eta\(^2\) als Effektschätzer zurück. Am Ende finde ich die Ausgabe etwas sperrig.
Das R Paket {afex} erlaubt es auch die Modelle nach der Schreibweise der linearen gemischten Modelle zu definieren. Das ist für mich eigentlich etwas angenehmer, da man sich dann in einer Formelwelt aufhält mehr dazu dann im Kapitel zu gemischten Modellen. Also die Ergebnisse sind die Gleichen. Ich brauche es nur als Referenz für mich zum nachschlagen.
R Code [zeigen / verbergen]
aov_4(jump_length ~ feeding * time_fct + (time_fct|.id),
data = mixed_fac2_tbl)Mit dem R Paket {MANOVA.RM} haben wir auch die Möglichkeit eine mixed ANOVA zu rechnen in dem wir die Funktion RM() nutzen. Das ist natürlich etwas verwirrend vom Namen, aber wenn wir die Optionen richtig mit den Variablen belegen, dann rechnen wir auch eine mixed ANOVA. Hier müssen wir dann angeben in welche Spalte die Beobachtungen sind und was der within Faktor ist. Bei uns ist das dann einfach. Wir nutzen die Variable .id sowie als Faktor der Messwiederholungen time_fct. Der Rest ist eigentlich nicht so wichtig. Wenn du die Funktion in deiner Analyse rechnest, dann würde ich die Option iter auf 1000 setzen.
Am Ende muss sich aber dann doch fragen, ob es eine noch hässlichere Ausgabe einer Funktion gibt. Ich glaube eher nicht, aber das soll die Leistung des Algorithmus nicht schmälern. Wir kämpfen uns da doch mal durch. Manche Informationen sind dann einfach auch wieder über, da wir je eigentlich nur wissen wollen, ob wir einen signifikanten Einfluss der Behandlung feeding vorliegen haben.
R Code [zeigen / verbergen]
RM(jump_length ~ feeding * time_fct, data = mixed_fac2_tbl,
subject = ".id", within = c("time_fct"), iter = 100,
resampling = "Perm", seed = 1234) |>
summary()Call:
jump_length ~ feeding * time_fct
A repeated measures analysis with 1 within-subject factor(s) ( time_fct ) and 1 between-subject factor(s).
Descriptive:
feeding time_fct n Means Lower 95 % CI Upper 95 % CI
1 ctrl t0 8 32.510 27.041 37.979
4 ctrl t1 8 28.544 22.315 34.773
7 ctrl t2 8 24.625 19.272 29.978
10 ctrl t3 8 21.853 17.032 26.673
13 ctrl t4 8 18.774 12.618 24.929
16 ctrl t5 8 5.239 0.156 10.321
2 blood t0 8 40.072 34.510 45.635
5 blood t1 8 28.651 24.154 33.148
8 blood t2 8 25.236 17.146 33.327
11 blood t3 8 18.062 11.397 24.728
14 blood t4 8 20.250 14.413 26.087
17 blood t5 8 17.351 11.442 23.260
3 ketchup t0 8 50.941 46.864 55.018
6 ketchup t1 8 45.538 41.039 50.036
9 ketchup t2 8 34.306 30.118 38.494
12 ketchup t3 8 29.146 23.946 34.346
15 ketchup t4 8 27.706 24.583 30.830
18 ketchup t5 8 31.434 25.615 37.252
Wald-Type Statistic (WTS):
Test statistic df p-value
feeding "396.501" "2" "<0.001"
time_fct "1009.823" "5" "<0.001"
feeding:time_fct "83.636" "10" "<0.001"
ANOVA-Type Statistic (ATS):
Test statistic df1 df2 p-value
feeding "214.133" "1.982" "72.609" "<0.001"
time_fct "115.543" "4.111" "Inf" "<0.001"
feeding:time_fct "8.488" "6.975" "Inf" "<0.001"
p-values resampling:
Perm (WTS)
feeding "<0.001"
time_fct "<0.001"
feeding:time_fct "0.02" Als erstes erhalten wir eine lange Tabelle mit der deskriptiven Statistik der Mittelwerte für jede Faktorkombination der Fütterungsart sowie dem Zeitpunkt. Dazu dann noch die passenden 95% Konfidenzintervalle. Dann kommen drei Tests nacheinander. Das ist auch immer das Problem, wenn du dann drei Ausgaben hast und nicht weißt, welche du nehmen sollst. Ich würde hier die Ergebnisse der ANOVA-Type Statistic (ATS) nehmen und den Rest liegen lassen. Die p-Werte aus dem Permutationstest p-values resampling haben zwar auch ihre Berechtigung bei nicht normalverteilten oder varianzheterogenen Daten aber ich würde dann hier auch eher das R Paket {WRS2} vorziehen.
R Code [zeigen / verbergen]
bwtrim(jump_length ~ feeding * time_fct, id = .id, data = mixed_fac2_tbl)Call:
bwtrim(formula = jump_length ~ feeding * time_fct, id = .id,
data = mixed_fac2_tbl)
value df1 df2 p.value
feeding 186.6213 2 9.6997 0.0000
time_fct 142.2153 5 9.7570 0.0000
feeding:time_fct 4.6665 10 9.4221 0.0133R Code [zeigen / verbergen]
mixed_fac2_tbl |>
ggplot(aes(x = time_fct, y = jump_length,
color = feeding, group = feeding)) +
theme_minimal() +
stat_summary(fun = "mean", geom = "point") +
stat_summary(fun = "mean", geom = "line") +
theme(legend.position = "top") +
scale_color_okabeito()
feeding und time_fct. Die Linien kreuzen sich leicht, es liegt eine Interaktion vor. Eine allgmeine Aussage über die Sprungweite und der Fütterung lässt sich nicht unabhängig vom Zeitpunkt treffen.
32.4.2 Dreifaktoriell
R Code [zeigen / verbergen]
aov(jump_length ~ feeding * workout * time_fct + Error(.id/time_fct),
data = mixed_fac3_tbl) |>
tidy() |>
na.omit() |>
select(term, df, statistic, p.value) |>
mutate(p.value = pvalue(p.value))# A tibble: 7 × 4
term df statistic p.value
<chr> <dbl> <dbl> <chr>
1 feeding 2 318. <0.001
2 workout 1 14.0 0.002
3 feeding:workout 2 79.0 <0.001
4 time_fct 5 113. <0.001
5 feeding:time_fct 10 4.03 <0.001
6 workout:time_fct 5 0.943 0.457
7 feeding:workout:time_fct 10 1.45 0.171 R Code [zeigen / verbergen]
rstatix_aov <- anova_test(data = mixed_fac3_tbl,
dv = jump_length, wid = .id,
between = c(feeding, workout),
within = time_fct)
get_anova_table(rstatix_aov)ANOVA Table (type II tests)
Effect DFn DFd F p p<.05 ges
1 feeding 2 18 318.133 9.03e-15 * 0.786
2 workout 1 18 13.980 2.00e-03 * 0.075
3 time_fct 5 90 112.777 3.34e-37 * 0.849
4 feeding:workout 2 18 79.042 1.22e-09 * 0.477
5 feeding:time_fct 10 90 4.028 1.40e-04 * 0.286
6 workout:time_fct 5 90 0.943 4.57e-01 0.045
7 feeding:workout:time_fct 10 90 1.451 1.71e-01 0.126R Code [zeigen / verbergen]
aov_car(jump_length ~ feeding * workout * time_fct + Error(.id/time_fct),
data = mixed_fac3_tbl)Anova Table (Type 3 tests)
Response: jump_length
Effect df MSE F ges p.value
1 feeding 2, 18 10.55 318.13 *** .786 <.001
2 workout 1, 18 10.55 13.98 ** .075 .002
3 feeding:workout 2, 18 10.55 79.04 *** .477 <.001
4 time_fct 3.88, 69.75 23.44 112.78 *** .849 <.001
5 feeding:time_fct 7.75, 69.75 23.44 4.03 *** .286 <.001
6 workout:time_fct 3.88, 69.75 23.44 0.94 .045 .442
7 feeding:workout:time_fct 7.75, 69.75 23.44 1.45 .126 .193
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Das R Paket {afex} erlaubt es auch die Modelle nach der Schreibweise der linearen gemischten Modelle zu definieren. Das ist für mich eigentlich etwas angenehmer, da man sich dann in einer Formelwelt aufhält mehr dazu dann im Kapitel zu gemischten Modellen. Also die Ergebnisse sind die Gleichen. Ich brauche es nur als Referenz für mich zum nachschlagen.
R Code [zeigen / verbergen]
aov_4(jump_length ~ feeding * workout * time_fct + (time_fct|.id),
data = mixed_fac3_tbl)R Code [zeigen / verbergen]
RM(jump_length ~ feeding * workout * time_fct, data = mixed_fac3_tbl,
subject = ".id", within = c("time_fct"), iter = 100,
resampling = "Perm", seed = 1234) |>
summary()Call:
jump_length ~ feeding * workout * time_fct
A repeated measures analysis with 1 within-subject factor(s) ( time_fct ) and 2 between-subject factor(s).
Descriptive:
feeding workout time_fct n Means Lower 95 % CI Upper 95 % CI
1 ctrl yes t0 4 37.718 19.752 55.683
7 ctrl yes t1 4 34.570 24.732 44.408
13 ctrl yes t2 4 31.290 26.489 36.091
19 ctrl yes t3 4 24.512 9.769 39.256
25 ctrl yes t4 4 26.165 14.107 38.223
31 ctrl yes t5 4 10.075 -4.546 24.696
4 ctrl no t0 4 33.405 21.056 45.754
10 ctrl no t1 4 29.212 14.480 43.945
16 ctrl no t2 4 21.595 8.099 35.091
22 ctrl no t3 4 19.015 1.083 36.947
28 ctrl no t4 4 16.195 8.369 24.021
34 ctrl no t5 4 7.775 -0.910 16.460
2 blood yes t0 4 46.035 35.756 56.314
8 blood yes t1 4 31.495 21.001 41.989
14 blood yes t2 4 25.635 7.979 43.291
20 blood yes t3 4 25.080 20.580 29.580
26 blood yes t4 4 21.587 2.378 40.797
32 blood yes t5 4 18.445 10.144 26.746
5 blood no t0 4 52.303 40.389 64.216
11 blood no t1 4 46.498 39.556 53.439
17 blood no t2 4 36.158 15.701 56.614
23 blood no t3 4 32.672 19.062 46.283
29 blood no t4 4 34.295 14.581 54.009
35 blood no t5 4 29.208 17.262 41.153
3 ketchup yes t0 4 56.428 41.856 70.999
9 ketchup yes t1 4 47.342 25.965 68.720
15 ketchup yes t2 4 36.908 26.465 47.350
21 ketchup yes t3 4 32.180 19.098 45.262
27 ketchup yes t4 4 32.422 26.616 38.229
33 ketchup yes t5 4 35.362 24.276 46.449
6 ketchup no t0 4 57.062 42.258 71.867
12 ketchup no t1 4 52.618 38.726 66.509
18 ketchup no t2 4 39.200 21.873 56.527
24 ketchup no t3 4 39.190 18.245 60.135
30 ketchup no t4 4 31.390 17.678 45.102
36 ketchup no t5 4 31.887 16.088 47.687
Wald-Type Statistic (WTS):
Test statistic df p-value
feeding "1042.469" "2" "<0.001"
workout "13.98" "1" "<0.001"
feeding:workout "133.546" "2" "<0.001"
time_fct "879.365" "5" "<0.001"
feeding:time_fct "47.128" "10" "<0.001"
workout:time_fct "5.119" "5" "0.402"
feeding:workout:time_fct "25.413" "10" "0.005"
ANOVA-Type Statistic (ATS):
Test statistic df1 df2 p-value
feeding "318.133" "1.725" "36.605" "<0.001"
workout "13.98" "1" "36.605" "0.001"
feeding:workout "79.042" "1.725" "36.605" "<0.001"
time_fct "112.777" "3.875" "Inf" "<0.001"
feeding:time_fct "4.028" "6.058" "Inf" "<0.001"
workout:time_fct "0.943" "3.875" "Inf" "0.435"
feeding:workout:time_fct "1.451" "6.058" "Inf" "0.19"
p-values resampling:
Perm (WTS)
feeding "<0.001"
workout "<0.001"
feeding:workout "<0.001"
time_fct "<0.001"
feeding:time_fct "0.16"
workout:time_fct "0.64"
feeding:workout:time_fct "0.39" R Code [zeigen / verbergen]
mixed_fac3_tbl |>
ggplot(aes(x = time_fct, y = jump_length,
color = feeding, linetype = workout,
group = interaction(feeding, workout))) +
theme_minimal() +
stat_summary(fun = "mean", geom = "point") +
stat_summary(fun = "mean", geom = "line") +
theme(legend.position = "top") +
scale_color_okabeito()
feeding sowie workout und time_fct. Die Linien sind in der Reichenfolge verdreht, es liegt eine Interaktion vor. Eine allgmeine Aussage über die Sprungweite und die Fütterung soweie dem Workout lässt sich nicht unabhängig vom Zeitpunkt treffen.
32.5 Repeated measurement ANOVA
Hier eine kurze Betrachtung. Mehr dazu in dem Kapitel zu gemischten Modellen
32.5.1 Einfaktoriell
R Code [zeigen / verbergen]
aov(jump_length ~ time_fct + Error(.id/time_fct), data = repeated_fac1_tbl)|>
tidy() |>
na.omit() |>
select(term, df, statistic, p.value) |>
mutate(p.value = pvalue(p.value))# A tibble: 1 × 4
term df statistic p.value
<chr> <dbl> <dbl> <chr>
1 time_fct 3 19.1 <0.001 R Code [zeigen / verbergen]
rstatix_aov <- anova_test(data = repeated_fac1_tbl,
dv = jump_length, wid = .id, within = time_fct)
get_anova_table(rstatix_aov)ANOVA Table (type III tests)
Effect DFn DFd F p p<.05 ges
1 time_fct 3 33 19.146 2.26e-07 * 0.465R Code [zeigen / verbergen]
aov_car(jump_length ~ time_fct + Error(.id/time_fct), data = repeated_fac1_tbl)Anova Table (Type 3 tests)
Response: jump_length
Effect df MSE F ges p.value
1 time_fct 2.68, 29.52 11.85 19.15 *** .465 <.001
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Das R Paket {afex} erlaubt es auch die Modelle nach der Schreibweise der linearen gemischten Modelle zu definieren. Das ist für mich eigentlich etwas angenehmer, da man sich dann in einer Formelwelt aufhält mehr dazu dann im Kapitel zu gemischten Modellen. Also die Ergebnisse sind die Gleichen. Ich brauche es nur als Referenz für mich zum nachschlagen.
R Code [zeigen / verbergen]
aov_4(jump_length ~ time_fct + (time_fct|.id), data = repeated_fac1_tbl)R Code [zeigen / verbergen]
RM(jump_length ~ time_fct, data = repeated_fac1_tbl,
subject = ".id", within = "time_fct", iter = 100,
resampling = "Perm", seed = 1234) |>
summary()Call:
jump_length ~ time_fct
A repeated measures analysis with 1 within-subject factor(s) ( time_fct ) and 0 between-subject factor(s).
Descriptive:
time_fct n Means Lower 95 % CI Upper 95 % CI
time_fct t1 12 24.081 22.212 25.950
time_fct.1 t2 12 30.642 28.189 33.096
time_fct.2 t3 12 28.107 24.974 31.239
time_fct.3 t4 12 21.410 18.982 23.838
Wald-Type Statistic (WTS):
Test statistic df p-value
time_fct "69.318" "3" "<0.001"
ANOVA-Type Statistic (ATS):
Test statistic df1 df2 p-value
time_fct "19.146" "2.684" "Inf" "<0.001"
p-values resampling:
Perm (WTS)
time_fct "<0.001" R Code [zeigen / verbergen]
with(repeated_fac1_tbl, rmanova(y = jump_length, groups = time_fct, block = .id))Call:
rmanova(y = jump_length, groups = time_fct, blocks = .id)
Test statistic: F = 12.6768
Degrees of freedom 1: 3
Degrees of freedom 2: 21
p-value: 6e-05 Effektschätzer
R Code [zeigen / verbergen]
aov_car(jump_length ~ time_fct + Error(.id/time_fct), data = repeated_fac1_tbl) |>
eta_squared()# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
time_fct | 0.64 | [0.44, 1.00]
- One-sided CIs: upper bound fixed at [1.00].32.5.2 Zweifaktoriell
R Code [zeigen / verbergen]
aov(jump_length ~ feeding * time_fct + Error(.id/(feeding * time_fct)),
data = repeated_fac2_tbl) |>
tidy() |>
na.omit() |>
select(term, df, statistic, p.value) |>
mutate(p.value = pvalue(p.value))# A tibble: 3 × 4
term df statistic p.value
<chr> <dbl> <dbl> <chr>
1 feeding 2 151. <0.001
2 time_fct 3 64.7 <0.001
3 feeding:time_fct 6 9.48 <0.001 R Code [zeigen / verbergen]
rstatix_aov <- anova_test(data = repeated_fac2_tbl,
dv = jump_length, wid = .id,
within = c(feeding, time_fct))
get_anova_table(rstatix_aov)ANOVA Table (type III tests)
Effect DFn DFd F p p<.05 ges
1 feeding 2 22 150.881 1.43e-13 * 0.667
2 time_fct 3 33 64.693 6.57e-14 * 0.459
3 feeding:time_fct 6 66 9.480 1.69e-07 * 0.216R Code [zeigen / verbergen]
aov_car(jump_length ~ feeding * time_fct + Error(.id/feeding * time_fct),
data = repeated_fac2_tbl)Anova Table (Type 3 tests)
Response: jump_length
Effect df MSE F ges p.value
1 feeding 1.84, 20.24 18.95 150.88 *** .667 <.001
2 time_fct 2.43, 26.71 14.23 64.69 *** .459 <.001
3 feeding:time_fct 3.94, 43.39 19.43 9.48 *** .216 <.001
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Das R Paket {afex} erlaubt es auch die Modelle nach der Schreibweise der linearen gemischten Modelle zu definieren. Das ist für mich eigentlich etwas angenehmer, da man sich dann in einer Formelwelt aufhält mehr dazu dann im Kapitel zu gemischten Modellen. Also die Ergebnisse sind die Gleichen. Ich brauche es nur als Referenz für mich zum nachschlagen.
R Code [zeigen / verbergen]
aov_4(jump_length ~ feeding * time_fct + (feeding * time_fct|.id),
data = repeated_fac2_tbl)R Code [zeigen / verbergen]
RM(jump_length ~ feeding * time_fct, data = repeated_fac2_tbl,
subject = ".id", within = c("feeding", "time_fct"), iter = 100,
resampling = "Perm", seed = 1234) |>
summary()Call:
jump_length ~ feeding * time_fct
A repeated measures analysis with 2 within-subject factor(s) ( feeding,time_fct ) and 0 between-subject factor(s).
Descriptive:
feeding time_fct n Means Lower 95 % CI Upper 95 % CI
1 ctrl t1 12 39.066 36.419 41.713
4 ctrl t2 12 43.392 40.083 46.702
7 ctrl t3 12 40.567 37.514 43.619
10 ctrl t4 12 31.857 29.138 34.576
2 blood t1 12 49.712 47.175 52.249
5 blood t2 12 54.938 52.383 57.494
8 blood t3 12 59.444 56.625 62.263
11 blood t4 12 46.676 44.635 48.717
3 ketchup t1 12 43.495 41.083 45.907
6 ketchup t2 12 47.423 44.383 50.464
9 ketchup t3 12 39.304 36.049 42.559
12 ketchup t4 12 35.608 32.582 38.635
Wald-Type Statistic (WTS):
Test statistic df p-value
feeding "242.241" "2" "<0.001"
time_fct "343.172" "3" "<0.001"
feeding:time_fct "77.296" "6" "<0.001"
ANOVA-Type Statistic (ATS):
Test statistic df1 df2 p-value
feeding "150.881" "1.84" "Inf" "<0.001"
time_fct "64.693" "2.428" "Inf" "<0.001"
feeding:time_fct "9.48" "3.944" "Inf" "<0.001"
p-values resampling:
Perm (WTS)
feeding "<0.001"
time_fct "<0.001"
feeding:time_fct "<0.001" 32.5.3 Dreifaktoriell
R Code [zeigen / verbergen]
aov(jump_length ~ feeding * workout * time_fct +
Error(.id/(feeding * workout * time_fct)),
data = repeated_fac3_tbl) |>
tidy() |>
na.omit() |>
select(term, df, statistic, p.value) |>
mutate(p.value = pvalue(p.value))# A tibble: 7 × 4
term df statistic p.value
<chr> <dbl> <dbl> <chr>
1 feeding 2 781. <0.001
2 workout 1 1000. <0.001
3 time_fct 3 84.2 <0.001
4 feeding:workout 2 2.87 0.078
5 feeding:time_fct 6 52.3 <0.001
6 workout:time_fct 3 0.472 0.704
7 feeding:workout:time_fct 6 0.873 0.520 R Code [zeigen / verbergen]
rstatix_aov <- anova_test(data = repeated_fac3_tbl,
dv = jump_length, wid = .id,
within = c(feeding, workout, time_fct))
get_anova_table(rstatix_aov)ANOVA Table (type III tests)
Effect DFn DFd F p p<.05 ges
1 feeding 2 22 781.468 3.69e-21 * 0.700
2 workout 1 11 999.674 3.76e-12 * 0.749
3 time_fct 3 33 84.199 1.52e-15 * 0.518
4 feeding:workout 2 22 2.874 7.80e-02 0.018
5 feeding:time_fct 6 66 52.306 3.41e-23 * 0.428
6 workout:time_fct 3 33 0.472 7.04e-01 0.003
7 feeding:workout:time_fct 6 66 0.873 5.20e-01 0.013R Code [zeigen / verbergen]
aov_car(jump_length ~ feeding * workout * time_fct +
Error(.id/feeding * workout * time_fct),
data = repeated_fac3_tbl)Anova Table (Type 3 tests)
Response: jump_length
Effect df MSE F ges p.value
1 feeding 1.65, 18.15 10.60 781.47 *** .700 <.001
2 workout 1, 11 17.54 999.67 *** .749 <.001
3 time_fct 2.68, 29.52 27.90 84.20 *** .518 <.001
4 feeding:workout 1.71, 18.84 21.94 2.87 + .018 .088
5 feeding:time_fct 3.96, 43.61 21.16 52.31 *** .428 <.001
6 workout:time_fct 2.54, 27.91 15.60 0.47 .003 .674
7 feeding:workout:time_fct 3.57, 39.22 24.44 0.87 .013 .479
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Das R Paket {afex} erlaubt es auch die Modelle nach der Schreibweise der linearen gemischten Modelle zu definieren. Das ist für mich eigentlich etwas angenehmer, da man sich dann in einer Formelwelt aufhält mehr dazu dann im Kapitel zu gemischten Modellen. Also die Ergebnisse sind die Gleichen. Ich brauche es nur als Referenz für mich zum nachschlagen.
R Code [zeigen / verbergen]
aov_4(jump_length ~ feeding * workout * time_fct +
(feeding * workout * time_fct|.id),
data = repeated_fac3_tbl)R Code [zeigen / verbergen]
RM(jump_length ~ feeding * workout * time_fct, data = repeated_fac3_tbl,
subject = ".id", within = c("feeding", "workout", "time_fct"), iter = 100,
resampling = "Perm", seed = 1234) |>
summary()Call:
jump_length ~ feeding * workout * time_fct
A repeated measures analysis with 3 within-subject factor(s) ( feeding,workout,time_fct ) and 0 between-subject factor(s).
Descriptive:
feeding workout time_fct n Means Lower 95 % CI Upper 95 % CI
1 blood yes t1 12 54.481 51.441 57.520
7 blood yes t2 12 61.332 58.027 64.636
13 blood yes t3 12 65.201 61.212 69.190
19 blood yes t4 12 59.774 56.497 63.052
4 blood no t1 12 38.211 35.494 40.928
10 blood no t2 12 44.127 42.030 46.223
16 blood no t3 12 50.966 48.882 53.050
22 blood no t4 12 41.111 38.492 43.730
2 ctrl yes t1 12 45.737 42.965 48.509
8 ctrl yes t2 12 49.628 46.124 53.131
14 ctrl yes t3 12 45.892 43.201 48.582
20 ctrl yes t4 12 29.461 26.105 32.817
5 ctrl no t1 12 32.153 30.342 33.963
11 ctrl no t2 12 34.882 32.284 37.479
17 ctrl no t3 12 31.833 28.145 35.522
23 ctrl no t4 12 16.332 13.308 19.357
3 ketchup yes t1 12 51.518 48.364 54.672
9 ketchup yes t2 12 54.800 51.450 58.150
15 ketchup yes t3 12 45.144 42.834 47.455
21 ketchup yes t4 12 39.668 36.577 42.759
6 ketchup no t1 12 33.944 31.187 36.701
12 ketchup no t2 12 39.954 37.452 42.456
18 ketchup no t3 12 29.016 26.209 31.823
24 ketchup no t4 12 22.848 19.392 26.303
Wald-Type Statistic (WTS):
Test statistic df p-value
feeding "2004.945" "2" "<0.001"
workout "999.674" "1" "<0.001"
feeding:workout "6.011" "2" "0.05"
time_fct "234.391" "3" "<0.001"
feeding:time_fct "511.924" "6" "<0.001"
workout:time_fct "2.2" "3" "0.532"
feeding:workout:time_fct "24.596" "6" "<0.001"
ANOVA-Type Statistic (ATS):
Test statistic df1 df2 p-value
feeding "781.468" "1.65" "Inf" "<0.001"
workout "999.674" "1" "Inf" "<0.001"
feeding:workout "2.874" "1.713" "Inf" "0.065"
time_fct "84.199" "2.684" "Inf" "<0.001"
feeding:time_fct "52.306" "3.965" "Inf" "<0.001"
workout:time_fct "0.472" "2.538" "Inf" "0.67"
feeding:workout:time_fct "0.873" "3.565" "Inf" "0.469"
p-values resampling:
Perm (WTS)
feeding "<0.001"
workout "<0.001"
feeding:workout "0.17"
time_fct "<0.001"
feeding:time_fct "<0.001"
workout:time_fct "0.61"
feeding:workout:time_fct "0.17"
Anwendungsbeispiel: foo (3-faktoriell)
R Code [zeigen / verbergen]
gurke_tbl <- read_excel("data/wachstum_gurke.xlsx") |>
clean_names() |>
select(group = versuchsgruppe, t1:t14) |>
mutate(group = as_factor(group)) | group | t1 | t2 | t3 | t4 |
|---|---|---|---|---|
| Proloog L | 5.5 | 6.1 | 7.4 | 8.9 |
| Proloog L | 4.6 | 5.1 | 6.4 | 5.7 |
| Proloog L | 5.3 | 5.8 | 6.8 | 8.3 |
| ... | ... | ... | ... | ... |
| Katrina D | 0.86 | 0.87 | 0.9 | 1.1 |
| Katrina D | 0.62 | 0.67 | 0.7 | 0.71 |
| Tageslänge | 13.98 | 14.05 | 14.12 | 14.18 |
R Code [zeigen / verbergen]
gurke_time_len_tbl <- gurke_tbl |>
filter(str_detect(group, "L$")) |>
rownames_to_column(".id") |>
mutate(group = factor(group, labels = c("Proloog", "Quarto", "Katrina"))) |>
pivot_longer(cols = t1:t14,
values_to = "length",
names_to = "time") |>
mutate(time_fct = as_factor(time),
time_num = as.numeric(time_fct)) |>
filter(length != 0)R Code [zeigen / verbergen]
ggplot(gurke_time_len_tbl,
aes(time_fct, length, color = group,
group = group)) +
theme_minimal() +
geom_point() +
stat_summary(fun = "mean", geom = "line") +
stat_summary(fun = "median", geom = "line", linetype = 2) +
scale_color_okabeito() +
theme(legend.position = "top")+
labs(x = "Zeitpunkte der Messung", y = "Gurkenlänge [cm]",
color = "Sorte")
R Code [zeigen / verbergen]
ggplot(gurke_time_len_tbl,
aes(time_fct, log(length), color = group,
group = .id)) +
theme_minimal() +
geom_point2() +
geom_line(linetype = 11) +
scale_color_okabeito() +
theme(legend.position = "top")+
labs(x = "Zeitpunkte der Messung", y = "Gurkenlänge [cm] (log-scale)",
color = "Sorte")
R Code [zeigen / verbergen]
aov_car(length ~ group * time_fct + Error(.id/time_fct),
data = gurke_time_len_tbl)Warning: Missing values for 2 ID(s), which were removed before analysis:
15, 2
Below the first few rows (in wide format) of the removed cases with missing data.
.id group t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14
# 7 15 Katrina 3.3 3.3 3.4 3.6 3.5 3.5 3.4 3.7 3.5 3.6 3.6 3.6 2.8 NA
# 11 2 Proloog 4.6 5.1 6.4 5.7 5.5 5.2 5.0 5.0 4.5 NA NA NA NA NAContrasts set to contr.sum for the following variables: groupAnova Table (Type 3 tests)
Response: length
Effect df MSE F ges p.value
1 group 2, 13 74.07 25.08 *** .704 <.001
2 time_fct 1.39, 18.11 33.11 92.69 *** .732 <.001
3 group:time_fct 2.79, 18.11 33.11 15.46 *** .477 <.001
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG R Code [zeigen / verbergen]
aov_4(length ~ group * time_fct + (time_fct|.id),
data = gurke_time_len_tbl)Warning: Missing values for 2 ID(s), which were removed before analysis:
15, 2
Below the first few rows (in wide format) of the removed cases with missing data.
.id group t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14
# 7 15 Katrina 3.3 3.3 3.4 3.6 3.5 3.5 3.4 3.7 3.5 3.6 3.6 3.6 2.8 NA
# 11 2 Proloog 4.6 5.1 6.4 5.7 5.5 5.2 5.0 5.0 4.5 NA NA NA NA NAContrasts set to contr.sum for the following variables: groupAnova Table (Type 3 tests)
Response: length
Effect df MSE F ges p.value
1 group 2, 13 74.07 25.08 *** .704 <.001
2 time_fct 1.39, 18.11 33.11 92.69 *** .732 <.001
3 group:time_fct 2.79, 18.11 33.11 15.46 *** .477 <.001
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Referenzen
Friedrich S, Konietschke F, Pauly M. 2019. Resampling-based analysis of multivariate data and repeated measures designs with the R package MANOVA. RM.
Gueorguieva R, Krystal JH. 2004. Move over anova: progress in analyzing repeated-measures data andits reflection in papers published in the archives of general psychiatry. Archives of general psychiatry 61: 310–317.