[R] Need advises on mixed-effect model ( a concrete example)
wphantomfr
wphantomfr at gmail.com
Thu Feb 3 20:13:59 CET 2011
Dear R-help members,
I'm trying to run LME model on some behavioral data and need
confirmations about what I'm doing...
Here's the story...
I have some behavioral reaction time (RT) data (participants have to
detect dome kind of auditory stimuli). the dependant variable is RT
measured in milliseconds. 61 participants were tested separated in 4 age
groups (unblanced groups, factor GROUP). Each participant (SUBJECT) was
tested in each of the 3 experimental condition (repeated measures,
factor COND).
Futhermore. In each condition there were 16 trials (and thus 16 measures
of RT) with 16 different stimuli. However some trials were rejected from
analysis because the participant did not detect the stimulus or because
RT were outliers (extremely high or low).
I also have a numeric acoustic parameter with 2 possible values (A1
and A2) for each stimuli that I may be interested to take into account.
Since the design was unbalanced, and after reading a lot here, I
decided to try to use linear mixed effect. I thus have the Pinheiro &
Bates. I have to admit that, as a psychologist, I did not understand
everything in the book. However I try to applied what I have read to my
data.
I wanted a confirmation on what I am doing.
I see two ways using my data with lme models depending on the granular
level of my dependant variable that give quite different results.
I tend to prefer the first approach for some reason but I need
confirmation the correctness of both approach.
My questions :
1) are both possibility correctly implemented in my analysis ?
2) Is one way better than the other ?
3) Is there better ways of testing the effects of group, condition and
acoustQ ?
4) If one of my factor's effect is significant, how can I test correctly
the difference between the modalities of this factor (since it seems
that the t-values in the model summary are not a good way of doing
this). Is there a specific post-hoc procedure ?
Thanks in advance for any help
Sylvain Clément
University of Lille Nord de France
######################################
1st APPROACH Using Mean RT of participants as my DV.
######################################
# first fit a model with single fixed effects of condition and groups
result.lme1<-lme(RT~CONDITION+GROUP,data= meandata,random=~1|SUBJECT,
method="ML")
plot(effet.lme) # residuals are goodlooking
#then I try to add an interaction term in the model
result.lme2<-lme(RT~CONDITION*GROUP,data= meandata,random=~1|SUBJECT,
method="ML")
#I then compare both models
anova(result.lme1,result.lme2)
#gives
Model df AIC BIC logLik Test L.Ratio p-value
effet.lme 1 8 3696.878 3726.822 -1840.439
effet.lme2 2 14 3705.821 3758.223 -1838.910 1 vs 2 3.056792 0.8017
# Thus I only use the 1st model without interaction term
# try another model adding an acoustic parameter as a fixed effect
result.lme3<-lme(RT~CONDITION+GROUP+ACOUSTQ,data=
meandata,random=~1|SUBJECT, method="ML")
anova(result.lme1,result.lme3) # gives :
Model df AIC BIC logLik Test L.Ratio p-value
effet.lme 1 8 3696.878 3726.822 -1840.439
effet.lme3 2 9 3668.755 3702.442 -1825.377 1 vs 2 30.12318 <.0001
# thus it seems interesting to add my accoustic parameter...
# I refit my 3rd model using the REML method and look at the results
result.lme3<-lme(RT~CONDITION+GROUP+ACOUSTQ,data=
meandata,random=~1|SUBJECT, method="REML")
summary(result.lme3)
Linear mixed-effects model fit by maximum likelihood
Data: meandata
AIC BIC logLik
3668.755 3702.442 -1825.377
Random effects:
Formula: ~1 | SUBJECT
(Intercept) Residual
StdDev: 54.41472 74.60438
Fixed effects: TR ~ CONDITION + GROUP + ACOUSTQ
Value Std.Error DF t-value p-value
(Intercept) 564.6296 11.908961 257 47.41216 0.0000
CONDITIONCond2 -23.1851 10.463814 257 -2.21574 0.0276
CONDITIONCond3 2.6871 10.463814 257 0.25680 0.7975
GROUP.L 119.2404 20.725884 48 5.75321 0.0000
GROUP.Q 43.7518 18.663080 48 2.34430 0.0232
GROUP.C -9.8656 16.341936 48 -0.60370 0.5489
ACOUSTQA2 -47.7384 8.543669 257 -5.58758 0.0000
Correlation:
(Intr) CONDITIONCond2CONDITIONCond1GROUP.L
GROUP.Q GROUP.C
CONDITIONCond2 -0.439
CONDITIONCond1 -0.439 0.500
GROUP.L 0.178 0.000 0.000
GROUP.Q 0.228 0.000 0.000 0.184
GROUP.C 0.004 0.000 0.000 0.180 0.169
ACOUSTQA2 -0.359 0.000 0.000 0.000
0.000 0.000
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.54282265 -0.59821033 -0.04686768 0.52178476 2.84475825
> anova(result.lme3)
numDF denDF F-value p-value
(Intercept) 1 257 3498.963 <.0001
CONDITION 2 257 3.696 0.0261
GROUP 3 48 12.769 <.0001
ACOUSTQ 1 257 31.221 <.0001
This gives me an significant effect my 3 fixed factor
######################################
2ND APPROACH : using RT at single trial level
######################################
> result.lme1<-lme(TR~CONDITION+GRPEAGE,data=
TOUTESCOND,random=~1|SUJET,method="ML")
> plot(result.lme1) #residuals OK
result.lme2<-lme(TR~CONDITION*GRPEAGE,data=
TOUTESCOND,random=~1|SUJET,method="ML")
> anova(result.lme1,result.lme2)
Model df AIC BIC logLik Test L.Ratio p-value
result.lme1 1 8 27464.80 27509.88 -13724.40
result.lme2 2 14 27472.97 27551.87 -13722.49 1 vs 2 3.828608 0.6999
> result.lme3<-lme(TR~CONDITION+GRPEAGE+ACOUSTQ,data=
TOUTESCOND,random=~1|SUJET,method="ML")
> anova(result.lme1,result.lme3)
Model df AIC BIC logLik Test L.Ratio p-value
result.lme1 1 8 27464.8 27509.88 -13724.40
result.lme3 2 9 27452.3 27503.02 -13717.15 1 vs 2 14.49784 1e-04
> result.lme3<-lme(TR~CONDITION+GRPEAGE+ACOUSTQ,data=
TOUTESCOND,random=~1|SUJET,method="REML")
Linear mixed-effects model fit by REML
Data: TOUTESCOND
AIC BIC logLik
27403.99 27454.68 -13693.00
Random effects:
Formula: ~1 | SUJET
(Intercept) Residual
StdDev: 60.8684 177.9567
Fixed effects: TR ~ CONDITION + GROUP + ACOUSTQ
Value Std.Error DF t-value p-value
(Intercept) 627.3987 19.596021 2001 32.01664 0.0000
CONDITIONCond2 -13.0783 9.575935 2001 -1.36574 0.1722
CONDITIONCond3 6.8634 9.637296 2001 0.71217 0.4764
GROUPE8ANS -58.5550 24.369628 63 -2.40279 0.0192
GROUPE9ANS -121.7949 23.811142 63 -5.11504 0.0000
GROUPE10ANS -138.3129 26.486620 63 -5.22199 0.0000
ACOUSTQ -30.3613 7.968108 2001 -3.81036 0.0001
Correlation:
(Intr) CONDITIONCnd2 CONDITIONCnd3 GROUP8
GROUP9 GROUP1
CONDITIONCond2 -0.251
CONDITIONCond3 -0.244 0.502
GROUP8ANS -0.717 0.004 0.003
GROUP9ANS -0.733 0.002 -0.001 0.590
GROUP10ANS -0.657 -0.004 -0.008 0.531 0.543
ACOUSTQ -0.166 0.018 0.006 -0.007
-0.006 -0.005
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.6253233 -0.6573193 -0.0620116 0.5675166 4.9655693
Number of Observations: 2071
Number of Groups: 67
# Finally
> anova(result.lme3)
numDF denDF F-value p-value
(Intercept) 1 2001 3979.373 <.0001
CONDITION 2 2001 2.079 0.1254
GROUP 3 63 12.648 <.0001
ACOUSTQ 1 2001 14.519 0.0001
This gives me no significant effect of CONDITION in this case.
More information about the R-help
mailing list