fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
class(fm1)
## Here the adjusted and unadjusted covariance matrices are identical,
## but that is not generally the case:
v1 <- vcov(fm1)
v2 <- vcovAdj(fm1, details=0)
v2 / v1
## For comparison, an alternative estimate of the variance-covariance
## matrix is based on parametric bootstrap (and this is easily
## parallelized):
if (FALSE) {
nsim <- 100
sim <- simulate(fm.ml, nsim)
B <- lapply(sim, function(newy) try(fixef(refit(fm.ml, newresp=newy))))
B <- do.call(rbind, B)
v3 <- cov.wt(B)$cov
v2/v1
v3/v1
}
Run the code above in your browser using DataLab