Dear Benedetta, I think you might want (1+T+Z|subject) as random effects rather than (1+T|subject) + (1 + Z|subject). The latter has two random intercepts per subject: a recipe for disaster.
Follow-up posts should only go to the mixed models mailing list which I'm cc'ing. Best regards, ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance Kliniekstraat 25 1070 Anderlecht Belgium + 32 2 525 02 51 + 32 54 43 61 85 thierry.onkel...@inbo.be www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -----Oorspronkelijk bericht----- Van: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] Namens Benedetta Cesqui Verzonden: maandag 25 november 2013 11:13 Aan: r-help@r-project.org Onderwerp: [R] lmer specification for random effects: contradictory reults Hi All, I was wondering if someone could help me to solve this issue with lmer. In order to understand the best mixed effects model to fit my data, I compared the following options according to the procedures specified in many papers (i.e. Baayen <http://www.google.it/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&ved=0CDsQFjAA &url=http%3A%2F%2Fwww.ualberta.ca%2F~baayen%2Fpublications%2FbaayenDavidsonB ates.pdf&ei=FhqTUoXuJKKV7Abds4GYBA&usg=AFQjCNFst7GT7mBX7w9lXItJTtELJSKWJg&si g2=KGA5MHxOvEGwDxf-Gcqi6g&bvm> R.H. et al 2008) Here, dT_purs is the response variable, T and Z are the fixed effects, and subject is the random effect. Random and fixed effects are crossed.: mod0 <- lmer(dT_purs ~ T + Z + (1|subject), data = x) mod1 <- lmer(dT_purs ~ T + Z + (1 +tempo| subject), data = x) mod2 <- lmer(dT_purs ~ T + Z + (1 +tempo| subject) + (1+ Z| subject), data = x) mod3 <- lmer(dT_purs ~ T * Z + (1 +tempo| subject) + (1+ Z| subject), data = x) mod4 <- lmer(dT_purs ~ T * Z + (1| subject), data = x) anova(mod0, mod1,mod2, mod3, mod4) Data: x Models: mod0: dT_purs ~ T + Z + (1 | subject) mod4: dT_purs ~ T * Z + (1 | subject ) mod1: dT_purs ~ T + Z + (1 + T| subject) mod2: dT_purs ~ T + Z + (1 + T| subject ) + (1 + Z | subject) mod3: dT_purs ~ T * Z + (1 + T| subject) + (1 + Z | subject) Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq) mod0 5 -689.81 -669.46 349.91 -699.81 mod4 6 -689.57 -665.14 350.78 -701.57 1.7532 1 0.185473 mod1 7 -689.12 -660.62 351.56 -703.12 1.5504 1 0.213070 mod2 10 -695.67 -654.97 357.84 -715.67 12.5563 3 0.005701 ** mod3 11 -695.83 -651.05 358.92 -717.83 2.1580 1 0.141825 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 It turns out that mod2 has the right level of complexity for this dataset. However when I looked at its summary, I got a correlation of -0.87 for the random effects relative to the T effect and -1 for the random effects relatively to the Z. summary(mod2) Linear mixed model fit by maximum likelihood ['lmerMod'] Formula: dT_purs ~T + Z + (1 + T | subject) + (1 + Z | subject) Data: x AIC BIC logLik deviance -695.6729 -654.9655 357.8364 -715.6729 Random effects: Groups Name Variance Std.Dev. Corr subject (Intercept) 0.0032063 0.05662 T 0.0117204 0.10826 -0.87 subject.1 (Intercept) 0.0005673 0.02382 Z 0.0025859 0.05085 1.00 Residual 0.0104551 0.10225 Number of obs: 433, groups: soggetto, 7 Fixed effects: Estimate Std. Error t value (Intercept) 0.02489 0.03833 0.650 T 0.52010 0.05905 8.808 Z -0.09019 0.02199 -4.101 Correlation of Fixed Effects: (Intr) tempo T -0.901 Z 0.218 -0.026 If I understand correctly what the correlation parameters reported in the table are, the correlation of 1 means that, for the Z effects the random intercept is perfectly collinear with the random slope. Thus, we fit the wrong model. A random intercept only model would have sufficed. Am I correct? If so, should I take mod1 (mod1 <- dT_purs ~ T + Z + (1 + T | subject ) instead of mod2 to fit my data? Why are these results contradictory? Finally is a correlation value of -0.87 a too high or an acceptable value ? Thanks for help me in advance! Best Benedetta --- Benedetta Cesqui, Ph.D. Laboratory of Neuromotor Physiology IRCCS Fondazione Santa Lucia Via Ardeatina 306 00179 Rome, Italy tel: (+39) 06-51501485 fax:(+39) 06-51501482 E_mail: b.ces...@hsantalucia.it [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. * * * * * * * * * * * * * D I S C L A I M E R * * * * * * * * * * * * * Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document. ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
