Dear colleagues,
I have a data from a repeated measures design that I'm analysing through
a mixed model. Nine independent sampling units (flasks with culture medium
with algae) were randomly divided into 3 groups ("c", "t1", "t2"). There is
no need for inclusion of the random effect of the intercept, because the
nine sample units are homogeneous among each other (samples taken from the
same culture). The algal concentration was measured every two days for 10
days. The goal is to test differences between treatments.
I estimated a model with only the intercept and the interaction Group * Day
to test which is the best:
library(nlme)
library(lme4)
Day = rep(c(0,2,4,6,8,10),each=9)
Group = rep(c("c","c","c","t1","t1","t1","t2","t2","t2"),6)
Individual = rep(1:9,6)
X = c(0.71,0.72,0.71,0.72,0.72,0.72,0.70,0.69,0.70,0.72,0.72,
0.71,0.72,0.72,0.71,0.71,0.70,0.71,0.73,0.73,0.69,0.74,
0.69,0.73,0.67,0.71,0.69,0.71,0.71,0.72,0.70,0.71,0.70,
0.52,0.64,0.60,0.70,0.73,0.73,0.67,0.66,0.71,0.47,0.56,
0.54,0.65,0.73,0.73,0.67,0.71,0.58,0.44,0.52,0.58)
xyplot(X~Day, groups=Group)
LME = lme(X ~ 1, random = ~Day|Individual)
Erro em lme.formula(X ~ 1, random = ~Day | Individual) :
nlminb problem, convergence error code = 1
message = iteration limit reached without convergence (10)
LME1 = lme(X ~ Group*Day, random = ~Day|Individual)
Erro em lme.formula(X ~ Group * Day, random = ~Day | Individual) :
nlminb problem, convergence error code = 1
message = iteration limit reached without convergence (10)
LMER = lmer(X ~ 1 + (Day|Individual))
LMER1 = lmer(X ~ Group*Day + (Day|Individual))
AIC(LMER)
[1] -179.0302
AIC(LMER1)
[1] -151.1938
anova(LMER,LMER1)
Data:
Models:
LMER: X ~ 1 + (Day | Individual)
LMER1: X ~ Group * Day + (Day | Individual)
Df AIC BIC logLik Chisq Chi Df Pr(>Chisq)
LMER 5 -187.2 -177.26 98.602
LMER1 10 -203.2 -183.31 111.600 25.996 5 8.939e-05 ***
xyplot(fitted(LMER)~Day, groups=Group)
xyplot(fitted(LMER1)~Day, groups=Group)
1st question: Why function lme4:lmer converge, but the nlme:lme
doesn't? The first is better than the second?
2nd question: Why does the "anova" give distinct values of AIC for the
two models. If we look the AIC value of each model, the best model is LMER,
but the "anova" says that LMER1 is the best.
3rd question: Why the fitted values of the model with only the
intercept (LMER) vary over time?
Thank you very much for any help
Diego PJ
[[alternative HTML version deleted]]