Letzte Änderung am 22. April 2025 um 18:59:00

“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.

Abbildung 32.1— Konzeptionelle Darstellung der mixed und repeated ANOVA in einem Flussdiagramm mit unterschiedlichen Individuen für einen Behandlungsfaktor mit drei Leveln als Arten der Behandlung. (A) In der mixed ANOVA werden pro Art der Behandlung ein Individum wiederholt gemessen. (B) In der repeated ANOVA wird ein Individum wiederholt für jede Art der Behandlung gemessen. [Zum Vergrößern anklicken]

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 Formelschreibweise y ~ 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 Funktion anova_test() die es erlaubt einfacher auf die Funktionalität der Funktion aov() 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.

Tabelle 32.1— Rohdaten der Sprungweiten in [cm] gemessen an sechs Zeitpunkten von drei Fütterungsarten im Wide-Format für die bessere Übersicht.
tinytable_dziyyaulolcr2qe0ohqh
.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.

Abbildung 32.2— Versuchsdesign der zweifaktoriellen mixed ANOVA. Der erste Faktor ist die Fütterungsart 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.

Abbildung 32.3— Darstellung des Zusammenhangs zwischen der Fütterung von Katzenflöhen und der gemessenen Sprungweite in [cm] über sechs Messzeitpunkte. Die Linie stellt den jeweiligen Mittelwert zu dem Zeitpunkt und der Fütterungsart dar. [Zum Vergrößern anklicken]

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.

Tabelle 32.2— Rohdaten der Sprungweiten in [cm] gemessen an sechs Zeitpunkten von drei Fütterungsarten sowie dem Status des Workouts im Wide-Format für die bessere Übersicht.
tinytable_got1x4cbpl5v9leiwp2g
.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.

Abbildung 32.4— Versuchsdesign der dreifaktoriellen mixed ANOVA. Der erste Faktor ist die Fütterungsart 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.

Abbildung 32.5— Darstellung des Zusammenhangs zwischen der Fütterung von Katzenflöhen sowie der Trainingszustand und der gemessenen Sprungweite in [cm] über sechs Messzeitpunkte. Die Linie stellt den jeweiligen Mittelwert zu dem Zeitpunkt und der Fütterungsart sowie des Trainings dar. [Zum Vergrößern anklicken]

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.

Tabelle 32.3— Rohdaten der Sprungweiten in [cm] gemessen an vier Zeitpunkten im Wide-Format für die bessere Übersicht.
tinytable_qpcdrf6i1r67y1i4apfz
.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.

Abbildung 32.6— Versuchsdesign der einfaktoriellen repeated ANOVA. Jeder Katzenfloh springt viermal und wird dabei gemessen. Keine weiteren Faktoren liegen vor. [Zum Vergrößern anklicken]

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.

Abbildung 32.7— Darstellung des Zusammenhangs zwischen der gemessenen Sprungweite in [cm] und den vier Messzeitpunkten. Die gestrichelte Linie verbindet die einzelnen individuellen Katzenflöhe. Die durchgezogene Linie stellt den Mittelwert zu dem jeweiligen Zeitpunkt dar. [Zum Vergrößern anklicken]

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.

Tabelle 32.4— Rohdaten der Sprungweiten in [cm] nach unterschiedlichen Fütterungsarten gemessen an vier Zeitpunkten im Wide-Format für die bessere Übersicht. Hier erhält jeder Floh jede Behandlung als Fütterung einmal.
tinytable_p23zeyiyzc4ozr31vagm
.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.

Abbildung 32.8— Versuchsdesign der zweifaktoriellen repeated ANOVA. Jeder Katzenfloh springt viermal unter jeder Fütterungsart und wird dabei gemessen. Damit messen wir einen Floh insgesamt zwölfmal. [Zum Vergrößern anklicken]

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.

Abbildung 32.9— Darstellung des Zusammenhangs zwischen der gemessenen Sprungweite in [cm] und den vier Messzeitpunkten sowie den Fütterungsarten. Die gestrichelte Linie verbindet die einzelnen individuellen Katzenflöhe. Die durchgezogene Linie stellt den Mittelwert zu dem jeweiligen Zeitpunkt dar. [Zum Vergrößern anklicken]

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.

Tabelle 32.5— Rohdaten der Sprungweiten in [cm] nach unterschiedlichen Fütterungsarten und Workoutstatus gemessen an vier Zeitpunkten im Wide-Format für die bessere Übersicht. Hier erhält jeder Floh jede Behandlungskombination aus Fütterung und Workout einmal.
tinytable_y6hms71jq1u0bukd0int
.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.

Abbildung 32.10— Versuchsdesign der dreifaktoriellen repeated ANOVA. Jeder Katzenfloh springt viermal unter jeder Kombination von Fütterungsart sowie Workoutstatus und wird dabei gemessen. Damit messen wir einen Floh insgesamt vierundzwanzigmal. [Zum Vergrößern anklicken]

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.

Abbildung 32.11— Darstellung des Zusammenhangs zwischen der gemessenen Sprungweite in [cm] und den vier Messzeitpunkten sowie den Fütterungsarten und Workoutstatus. Die Linienarten stellen den Workoutstatus dar. [Zum Vergrößern anklicken]

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.408

Numerisch 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.

Kommen wir nun zum letzten Paket, was ich hier noch für die zweifaktorielle mixed ANOVA zeigen möchte. Es handelt sich um {WRS2}, welches dann auch erlaubt die zweifaktorielle mixed ANOVA zu rechnen, wenn wir die Annahme der Normalverteilung an den Messwert oder keine Varianzhomogenität in den Gruppen vorliegen haben. Es geht dann auch beides. Das ist dann auch der große Vorteil von {WRS2}. Das Modell sieht genauso aus wie auch bei den anderen Anwendungen. Wir müssen hier nur einmal in der Option id angeben, welche Spalte unsere ID ist.

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.0133

Die Ergebnisse sind dann vom p-Wert inhaltlich gleich wie bei den anderen Anwendungen. Das habe ich hier auch erwartet, denn unsere Sprungweiten habe ich so gebaut, dass die Sprungweiten normalverteilt sind. Nicht immer können wir aber in unseren Daten davon ausgehen, dass wir normalverteilte Messwerte vorliegen haben. Da ist diese Alternative hier sehr gut.

Wir haben ja in den obigen Anwendungen der mixed ANOVA immer eine signifikante Interaktion gefunden. Daher wollen wir uns einmal in einem Liniendiagramm der Mittelwerte die Sachlage näher anschauen. Wenn wir eine Interaktion erwarten, dann sollten sich die Linien der Mittelwerte kreuzen. Wir haben ja auch in allen Analysen eine signifkante Interaktion voriegen. Wenn wir jetzt auf die Abbildung schauen, dann sehen wir auch, dass sich die Lienien für die Fütterung mit Blut und der Kontrolle kreuzen. Wir haben hier eine Interaktion vorliegen. Auch sind die Unterschiede zwischend den Fütterungen nicht immer parallel zueinander mit dem gleichen Effekt oder eben Abstand.

R 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()
Abbildung 32.12— Liniendigramm mit den Mittelwerten aller Faktorkombinationen der beiden Faktoren feeding und time_fct. Die Linien kreuzen sich leicht, es liegt eine Interaktion vor. Eine allgemeine Aussage über die Sprungweite und der Fütterung lässt sich nicht unabhängig vom Zeitpunkt treffen. Die Interaktion ist aber nicht sehr stark.

32.4.2 Dreifaktoriell

Jetzt betrachten wir noch die dreifaktorielle mixed ANOVA. Dieses Design kommt dann vor, wenn du zwei Behandlungsfaktoren sowie einen Zeitfaktor in deinen Daten vorliegen hast. In unserem Fall ist der erste Faktor die Fütterungsart mit den Gruppen Kontrolle, Blut und Ketchup sowie als zweiten Faktor der Status des Workouts des einzlenen Flohs. Wir schauen hier nur Katzenflöhe an. Dann messen wir die Sprungweite jedes einzelnen Flohs nochmal zu sechs Zeitpunkten. Hier ist es wieder wichtig darauf zu achten, was mit einem einzelnen Floh passiert und gemessen wird. Wir messen hier jede Faktorkombination und Wiederholung mit einem neuen Floh. Schaue dir dazu nochmal die Übersicht des Designs in dem Abschnitt zu den Daten in der Abbildung 32.4 an. Auch hier ist es dann manchmal von Vorteil zu wissen, wechler Faktor einem eigentlich am meisten interessiert. Wir betrachten hier jetzt alle drei Faktoren und deren Interaktionen als gleichwertig. Wir können die dreifaktorielle mixed ANOVA nicht mehr in dem R Paket {WRS2} rechnen. Das Paket biette keine Funktion dafür an.

Beginnen wir wieder mit der Standardfunktion in {base} mit aov() und deren Anwendung. Hier kommt dann die verkürzende Formelschreibweise mit dem Stern * zum Einsatz. Ich kann mir nicht merken was alle Faktorkombinationen zusammen mit allen Interaktionen sind. Daher nutze ich dann einmal die Schreibweise, die mir diese Überlegung abnimmt. Dann muss ich noch in dem Fehlerterm mitteilen wer .id in welchen Faktor time_fct wiederholt wird. Damit ist die dreifaktorielle mixed ANOVA der zweifaktoriellen mixed ANOVA sehr ähnlich. Wir erhalten dann folgende Ausagen, die ich nochmal aufgereinigt habe, da ich nicht alle Spalten und Informationen brauche. Auch mag ich den p-Wert gerne etwas übersichtlicher.

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  

Wir sehen sofort, dass die beiden Faktoren der Fütterungsart sowie des Workout signifikant sind. Auch haben wir eine signifikante Interaktion zwischen den beiden Faktoren vorliegen. Die Gruppen der Fütterung verhalten sich in den beiden Gruppen des Workouts nicht gleich. Auch ist die Zeit signifkant, es passiert also etwas mit der Sprungweite über die Zeit. Im Weiteren haben wir eine signifikante Interaktion zwischen der Fütterungsart und der Zeit, so dass sich hier auch über die Zeit die Gruppen der Fütterung anders verhalten. Die anderen Interaktionen sind nicht signifikant.

Jetzt kürze ich ein wenig ab, da die Ergebnisse ja gleich sind in {rstatix} zu der mixed ANOVA mit aov(). Die Funktionen sind ja identisch nur das sich die Eingabe unterscheidet. Hier sieht man dann eventuell besser, dass wir zwischen den Fütterungsarten und den Workouts vergleichen wollen und wir als Wiederholung die .id der einzelnen Flöhe vorliegen haben. Wir erhalten hier nochmal als Effektschätzer das Eta\(^2\) wiedergegeben, was wirklich ein Vorteil ist. Ansonsten würde ich fast die aov() Funktion vorziehen. Am Ende filtere ich wieder um dann eine etwas ausgeräumte Ausgabe zu haben.

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) |> 
  as_tibble() |> 
  select(Effect, F, p, ges) |> 
  mutate(p = pvalue(p))
# A tibble: 7 × 4
  Effect                         F p        ges
  <chr>                      <dbl> <chr>  <dbl>
1 feeding                  318.    <0.001 0.786
2 workout                   14.0   0.002  0.075
3 time_fct                 113.    <0.001 0.849
4 feeding:workout           79.0   <0.001 0.477
5 feeding:time_fct           4.03  <0.001 0.286
6 workout:time_fct           0.943 0.457  0.045
7 feeding:workout:time_fct   1.45  0.171  0.126

Bei dem R Paket {afex} ist es ähnlich wie schon bei dem R Paket {rstatix}. Hier ist dann die Formelschreibweise ähnlich wie schon bei aov(). Das macht die Nutzung angenehmer für mich. Als Besonderheit können wir dann auch {afex} nutzen wie ein lineares gemischtes Modell. Intern ruft dann aber auch {afex} das Paket {car} auf und nutzt dessen Algorithmus. Daher ist hier die Nutzung der Type II ANOVA und Type III ANOVA einfacher umzusetzen. Ich mag das Paket gerne, da es auch speziell für Messwiederholungen gebaut wurde.

R 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)

Am Ende stelle ich noch das R Paket {MANOVA.RM} vor. Auch hier haben wir dann wieder zu viel Ausgabe. Ich finde hier wird dann deutlich, dass die deskriptive Statistik über ist. Ich will doch nur eine mixed ANOVA rechnen. In der Funktion RM() spezifizieren wir ja nur den within Faktor, der sich hier ja auch bei einer dreifaktoriellen mixed ANOVA nicht ändert. Wir haben auch hier wieder drei Tests, wobei ich dann auch nur die Ergebnisse der ANOVA-Type Statistic (ATS) nehmen und den Rest liegen lassen würde.

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"    

Nachdem wir ja unsere mixed ANOVA mit einem der vorgeschlagenen Algorithmen gerechnet haben, sehen wir dass wir eine signifikante Interaktion vorliegen haben. Zum einen ist der Interaktionsterm feeding:workout signifikant sowie der Interaktionsterm feeding:time_fct auch. Damit haben wir andere Mittelwerte der Sprungweite für die Fütterungsarten zu den beiden Workouts sowie unterschiedliche Mittelwerte der Sprungweiten über die Zeit und Fütterungsarten. Ja und das sehen wir dann auch in der folgenden Abbildung. Bei der Kontrollbehandlung und der Ketchupfütterung hat das Workout exakt einen anderen Effekt. Bei der Kontrolle springen die nicht trainierten Flöhe weniger weit. Bei der Ernährung mit Ketchup haben das Training einen negativen Effekt. Eine Drehung des Effekts des Workouts nach der Fütterung. Auch kreuzen die Mittelwerte der Sprungweiten über die Zeit nach der Fütterungsart.

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()
Abbildung 32.13— Liniendigramm mit den Mittelwerten aller Faktorkombinationen der beiden Faktoren feeding sowie workout und time_fct. Die Linien sind in der Reichenfolge verdreht, es liegt eine Interaktion vor. Auch kruezen sich die Fütterungsarten über die Zeit. 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

Achtung, bitte beachten!

Willst du wirklich eine repeated ANOVA rechnen oder doch eher eine mixed ANOVA? Der Unterschied ist subtil und hängt davon ab, welche Behandlungskombinationen deine einzelnen Beobachtungen erhalten haben. Schaue im Zweifel nochmal oben in der Einleitung um dich zu versichern. Häufig möchtest du nämlich eine mixed ANOVA rechnen.

Schauen wir also einmal was die repeated ANOVA einmal macht. Hier nochmal ganz wichtig zu Anfang. Wir werden jetzte jede Beobachtung mit jeder Faktorkombination zu jedem Zeitpunkt messen. Das heißt in unserem Konkreten Fall, dass wir bei jeden Katzenfloh die Sprungweite zu jedem Zeitpunkt messen. Das macht ja auch erstmal Sinn. Das wäre dann die einfaktorielle repated ANOVA. Wir haben sonst auch keinen weiteren Faktor im Experiment. Wenn wir jetzt noch die Fütterungsart als Faktor hinzunehmen, dann füttern wir jeden Floh mit allen drei Arten der Fütterung. Jeder Floh erhält Zuckerwasser, springt dann mehrfach und wir messen die Sprungweite. Das machen wir dann mit dem selben Floh nochmal für die Ernährungsart Blut und Ketchup. Wenn dann noch ein zweiter Behandlungsfaktor hinzukommt wie das Workout, dann müssen wir unter jeder Workoutbedinung einmal alle drei Fütterungsarten die Sprungweiten für einen Floh durchmessen. Das heißt wir messen dann einen einzelnen Floh ganz schön oft. Im Bereich der Agrarwissenschaften mag diese Art der wiederholten Behandlung eher selten vorkommen. Daher hier auch ausdrücklich die Warnung, sich nochmal klar zu machen, wie du die einzelnen Beobachtungen nun bemisst. Ich habe auch erst ein wenig gebraucht um die repeated ANOVA von Design her zu verstehen.

32.5.1 Einfaktoriell

Beginnen wir mit der einfaktoriellen repeated ANOVA. Hier haben wir die Fragestellung vorliegen, ob die Zeit einen Einfluss auf die Sprungweiten der Katzenflöhe hat. Wir haben sonst keinen weiteren Behandlungsfaktor. Wir messen einfach nur die Sprungweite der Katzenflöhe wiederholt und schauen, ob sich die Sprungweiten über die Zeit ändern. Das Experiment ist eher ungewöhnlich, da wir hier ja nur wiederholt messen und sonst nichts tun. Faktisch ist es ein proof-of-principle, ob wir überhaupt einen Effekt der Zeit auf die Sprungweite haben.

Jetzt ist im Prinzip die Schreiweise das Besondere. Wir haben hier jetzt nur die Zeit als einen Einflussfaktor. Den packen wir jetzt einmal vorne in unser Modell und dann einmal hinten in den Fehlerterm. Dann müssen wir auch hier sagen, dass wir die Beobachtungen wiederholt im Faktor Zeit gemessen haben. Dann können wir auch schon die ANOVA rechnen.

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 

Wir sehen, dass der p-Wert kleiner ist als das Signifikanzniveau \(\alpha\) gleich 5% und damit haben wir einen signifkanten Einfluss der Zeit auf die Sprungweite der Katzenflöhe. Mehr testen wir hier auch nicht. Wir haben ja auch nur die Zeit als Faktor in unserem Modell.

Wenn du dann noch einen Effektschätzer brauchst, dann kannst du das R Paket {effectsize} mit der Funktion eta_squared() nutzen. Wir wollen dann noch das generalisierte Eta\(^2\) nutzen, da wir die Zeitpunkte nur beobachten und nicht als experimentellen Faktor festgelegt haben. Dazu dann auch mehr in der Arbeit von Olejnik und Algina (2003) mit dem Titel Generalized eta and omega squared statistics: measures of effect size for some common research designs.

R Code [zeigen / verbergen]
aov(jump_length ~ time_fct + Error(.id/time_fct), data = repeated_fac1_tbl) |> 
  eta_squared(generalized = "time_fct")
# Effect Size for ANOVA (Type I)

Group        | Parameter | Eta2 (generalized) |       95% CI
------------------------------------------------------------
.id:time_fct |  time_fct |               0.46 | [0.22, 1.00]

- Observed variables: time_fct
- One-sided CIs: upper bound fixed at [1.00].

Wir können auch das R Paket {rstatix} nutzen. Hier rufen wir am Ende dann doch auch nur die aov() Funktion auf, schreiben aber die Variablen in unseren Daten etwas anders in die Funkion anova_test(). Du musst eben wissen, dass hier die abhängige Variable dv die Sprungweite ist, du dann die Beobachtungen in wid schreibst und die Zeit in die Option within. Am Ende ist es dann auch hier eine Geschmackssache, welche Funktion dir besser gefällt und passt.

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.465

Es kommt das Gleiche numerisch raus, wie schon bei aov(). Wir erhalten hier aber zusätzlich noch einen Effektschätzer in der Spalte ges, was dem Eta\(^2\) entspricht. Demensprechend ein Vorteil der Funktion.

Das R Paket {afex} kommt bei einer so simplen einfaktoriellen repeated ANOVA dann auf die gleichen Zahlen wie schon die beiden vorherigen Algorithmen. Auch erhälst du bei der Funktion aov_car() einen Wert für Eta\(^2\) in der Spalte ges wiedergegeben. Ich bevorzuge dann ja die Schreibweise hier von {afex}, da ich diese hier etwas klarer finde als bei {rstatix}.

R 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)

Am Ende dann noch die Funkton RM() aus dem R Paket. Erstmal erhalten wir auch hier weder drei statistische Tests, obwohl wir nur nach einem gefragt haben. Ich würde auch hier nur die Ausageb von ANOVA-Type Statistic (ATS) betrachten. Auch hier sehen wir, dass wir in einem einfaktoriellen Fall die gleichen Werte wie in den anderen R Paketen erhalten. Hir haben wir dann aber noch mehr Informatione, die wir dann brauchen oder auch nicht.

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"  

32.5.2 Zweifaktoriell

Die zweifaktorielle repeated ANOVA wird jetzt schon interessanter. Wir haben jetzt nämlich neben der Zeit als Faktor noch eine Behandlung vorliegen. In unserem Beispiel erhält jeder Floh drei Arten der Fütterung wiederholt gefüttert. Zum einen erhält jeder Floh ein Zuckerwasser, dann Blut und zum Schluss Ketchup als Nahrung. Die Reihenfolge sollte sich dann von Floh Floh noch unterscheiden und die Zeiten zwischen den Behandlungen dann so lang gewählt werden, dass die Ernährung förmlich ausgewaschen ist. Die Idee ist, dass die vorherige Ernährungsform keinen Einfluss auf die neue Ernährung hat. Dann müssen wir das Ganze ja auch noch wiederholt über vier Zeitpunkte die Sprungweite messen. Daher misst du hier ganz schön viel. Häufig haben wir dieses Design nicht in den Agrarwissenschaften vorliegen. Wir können aber im Bereich der Nutztierwissenschaften schon das ein oder andere Mal so einen Fall vorliegen haben, wo wir Tiere wiederholt mit einer Behandlung verarzten.

Die Frage ist jetzt immer, wie schreibe ich jetzt das Modell richtig? Der koplizierte Teil ist eigentlich der Fehlerterm Error() aber auch der Teil davor. Ich nutze hier wieder die Kurzschreibweise mit dem Stern * um alle Faktorkombinationen plus Interaktion im Modell zu haben. Das wird eben dann ab zwei Faktoren etwas groß. So ist es dann weniger anfällig für Fehler. Dann musst du noch darauf achten, dass du in dem Fehlerterm die Klammern richtig setzt, sonst funktioniert es dort auch nicht. Erst kommt der Spaltenname für die einzelnen Flöhe .id und dahinter nach dem Slash / dann die Faktorkombinationen der Behandlungen und der Zeit. Dann räume ich die Ausgabe nochmal vollkommen auf, damit ich auch was sehe.

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 

Interpreation fehlt

Dann können wir uns auch den Effektschätzer Eta^\(^2\) berechnen lassen. Hier müssen wir dann das ANOVA Modell einmal in die Funktion eta_squared() stecken und darauf achten, dass wir angeben, dass wir die Zeitpunkte nur beobachtet haben und nicht experimentell festgelegt. Deshalb kommt time_fct noch in die Option generalized mit rein.

R Code [zeigen / verbergen]
aov(jump_length ~ feeding * time_fct + Error(.id/(feeding * time_fct)), 
    data = repeated_fac2_tbl) |> 
  eta_squared(generalized = c("time_fct"))
# Effect Size for ANOVA (Type I)

Group                |        Parameter | Eta2 (generalized) |       95% CI
---------------------------------------------------------------------------
.id:feeding          |          feeding |               0.48 | [0.20, 1.00]
.id:time_fct         |         time_fct |               0.40 | [0.15, 1.00]
.id:feeding:time_fct | feeding:time_fct |               0.13 | [0.00, 1.00]

- Observed variables: time_fct
- One-sided CIs: upper bound fixed at [1.00].

Das R Paket {rstatix} hat eine andere Schreibweise des Modells. Im Prinzip schreiben wir hier wieder nicht das Modell auf, sondern weisen die Spaltennamen dann den jeweiligen Optionen zu. Hier werden dann auch alle Interaktionen gerechnet, du kannst die Interaktionen aller Faktorkombinationen nicht abstellen. Wir setzten dann unseren Endpunkt jump_length in die Option dv und die einzelnen Beobachtungen in die Option wid. Dann müssen wir noch die beiden Faktoren der Zeit und der Beahndlung in die Option within einfügen. Dann haben wir auch schon alles zusammen. Was hier etwas doof ist, ist das wir hier die Formelschreibweise in R verlassen, die wir dann aber überall sonst nutzen. Dafür erhalten wir dann aber das Eta\(^2\) gleich als ges wieder aus der Funktion zurück.

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.216

Wir erhalten die gleichen Ergebnisse wie schon bei der Funktion aov() wiedergegeben. Das überrascht auch nicht, den die Berechnung ist ja die gleiche Berechnung wie auch in aov() nur eben dann mit einem anderen Funktionsaufruf.

Das R Paket {afex} nutzt die Formelschreibweise und kann auch genau so bedient werden wie schon die Funktion aov(). Dafür können wir aber mit der Option type auch anderen Typen der ANOVA nutzen und nicht nur die Type I ANOVA in aov(). Das ist eigentlich die große Stärke von dem R Paket {afex}, dass wir hier flexibel bei der Auswahl der Typen sind. Bei der einfaktoriellen ANOVA war es noch egal, aber hier macht es mit zwei Faktoren dann schon einen Unterschied welche Art der ANOVA du rechnen möchtest. Ich nutze hier den Standard und damit die Type III ANOVA. Auch erhälst du bei der Funktion aov_car() einen Wert für Eta\(^2\) in der Spalte ges wiedergegeben.

R 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 

Die Ergebnisse sind ähnlich wie schon bei den beiden anderen Algorithmen, aber hier macht es dann einen Unterschied welche ANOVA gerechnet wird???

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.013
R 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)) 
Tabelle 32.6— Datensatz zu dem Längen- und Dickenwachstum von Gurken.
tinytable_cvvkd5gjs1behry8boh4
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")
Abbildung 32.14— Scatterplot des Längenwachstums der drei Gurkensorten über vierzehn Wochen. Die gestrichtelten Linien stellen den Median und die durchgezogene Line den Mittelwert der Gruppen dar.
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")
Abbildung 32.15— Scatterplot des Längenwachstums der drei Gurkensorten über vierzehn Wochen. Die gestrichtelten Linien stellen den Median und die durchgezogene Line den Mittelwert der Gruppen dar.
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  NA
Contrasts set to contr.sum for the following variables: group
Anova 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  NA
Contrasts set to contr.sum for the following variables: group
Anova 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.
Olejnik S, Algina J. 2003. Generalized eta and omega squared statistics: measures of effect size for some common research designs. Psychological methods 8: 434.