Mixed–effects modeling
Generalising to the universe with random item and subject
selection
Davide Crepaldi
MoMo Lab, Department of Psychology
University of Milano Bicocca, Italy
www.davidecrepaldi.net
[email protected]
Spring 2013
Davide Crepaldi
1 / 123
Installation
R objects
Working with R objects
Mixed models
Part III
A practical introduction using R
Davide Crepaldi
99 / 148
Installation
R objects
Working with R objects
Mixed models
Outline of part III
Installation
R objects
Vectors
Matrices
Lists and data frame
9 Working with R objects
Graphics
Session management
Sample sessions
Import–export
Sample sessions
10 Mixed models
Preliminary thoughts
Continuous Y
Sample sessions
Binary Y
Sample sessions
7
8
Davide Crepaldi
100 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Vectors
Vectors: x <- c(10.4, 5.6, 3.1, 6.4, 21.7); c(7.3,
4.0, 6.2, 19.3, 13.4) -> y; x = c(10.4, 5.6, 3.1,
6.4, 21.7)
Operations: 5+4; x+y; z <- x+y
Logical and character vectors: z <- x>5.7; z <c("lab1","lab2","lab1","lab2","lab1")
Missing values: z <- c(1:3,NA)
Indexing: x[3]; x[x>5.7]
Davide Crepaldi
101 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Vectors in use
state <- c("tas", "sa", "qld", "nsw", "nsw",
"nt", "wa", "wa", "qld", "vic", "nsw", "vic",
"qld", "qld", "sa", "tas", "sa", "nt", "wa",
"vic", "qld", "nsw", "nsw", "wa", "sa", "act",
"nsw", "vic", "vic", "act")
incomes <- c(60, 49, 40, 61, 64, 60, 59, 54, 62,
69, 70, 42, 56, 61, 61, 61, 58, 51, 48, 65, 49,
49, 41, 48, 52, 46, 59, 46, 58, 43)
sort(incomes)
cut(incomes, breaks = 35+10*(0:7)) -> incomef
Davide Crepaldi
102 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Vectors – practice
Create un vettore con il nome dei vostri cugini; chiamatelo nomi.
Create un secondo vettore con l’età dei vostri cugini; chiamatelo
eta. Calcolate:
l’età media dei vostri cugini
l’età media che avranno tra 8 anni
l’età media dei cugini con più di 30 anni
l’età media dei 3 cugini più vecchi
Davide Crepaldi
103 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Matrices
x <- matrix(1:20, nrow=4, byrow=F); y <array(1:20, dim=c(4,5))
Operations: x+2; x+y; x*y
Indexing: x[2,3]; x[,3]; x[2,3] <- 54
Dimensions: dim(x); nrow(x); ncol(x)
Davide Crepaldi
104 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Lists and data frame
list1 <- list(name="Fred", wife="Mary",
no.children=3, child.ages=c(4,7,9));
list2 <- list(company="robfromtheriches",
ceo="robinhood", staff=c("littlejohn",
"muchthemillersson","willscarlett"))
Indexing: list1[[4]]; list1$wife
Operations: superlist <- c(list1,list2);
length(superlist)
Data frames: accountants <- data.frame(home=state,
loot=incomes, shot=incomef)
Davide Crepaldi
105 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Working with data frames
attach(accountants); detach(accountants);
summary(accountants); head(accountants, n=20)
accountants$state <- as.character(accountants$state);
accountants$state <- as.factor(accountants$state)
table(home); table(home, shot)
tapply(loot, home, mean); tapply(loot, home, sd);
tapply(shot, home, mean)
aggregate(loot, list(home), mean); aggregate(loot,
list(home,shot), mean)
Davide Crepaldi
106 / 148
Installation
R objects
Working with R objects
Mixed models
Vectors
Matrices
Lists and data frame
Data frame – practice
Con i vettori nomi e eta costruite un dataframe.
Controllate di aver creato il dataframe correttamente facendovi dare
il numero di righe e il numero di colonne
Visualizzate a schermo le prime 3 righe
Aggiungete una colonna al dataframe (ramo) in cui scrivete se il
cugino viene da parte di mamma e di papà
Controllate che tutto sia andato per il verso giusto con un comando
appropriato
Calcolate l’età media dei cugini, separando il ramo materno da
quello paterno
Davide Crepaldi
107 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Graphics
plot(loot); hist(loot); plot(density(loot))
Overlapping graphs: h <- hist(loot, freq=F, plot=F);
d <- density(loot); xlimit <- range(h$breaks,
d$x); ylimit <- range(0, h$density, d$y);
hist(loot, freq=FALSE, xlim=xlimit, ylim=ylimit,
main="", xlab="log RT", ylab="", col="lightgrey",
border="darkgrey"); lines(d$x, d$y)
Davide Crepaldi
108 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Graphics
barplot(table(home)); mosaicplot(table(home,shot))
boxplot(loot ∼ home); plot(1:nrow(accountants),
loot)
Devices: jpeg("plot1.jpg", res=200, height=1500,
width=1500); pdf(); dev.off()
Davide Crepaldi
109 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Data frame and graphics – practice
Costruite un dataframe con i dati dei tempi di uscita dal labirinto dei
topini del primo incontro (ricordate che ogni caso statistico – cella della
tabella – deve stare su una riga a sé stante)
Controllate di aver creato il dataframe correttamente facendovi dare
il numero di righe, il numero di colonne e le statistiche descrittive di
tutte le variabili; stampate inoltre le prime 12 righe per controllare
l’allineamento
Controllate la distribuzione incrociata dei topini nei due diversi tipi
di labirinto attraverso table() e mosaicplot()
Calcolate il tempo medio di uscita dai sei labirinti per ciascun topino
Calcolate il tempo medio di uscita dei sei topini per ciascun labirinto
Confrontate graficamente i tempi di uscita nei due tipi di labirinto
Davide Crepaldi
110 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Data frame and graphics – practice
Costruite un dataframe con i dati sui tempi di lettura nell’esperimento di
priming del secondo incontro (ricordate che ogni caso statistico – cella
della tabella – deve stare su una riga a sé stante)
Controllate di aver creato il dataframe correttamente
Controllate la distribuzione di frequenza dei partecipanti, e la loro
distribuzione incrociata con la variabile target e la variabile soa
Confrontate graficamente i tempi di lettura dei partecipanti,
indipendentemente dalle altre variabili
Confrontate graficamente i tempi di lettura dei partecipanti,
separatamente per l’SOA lungo e l’SOA corto
Fate lo stesso tipo di operazione per le parole
Confrontate graficamente i tempi di risposta nell’SOA lungo e
nell’SOA breve; elaborate il grafico; e salvatelo in un file .jpg con un
nome a vostro piacimento
Davide Crepaldi
111 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session management
Objects: ls(); rm(OldStyleAnova)
History: savehistory("TodaySession.RHistory")
Data saving: save.image("TodaySession.RData")
From the outside: source("Commands.R");
sink("OldStyleAnova.txt")
Help: ?mean; ??deviation
Davide Crepaldi
112 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session management – practice
Controllate le variabili che avete nel vostro workspace
Eliminate tutte le variabili tranne i dataframes
Salvate lo storico dei vostri comandi
Salvate il workspace
Davide Crepaldi
113 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session A1
Lavoreremo sul dataframe lexdec, inserito nel pacchetto languageR
Esploriamo il dataframe (quante variabili, di quale tipo, con quali caratteristiche,
quanti casi statistici, . . . )
Controlliamo la distribuzione di PrevType, quella di NativeLanguage, e quella
bivariata di entrambe le variabili
Controlliamo che le parole target in ciascuna delle quattro condizioni create
dall’incrocio di PrevType e NativeLanguage siano bilanciate per Frequency
RT cambia a seconda di PrevType e NativeLanguage?
Elaboriamo il grafico
Come è la distribuzione di RT? Essa cambia a seconda di PrevType e
NativeLanguage?
Elaboriamo il grafico
Davide Crepaldi
114 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session A2
Lavoreremo sul dataframe lexdec, inserito nel pacchetto
languageR
Controlliamo la distribuzione di Complex, quella di Class, e
quella bivariata di entrambe le variabili
Studiamo come varia la distribuzione di DerivEntropy al
variare di Class
RT e DerivEntropy sono relate?
Studiamo come varia la distribuzione di RT al variare di Class
e Complex
Elaboriamo il grafico e salviamolo in un file .pdf
Davide Crepaldi
115 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session A3
Lavoreremo sul dataframe spanishMeta, inserito nel pacchetto languageR
Calcolate quanti testi diversi (TextName) ci sono per ciascun autore (FullName)
nel dataframe
Calcolate la data media di pubblicazione (TextName) e l’arco temporale di
pubblicazione per ciascun autore.
Estraete dal dataframe tutte le pubblicazioni uscite prima del 1980. Qual é la
media dell’anno di pubblicazione per questi testi?
Studiate la distribuzione del numero di parole autore per autore, anche
graficamente. Salvate il grafico in formato .jpg dopo aver settato
appropriatamente i titoli dell’asse x, dell’asse y e il titolo generale, e dopo aver
scritto nel grafico quanti testi sono stati considerati per ciascun autore
Davide Crepaldi
116 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Excel to R
read.table(...); write.table(...)
File format: read.table(‘‘data.txt’’);
read.csv(‘‘data.csv’’)
Column separator: read.table(‘‘data.txt’’,
sep=‘‘,’’), read.table(‘‘data.txt’’, sep=‘‘’’),
read.table(‘‘data.txt’’, sep=‘‘\t’’)
Decimals: read.table(‘‘data.txt’’, dec=‘‘,’’);
read.table(‘‘data.txt’’, dec=‘‘.’’)
Variable names: read.table(‘‘data.txt’’, header=T)
Davide Crepaldi
117 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Excel to R – practice
Esportate il data frame accountants nel formato .txt
usando il punto e virgola come separatore di colonna
Esportate di nuovo il data frame accountants questa volta
usando il formato .csv
Generate una nuova colonna nel data frame dei cugini con
l’età media del ramo cui ciascun cugino appartiene. Quindi
esportate il data frame in un formato a vostro piacimento
Aprite il file in Office e controllate il separatore dei decimali.
Quindi provate a re–importare il data frame in R con
read.table()
Davide Crepaldi
118 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session B1
Importiamo il data frame heid1.txt
Esploriamo il data frame (quante variabili, di che tipo, con
quali caratteristiche, . . . )
Generiamo un nuovo data frame con il tempo di riposta medio
per ciascun soggetto:
Nominiamo le colonne in modo appropriato
Quanti soggetti hanno partecipato all’esperimento?
C’è qualche partecipante che sarebbe bene escludere dalle
analisi?
Davide Crepaldi
119 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session B2
Importiamo il data frame heid2.txt
Generiamo un nuovo data frame con il tempo di riposta medio
per ciascuna parola:
Nominiamo le colonne in modo appropriato
Attacchiamo per ciascuna parola il relativo valore della
frequenza
Esploriamo i dati cosı̀ da capire se il tempo di risposta dipende
dalla frequenza della radice della parola
Sia il tempo di risposta che la frequenza della base sono
trasformati secondo la funzione logaritmo. Studiamo la
distribuzione di queste variabili, e la loro distribuzione in
versione non trasformata.
Davide Crepaldi
120 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session B3
Importiamo il data frame scramble.txt
Esploriamo il data frame (quante variabili, di che tipo, con
quali caratteristiche, . . . )
Aggiustiamo il formato (e le etichette dei livelli) delle singole
variabili nel data frame
Controlliamo che le diverse rotation dell’esperimento non
abbiano introdotto distorsioni nei dati. Facciamo lo stesso per
computer
Confrontiamo graficamente le distribuzioni dei RT nelle diverse
rotazioni e computer
Calcoliamo le medie per condizione sperimentale, rotazione per
rotazione (e computer per computer)
Davide Crepaldi
121 / 148
Installation
R objects
Working with R objects
Mixed models
Graphics
Session management
Import–export
Session B3
Calcoliamo le medie per soggetto, condizione per condizione,
e vediamo se qualche soggetto dovrebbe essere escluso dalle
analisi
Calcoliamo le medie per target, condizione per condizione, e
vediamo se qualche target dovrebbe essere escluso dalle analisi
Sul dataset per soggetto, esploriamo con boxplot e
aggregate la possibilit di un’interazione tra posizione e
morfologia.
Facciamo lo stesso sul dataset per target.
Davide Crepaldi
122 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
When to use mixed effects
There are random effects, typically, multiple individuals tested
over multiple trials
Standard analysis has problems (unbalanced designs, power
issues, obvious non–matched covariates, binomial Y , unclear
F1 and/or F2 results)
If standard analysis is OK, no point in more complex analyses!
Davide Crepaldi
123 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
We will be working on the data frame scramble; we need
functions within languageR
qqmath(∼ rt|factor(sbj.id), data=subset(scramble,
rt>0))
Davide Crepaldi
124 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
Davide Crepaldi
125 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
xylowess.fnc(rt ∼ item.position|factor(sbj.id),
data=subset(scramble, rt>0))
Davide Crepaldi
126 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
Davide Crepaldi
127 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
par(mfrow=c(1,3));
qqnorm(rt[rt>0]);
qqnorm(log(rt[rt>0]));
qqnorm(-1000/rt[rt>0]);
Davide Crepaldi
128 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Preliminary data exploration
Davide Crepaldi
129 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
scramble2 <- subset(scramble, rt>0);
m1 <- lmer(-1000/rt ∼ item.position + (1|sbj.id) +
(1|item.id), data=scramble2);
print(m1, cor=F);
ranef(m1)$sbj.id;
var(ranef(m1)$sbj.id);
head(fitted(m1));
-1.502e+00 + ranef(m1)$sbj.id["1",] +
ranef(m1)$item.id["105",] + -3.222e-04*38;
Davide Crepaldi
130 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
2*(1-pt(abs(-4.03), 2937-2));
pvals.fnc(m1)$fixed;
m2 <- lmer(-1000/rt ∼ item.position + (1 +
item.position|sbj.id) + (1|item.id), data=scramble2);
anova(m1,m2)
print(m2, cor=F);
Davide Crepaldi
131 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
m3 <- lmer(-1000/rt ∼ item.position + rt.preceding
+ acc.preceding + (1 + item.position|sbj.id) +
(1|item.id), data=scramble2);
anova(m3,m2)
print(m3, cor=F);
m4 <- lmer(-1000/rt ∼ item.position + rt.preceding +
acc.preceding + (1 + item.position + rt.preceding +
acc.preceding|sbj.id) + (1|item.id), data=scramble2);
anova(m4,m3)
print(m4, cor=F);
Davide Crepaldi
132 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
m5 <- lmer(-1000/rt ∼ item.position + rt.preceding
+ acc.preceding + X1.posit*X2.morph + (1 +
item.position + rt.preceding + acc.preceding|sbj.id)
+ (1|item.id), data=scramble2);
anova(m5,m4)
print(m5, cor=F);
m6 <- lmer(-1000/rt ∼ item.position + rt.preceding +
acc.preceding + X1.posit*X2.morph + (1 +
item.position + rt.preceding + acc.preceding|sbj.id)
+ (0+X1.posit*X2.morph|sbj.id) + (1|item.id),
data=scramble2);
anova(m6,m5)
print(m5, cor=F);
Davide Crepaldi
133 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
head(fitted(m6));
head(fitted(m1));
head(coef(m6)$sbj.id);
m6.trimmed <- lmer(-1000/rt ∼ item.position +
rt.preceding + acc.preceding + X1.posit*X2.morph + (1
+ item.position + rt.preceding +
acc.preceding|sbj.id) + (0+X1.posit*X2.morph|sbj.id)
+ (1|item.id), data=scramble2, subset =
abs(scale(resid(m6)))<2.5);
print(m6.trimmed, cor=F);
Davide Crepaldi
134 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
jpeg("Trimming.jpg", res=200, height=1500,
width=2000);
par(mfrow=c(1,2));
qqnorm(residuals(m6),main="Base model");
qqline(residuals(m6));
qqnorm(residuals(m6.trimmed),main="Trimmed model");
qqline(residuals(m6.trimmed));
dev.off();
Davide Crepaldi
135 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
Davide Crepaldi
136 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Mixed models
anova(m6.trimmed);
anova(m6.trimmed)$F;
anova(m6.trimmed)$Df;
length(fitted(m6.trimmed));
sum(nrow(anova(m6.trimmed)));
df(anova(m6.trimmed)$F, anova(m6.trimmed)$Df,
length(fitted(m6.trimmed))-sum(nrow(anova(m6.trimmed))));
Davide Crepaldi
137 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session C1
Lavoreremo sul dataframe lexdec, inserito nel pacchetto languageR
Rappresentiamo graficamente la possibilità che i partecipanti mostrino effetti di
affaticamento o di apprendimento del compito
Costruite il primo modello con l’effetto fisso della posizione del trial nella lista
sperimentale (Trial), e le random intercept di Subject e Word. Stampate il
modello. Quale effetto random spiega più varianza?
Calcolate la varianza della random intercept di Word usando la funzione ranef
Aggiungete al modello una random slope di Trial per i soggetti. L’aggiunta di
questo nuovo parametro è giustificata?
Aggiungete quindi al modello anche l’effetto fisso di NativeLanguage.
L’aggiunta di questo nuovo parametro è giustificata? E quella di Frequency? E
quella di FamilySize? E quella di DerivEntropy?
Una volta stabilito il modello migliore, calcolate la significatività dei fattori fissi
inclusi sia con la funzione pvals.fnc che attraverso il metodo analitico basato
su pt.
Davide Crepaldi
138 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session C2, part I
Lavoreremo sul dataframe FakeIrregulars, da importare dal file .txt
Trasformiamo in fattori le variabili item.id, target.id, sbj.id e
accuracy
Rappresentiamo graficamente la possibilità che i partecipanti mostrino
effetti di affaticamento o di apprendimento del compito
Usando la funzione qqnorm, decidere se usare i tempi di reazione grezzi
(rt) nelle analisi, oppure la loro trasformaziona logaritmica, oppure
ancora la loro trasformazione inversa.
Costruite il primo modello inserendo l’effetto fisso di tutte le covariate
(target.length, target.wf, target.sf e target.n), e le random
intercept di sbj.id e target.it. Stampate il modello.
Davide Crepaldi
139 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session C2, part II
Eliminate via via, uno per volta, gli effetti fissi meno significativi (con il t
più basso). Ad ogni passaggio controllate che non ci sia una perdita
significativa di goodness of fit attraverso la funzione anova; al primo
passaggio in cui anova dà un risultati significativo, fermatevi: quello é il
modello migliore.
Aggiungete al modello una random slope di target.wf per i soggetti.
L’aggiunta di questo nuovo parametro è giustificata?
Calcolate la significatività dei fattori fissi con la funzione pvals.fnc.
Calcolate la correlazione tra target.wf e target.sf. Rifittate quindi il
modello usando invece che target.sf il suo rapporto con target.wf; la
sintassi da usare I(target.sf/target.wf). Confrontate quindi i
parametri, e in particolare i t value, di questo nuovo modello con quello
in cui invece target.sf era usato senza alcuna trasformazione.
Davide Crepaldi
140 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Generalised mixed models
We will be working on the data frame scramble;
m1 <- lmer(acc ∼ item.position + (1|sbj.id) +
(1|item.id), data=scramble, family="binomial");
print(m1, cor=F);
ranef(m1)$sbj.id;
head(fitted(m1));
3.503162 + ranef(m1)$sbj.id["1",] +
ranef(m1)$item.id["105",] + .001787*38;
exp( 3.503162 + ranef(m1)$sbj.id["1",] +
ranef(m1)$item.id["105",] + .001787*38 );
Davide Crepaldi
141 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Generalised mixed models
m1 <- lmer(acc ∼ acc.preceding + rt.preceding +
(1|sbj.id) + (1|item.id), data=scramble,
family="binomial");
print(m1, cor=F);
m1 <- lmer(acc ∼ acc.preceding + rt.preceding + (1 +
acc.preceding + rt.preceding|sbj.id) + (1|item.id),
data=scramble, family="binomial");
print(m1, cor=F);
Davide Crepaldi
142 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Generalised mixed models
m2 <- lmer(acc ∼ rt.preceding + (1 +
rt.preceding|sbj.id) + (1|item.id), data=scramble,
family="binomial");
anova(m1,m2);
m2 <- lmer(acc ∼ acc.preceding + rt.preceding +
X1.posit*X2.morph + (1 + acc.preceding +
rt.preceding|sbj.id) + (1|item.id), data=scramble,
family="binomial");
anova(m1,m2);
Davide Crepaldi
143 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session D1, part I
Lavoreremo sul data frame dative, inserito nel pacchetto languageR
Controllate la distribuzione della variabile Verb
Adottiamo un approccio backward: costruiamo quindi il primo modello
con tutti i potenziali predittori, togliendo poi uno alla volta, partendo dal
meno significativo, tutti quelli che non sono giustificati nel modello. La
variabile dipendente é RealizationOfRecipient.
Gli effetti fissi sono AccessOfTheme, AccessOfRec, LengthOfRecipient,
AnimacyOfRec, AnimacyOfTheme, PronomOfTheme, DefinOfTheme,
LengthOfTheme, SemanticClass, and Modality.
Un solo effetto random, l’intercetta di Verb.
Davide Crepaldi
144 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session D1, part II
Notate quanti parametri (cioé quante righe nella tabella degli effetti fissi)
pertengano a ciascun predittore. Cercate di capire come mai il numero di
questi parametri possa variare, e a che cosa sia uguale. Un suggerimento:
ha qualcosa a che fare con length(levels(...)).
Estraete i primi 6 valori predetti con head(fitted(...)). Calcolate
analiticamente il primo, usando ranef e i parametri del modelli estratti
con print.
Aggiungere una random slope per LengthOfRecipient e controllate se
sia giustificata la sua presenza nel modello con la funzione anova.
Fate lo stesso per AnimacyOfTheme, AccessOfRec e Modality.
Davide Crepaldi
145 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session D2, part I
Lavoreremo sul dataframe FakeIrregulars, da importare dal file .txt
Trasformiamo in fattori le variabili item.id, target.id, sbj.id e
accuracy
Studiamo la distribuzione di risposte corrette ed errate a seconda di
prime.type, morph e del loro incrocio.
Costruiamo un primo modello con morph, prime.type e la loro
interazione come effetti fissi – la sintassi morph*prime.type specifica in
un solo colpo questi tre effetti –, una random intercept per sbj.id e una
random intercept per target.id. Stampiamo il modello e cerchiamo di
intuire quali sono i parametri dell’interazione.
Davide Crepaldi
146 / 148
Installation
R objects
Working with R objects
Mixed models
Preliminary thoughts
Continuous Y
Binary Y
Session D2, part II
Proviamo a togliere la random intercept per target.id e controlliamo se
la goodness of fit cambia significativamente. Facciamo lo stesso con la
random intercept per sbj.id.
Aggiungiamo ora, una per volta, le random slope per sbj.id di morph e
prime.type. La loro presenza é giustificata nel modello? E migliora la
significatività dei parametri corrispondenti tra gli effetti fissi?
Aggiungiamo ora, una per volta, gli effetti delle variabili concomitanti
target.wf, target.length e target.sf, ciascuna in interazione con le
due variabili già presenti tra gli effetti fissi del modello – la sintassi
corretta é morph*prime.type*target.wf. La loro presenza é giustificata
nel modello? E migliora la significatività dei parametri delle variabili
morph e prime.type?
Davide Crepaldi
147 / 148
Contacts
Contacts
www.davidecrepaldi.net
[email protected]
Davide Crepaldi
148 / 148
Scarica

Mixed–effects modeling - Generalising to the universe with random