Mixed models Jonathan Harrington Die RBefehle mixed txt

Mixed models Jonathan Harrington Die R-Befehle: mixed. txt pfad = "Das Verzeichnis, wo die Daten gespeichert sind" attach(paste(pfad, "anova 1"), sep="") library(car) library(language. R) library(multcomp) alc = read. table(paste(pfad, "alcdata. txt", sep="")) soa = read. table(paste(pfad, "soa. txt", sep="")) vot = read. table(paste(pfad, "vot. txt", sep=""))

Mixed models Baayen, R. H. (in press) Analyzing Linguistic Data: A practical introduction to Statistics. Kapitel 7 http: //www. ualberta. ca/~baayen/publications/baayen. CUPstats. pdf Artikel in einem Special Issue im Journal of Memory and Language, Vol. 59. insbesondere: Baayen, Davidson & Bates (2008); Quene & van den Bergh (2008); Jaeger (2008). 2 Präsentationen hier vorhanden http: //hlplab. wordpress. com/2009/05/03/multilevel-model-tutorial/ Levy & Jaeger (2009) A Brief and Friendly Introduction to Mixed-Effects Models in Psycholinguistics Frank & Jaeger (April, 2009) Post hoc comparisons Additional Issues: Random effects diagnostics, multiple comparisons Erste Veröffentlichung: Pinheiro & Bates (2000). http: //www. amazon. com/Mixed-Effects-Models-S-S-Plus/dp/0387989579

Subject und Items als Random Factors 3 Vpn ( produzierten 3 Wörter (Bart, Pfad, Start) jeweils 2 Mal: einmal phrasenmedial, einmal phrasenfinal. Unterscheiden sich phrasenmedial und –final in F 1? soa: eine modifizierte Datei von Baayen et al (2008) (selbe Werte, andere Namen) F 1 (F 1 -Werte im Vokal); Vpn (s 1, s 2, s 3), W (Bart, Pfad, Start), Pfinal (medial, final) Faktor W: between/within: within Faktor Pfinal: between/within: within Zuerst eine Abbildung

Zuerst eine Abbildung attach(soa) boxplot(F 1 ~ Pfinal * W) Signifikanter Unterschied in Pfinal (zwischen long und short)? Wort * Pfinal Interaktion?

soa. t = Anova. prepare(soa, c("s", "w", "d")) soa. lm = lm(soa. t$d ~ 1) Anova(soa. lm, idata=soa. t$w, idesign = ~ W * Pfinal) Type III Repeated Measures MANOVA Tests: Pillai test statistic Df test stat approx F num Df den Df Pr(>F) (Intercept) 1 1. 00 1509. 78 1 2 0. 0006617 *** W 1 0. 98 26. 67 2 1 0. 1356518 Pfinal 1 0. 78 7. 24 1 2 0. 1147994 W: Pfinal 1 0. 62 0. 81 2 1 0. 6186755

Subject und Items als Random Factors Jedoch hat Pfinal eindeutig denselben Einfluss pro Wort. d. h. wir sollten hier die Variation zwischen den Wörtern ausklammern. (Dass der Mittelwert von Bart kleiner ist im Vgl. zu Pfad oder Start ist für diese Fragestellung uninteressant genauso uninteressant wie die Variation zwischen Sprechern).

Subject und Items als Random Factors Wir wollen also 2 Sorten von Variation gleichzeitig ausklammern: wegen der Sprecher (Subject as a random factor) wegen der Wörter (Item as a random factor) Dadurch lösen wir gleichzeitig ein großes Problem in der Statistik (Clark, 1973): 'language as fixed effect fallacy'.

Subject und Items als Random Factors soa. t = Anova. prepare(soa, c("s", "w", "d")) soa. lm = lm(soa. t$d ~ 1) Anova(soa. lm, idata=soa. t$w, idesign = ~ W * Pfinal) Random: Subject bedeutet: signifikante Ergebnisse sind nicht nur für diese Vpn sondern allgemein für ähnliche Sprecher gültig. Fixed: W, Pfinal (beide within) signifikante Ergebnisse in Pfinal sind nur für W gültig (Bart, Pfad, Start): d. h. damit die Ergebnisse allgemein für ähnliche Wörter gültig sind, müsste noch einen RM-Manova durchgeführt werden mit Wort als Random Faktor.

Subject und Items als Random Factors ein by-subject RM-Manova (Subject als Random Faktor) ein by-item RM-Manova (Wort als Random Faktor) Die Ergebnisse werden kombiniert in einer Statistik min. F' Das Verfahren ist schrecklich (siehe Johnson, 2008) Keine solchen Komplikationen in einem Mixed-Model

Weitere Vorteile von einem Mixed-Model Es muss nicht gemittelt werden pro Stufen-Kombinationen Die Zellen müssen nicht vollständig pro Vpn sein. Sprache engl. oder span. between within Vpn Sprechtempo Vokal lang. i e a w 1. 1 w 2 w 3 w 1. 2 Mittelwert w 1. 3. . . w 1. 10 schnell i e a w 4 w 5 w 6

Nachteile von einem Mixed-Model library(lme 4) und mixed modelling (MM) überhaupt sind noch in der Entwicklungsphase. Daher bugs, häufige code Änderungen und einige Teile des Verfahrens sind in R noch nicht ganz vollständig. Siehe: http: //wiki. r-project. org/rwiki/doku. php? id=guides: lmer-tests Mit MM können zwar Werte aus der t- und F-Verteilung berechnet werden, aber diese lassen sich nur schwierig und vielleicht sogar ungenau in Wahrscheinlichkeiten umsetzen – weil die Freiheitsgrade nicht eindeutig berechnet werden können.

lmer() ist die Funktion für mixed modelling. berechnet ein lineares Modell, in dem der Abstand zwischen eingeschätzen und tatsächlichen Werten minimiert wird minimiert den Abstand durch REML (restricted maximum likelihood) muss mindestens einen Random-Faktor enthalten vereinheitlicht RM-Anova und logistische Regression 3. by-subject intercept adjustment 1. Neigung Eingeschätzte Werte ^y = bx + k + e 2. Intercept Vpn + e. W 4. by-item intercept adjustment x ist ein Faktor-Code* (z. B 0 oder 1 für 2 Stufen) *siehe: http: //www. ats. ucla. edu/stat/r/library/contrast_coding. htm

lmer() abhängige Variable Vpn und W als Random Factors unabhängige Variable(n) 1. by-subject and by-item intercept adjustment F 1. lmer = lmer(F 1 ~ Pfinal + (1 | Vpn) + (1 | W)) = Sprecher- und Wortvariationen werden für Pfinal (ohne zwischen den Stufen zu differenzieren) ausgeklammert. 2. by-subject and by-item intercept and slope adjustment F 1. lmer = lmer(F 1 ~ Pfinal + (1 + Pfinal | Vpn) + (1 +Pfinal | W)) = Sprecher- und Wortvariationen werden getrennt für die Stufen (long, short) von Pfinal berechnet und ausgeklammert.

lmer() ^y = bx + k + e print(F 1. lmer, corr=F) Linear mixed model fit by REML Fixed effects: Estimate Std. Error t value (Intercept) 522. 111 19. 604 26. 633 Pfinalshort -18. 889 5. 525 -3. 419 + e. W Fixed effects (bezieht sich daher auf Pfinal) Random effects (bezieht sich daher auf Vpn und W) ranef(F 1. lmer) $Vpn (Intercept) s 1 -20. 557646 s 2 22. 948070 s 3 -2. 390424 Vpn $W (Intercept) Bart -27. 94979 Pfad 14. 13553 Start 13. 81427 fitted(F 1. lmer) [1] 473. 6037 515. 6890 515. 3677 454. 7148. . . auch BLUPS genannt (best linear unbiased predictor). Ein BLUP pro Vpn und pro Wort summieren auf 0. Eingeschätzte Werte

lmer() Fixed $Vpn (Intercept) s 1 -20. 557646 $W (Intercept) Bart -27. 94979 Linear mixed model fit by REML Fixed effects: Estimate Std. Error t value (Intercept) 522. 111 19. 604 26. 633 Pfinalshort -18. 889 5. 525 -3. 419 ^y = bx + k + e Random Vpn + e. W Eingeschätzer Wert für phrasenmedialer (Stufe, short) Bart, Vpn. s 1 -18. 889 * 1 + 522. 111 -20. 557646 -27. 94979 1 [1] 454. 7146 fitted(F 1. lmer)[4] [1] 454. 7148 contrast(Pfinal) short long 0 short 1

lmer(), t-Werte und Basis-Kodierung Linear mixed model fit by REML Fixed effects: Estimate Std. Error t value (Intercept) 522. 111 19. 604 26. 633 Pfinalshort -18. 889 5. 525 -3. 419 Die Stufe short vom Faktor Pfinal ist 3. 419 Standard. Abweichung von der Basis-Stufe desselben Faktors entfernt. Basis-Stufe levels(Pfinal) [1] "long" "short" (Also: der absolute Abstand zwischen long und short ist 3. 419 Standard-Abweichungen). Die Basis-Stufe kann man mit relevel() auswählen.

lmer(), t-Werte und Basis-Kodierung Linear mixed model fit by REML Fixed effects: Estimate Std. Error t value (Intercept) 522. 111 19. 604 26. 633 Pfinalshort -18. 889 5. 525 -3. 419 Die Basis-Stufe kann man mit relevel() auswählen. Pfinal 2 = relevel(Pfinal, "short") F 1 b. lmer = lmer(F 1 ~ Pfinal 2 + (1 | Vpn) + (1 | W)) print(F 1 b. lmer, corr=F) Fixed effects: Estimate Std. Error t value (Intercept) 503. 222 19. 604 25. 670 Pfinal 2 long 18. 889 5. 525 3. 419 Kein Freiheitsgrad, keine p-Werte

lmer(), t-Werte, p-Werte Linear mixed model fit by REML Fixed effects: Estimate Std. Error t value (Intercept) 522. 111 19. 604 26. 633 Pfinalshort -18. 889 5. 525 -3. 419 Die Wahrscheinlichkeit könnte mit der höchst möglichen Anzahl der Freiheitsgrade von der t-Verteilung berechnet werden. Freiheistgradanzahl = Anzahl der Werte - (Anzahl der Stufen) = 16 2 * ( 1 - pt(3. 419, 16)) [1] 0. 003516309 Aber der p-Wert ist anti-konservativ (ggf. zu niedrig) und daher nicht zuverlässig.

lmer() und MCMC sampling Um genauere Wahrscheinlichkeiten der im lmer() entstandenen t-Werte zu berechnen, gibt es Markov Chain Monte Carlo sampling (MCMC). kann z. Zt. für solche Modelle angewandt werden F 1. lmer = lmer(F 1 ~ Pfinal + (1 | Vpn) + (1 | W)) nicht für diese F 1. lmer = lmer(F 1 ~ Pfinal + (1 + Pfinal | Vpn) + (1 +Pfinal | W))

lmer(), t-Werte und Basis-Kodierung Bitte beachten: es gibt einen Bug in pvals. fnc(), sodass die Werte von F 1. lmer (!!) nach der Durchführung der Funktion geändert werden. Daher bitte immer gleich nach pvals. fnc() nochmal die lmer() Funktion duchführen! F 1. lmer = lmer(F 1 ~ Pfinal + (1 | Vpn) + (1 | W)) F 1. fnc = pvals. fnc(F 1. lmer) F 1. lmer = lmer(F 1 ~ Pfinal + (1 | Vpn) + (1 | W))

lmer() und MCMC sampling F 1. fnc = pvals. fnc(F 1. lmer) F 1. fnc$fixed Highest Posterior Density Estimate MCMCmean HPD 95 lower HPD 95 upper p. MCMC Pr(>|t|) (Intercept) 522. 11 522. 13 491. 88 554. 614 0. 0001 0. 0000 Pfinalshort -18. 89 -18. 88 -35. 94 -1. 013 0. 0356 0. 0035 2 * ( 1 - pt(3. 419, 16)) Die Daten wurden mit einem Mixed Model mit Subject und Word als Random Faktoren und phrasenfinale Längung als unabhängige Variable analysiert. Alle Wahrscheinlichkeiten wurden mit Markov. Chain-Monte-Carlo sampling (MCMC) berechnet (Baayen et al, 2008). Die phrasenfinale Längung hatte einen signifikanten Einfluss auf F 1 (t = 3. 42, MCMC: p < 0. 05).