Interpretation of various output of "lmer" function in R

[Cynderella]

  library(lme4)
  fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)

The notation `(Days | Subject)` says to allow the intercept and `Days` to vary randomly for each level of `Subject` .

Can you please explain me the result of the following commands ?

    attr(summary(fm1)$varcor$Subject,"stddev")
    (Intercept)        Days 
     24.740448    5.922133 

    c(sd(ranef(fm1)$Subject[,1]),sd(ranef(fm1)$Subject[,2]))
    [1] 21.595943  5.455217

    summary(fm1)$sigma
    [1] 25.59182

    residuals(summary(fm1))
  
    sd(residuals(summary(fm1)))
    [1] 0.9183965

What is the INTERPRETATION of the results found from various commands?

That is , if one asks me what is the meaning of the results that you have found from "sd(ranef(fm1)$Subject[,1])" and "attr(summary(fm1)$varcor$Subject,"stddev")[1]" ? Both are standard deviation of "Intercept" but of course there is difference between these two results . But I don't know what is this ?

In "?getMe" , it is said that from "summary(fm1)$sigma" , we found residual standard error . But why doesn't the result match with "sd(residuals(summary(fm1)))" ?

Also , In "summary(fm1)$varcor" there is value 0.066 under the column "Corr" . Does it mean correlation between two random effects "(Intercept)" and "Days " is 0.066 ?


Any help is appreciated . Thank you .

 

[Jake]

    attr(summary(fm1)$varcor$Subject,"stddev")
    (Intercept)        Days 
     24.740448    5.922133 

    c(sd(ranef(fm1)$Subject[,1]),sd(ranef(fm1)$Subject[,2]))
    [1] 21.595943  5.455217

The difference between these two pairs of quantities is subtle but conceptually important. The first pair are the actual parameter estimates: they are our best guess about the standard deviation of the intercepts and the standard deviation of the slopes in the population of Subjects. The second pair are the sample standard deviations of the BLUPs for the Subjects in our particular study. This is typically (or maybe always, but I hesitate to make the general statement) lower than the estimated standard deviation for the population because the BLUPs have been subjected to some degree of shrinkage (some info here). Because of this shrinkage, the standard deviation of our predictions for a particular sample does not equal our estimate of the standard deviation in the population.

    summary(fm1)$sigma
    [1] 25.59182

This is the estimate of the standard deviation of the errors. It is similar conceptually to our estimates of the standard deviations of the random effects.

    residuals(summary(fm1))
  
    sd(residuals(summary(fm1)))
    [1] 0.9183965

I actually didn't know before now that you could call residuals() on a summary.merMod object. But apparently you can, and the results are the "scaled residuals", i.e., the raw residuals (which you would obtain with residuals(fm1) -- omitting the summary() command) divided by sigma. In other words:

identical(resid(summary(fm1)), resid(fm1)/summary(fm1)$sigma)
# [1] TRUE

 

[TheEcologist]

I actually didn't know before now that you could call residuals() on a summary.merMod object. But apparently you can, and the results are the "scaled residuals", i.e., the raw residuals (which you would obtain with residuals(fm1) -- omitting the summary() command) divided by sigma. In other words:

identical(resid(summary(fm1)), resid(fm1)/summary(fm1)$sigma)
# [1] TRUE

Yeah, lme4 is constantly evolving and adding more s3/s4 methods to the packages. Two years ago predict didn't work on a merMod object, now it does. It helps to read the the release notes. It's pain because it also invariably breaks legacy code, in each version they seem to change the slotNames or input forms (re.form, ReForm, REForm, REform) - but it's still the best package for hierarchical work short of going Bayesian.

 

[Cynderella]

    attr(summary(fm1)$varcor$Subject,"stddev")
    (Intercept)        Days 
     24.740448    5.922133 

    c(sd(ranef(fm1)$Subject[,1]),sd(ranef(fm1)$Subject[,2]))
    [1] 21.595943  5.455217

The first pair are the actual parameter estimates: they are our best guess about the standard deviation of the intercepts and the standard deviation of the slopes in the population of Subjects.

Does it mean the first pair are the point estimates of variance components of random effects ? This link https://stat.ethz.ch/pipermail/r-sig-mixed-models/2014q3/022647.html tells attr(summary(fm1)$varcor$Subject,"stddev") is the point estimates of variance components of random effects . Is it ?

Thank you again . Regards .

 

[Cynderella]

    attr(summary(fm1)$varcor$Subject,"stddev")
    (Intercept)        Days 
     24.740448    5.922133 

    c(sd(ranef(fm1)$Subject[,1]),sd(ranef(fm1)$Subject[,2]))
    [1] 21.595943  5.455217

The first pair are the actual parameter estimates: they are our best guess about the standard deviation of the intercepts and the standard deviation of the slopes in the population of Subjects.

Does it mean the first pair are the point estimates of variance components of random effects ? This link https://stat.ethz.ch/pipermail/r-sig-mixed-models/2014q3/022647.html tells attr(summary(fm1)$varcor$Subject,"stddev") is the point estimates of variance components of random effects . Is it ?

Thank you again . Regards .

 

[Jake]

Yes, exactly.