-----Original Message-----
From: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf
Of David Afshartous
Sent: Friday, August 08, 2008 11:06 AM
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] Different versions of lme4 and covariance
of randomeffects
All,
I recently re-estimated a model after upgrading to Rv2.7.1
that had been estimated Rv2.6.2 previously and was surprised
to see the estimated correlation between random effects in
the model change from 0.3 to 1.0.
It appears that the reason has more to do with the version of
lme4 since when I download the latest possible lme4 to
Rv2.6.2 the results agree.
Below is a reproducible example where the difference in
results is not as dramatic. Of course, for my data the
initial results with the correlation of .3 seem a lot more
plausible than that of 1.0 so I'm inclined to trust those
results more than the newer ones. It seems strange to get
the correlation of 1.
Cheers,
David
library("lme4")
set.seed(500)
n.timepoints <- 4 ## change for shorter examples
n.subj.per.tx <- 20 sd.d <- 5; sd.p <- 2; sd.res <- 1.3 drug
<- factor(rep(c("D", "P"), each = n.timepoints, times =
n.subj.per.tx))
drug.baseline <- rep( c(0,5), each=n.timepoints,
times=n.subj.per.tx ) Patient <- rep(1:(n.subj.per.tx*2),
each = n.timepoints) Patient.baseline <- rep( rnorm(
n.subj.per.tx*2, sd=c(sd.d, sd.p) ), each=n.timepoints ) time
<- factor(paste("Time-", rep(1:n.timepoints, n.subj.per.tx*2),
sep=""))
time.baseline <-
rep(1:n.timepoints,n.subj.per.tx*2)*as.numeric(drug=="D")
dv <- rnorm( n.subj.per.tx*n.timepoints*2,
mean=time.baseline+Patient.baseline+drug.baseline, sd=sd.res
) dat.new <- data.frame(time, drug, dv, Patient)
dat.new$Patient.cross <- rep(1:(n.subj.per.tx), each = 2*n.timepoints)
dat.new$Dind <- as.numeric(dat.new$drug == "D") dat.new$Pind
<- as.numeric(dat.new$drug == "P") dat.new$time.num =
rep(1:n.timepoints, n.subj.per.tx*2)
##################################################################
R version 2.6.2 (2008-02-08)
i386-pc-mingw32
locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United
States.1252;LC_MONETARY=English_United
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] lme4_0.99875-9 Matrix_0.999375-5 lattice_0.17-4
loaded via a namespace (and not attached):
[1] grid_2.6.2
( fm.het.3 <- lmer( dv ~ time.num*drug -1 + ( 0 + Dind + Pind |
Patient.cross ), data=dat.new ) )
Linear mixed-effects model fit by REML
Formula: dv ~ time.num * drug - 1 + (0 + Dind + Pind | Patient.cross)
Data: dat.new
AIC BIC logLik MLdeviance REMLdeviance
684.7 706.3 -335.4 669.2 670.7
Random effects:
Groups Name Variance Std.Dev. Corr
Patient.cross Dind 26.8380 5.1805
Pind 7.7623 2.7861 0.060
Residual 1.5906 1.2612
number of obs: 160, groups: Patient.cross, 20
Fixed effects:
Estimate Std. Error t value
time.num 0.9101 0.1261 7.216
drugD -0.2637 1.2088 -0.218
drugP 5.0330 0.7123 7.066
time.num:drugP -0.9476 0.1784 -5.313
Correlation of Fixed Effects:
tim.nm drugD drugP
drugD -0.261
drugP 0.000 0.050
tim.nm:drgP -0.707 0.184 -0.313
#####################################################################
R version 2.7.1 (2008-06-23)
i386-pc-mingw32
locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United
States.1252;LC_MONETARY=English_United
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] lme4_0.999375-24 Matrix_0.999375-11 lattice_0.17-8
loaded via a namespace (and not attached):
[1] grid_2.7.1
( fm.het.3 <- lmer( dv ~ time.num*drug -1 + ( 0 + Dind + Pind |
Patient.cross ), data=dat.new ) )
Linear mixed model fit by REML
Formula: dv ~ time.num * drug - 1 + (0 + Dind + Pind | Patient.cross)
Data: dat.new
AIC BIC logLik deviance REMLdev
705 729.7 -344.5 687.6 689
Random effects:
Groups Name Variance Std.Dev. Corr
Patient.cross Dind 29.3378 5.4164
Pind 4.8857 2.2104 0.169
Residual 1.9651 1.4018
Number of obs: 160, groups: Patient.cross, 20
Fixed effects:
Estimate Std. Error t value
time.num 0.8175 0.1402 5.832
drugD -0.8505 1.2705 -0.669
drugP 5.9720 0.6258 9.543
time.num:drugP -1.1262 0.1982 -5.681
Correlation of Fixed Effects:
tim.nm drugD drugP
drugD -0.276
drugP 0.000 0.127
tim.nm:drgP -0.707 0.195 -0.396