[R] Getting residual term out of lmer summary table
Ben Bolker
bbolker at gmail.com
Tue Nov 12 16:34:34 CET 2013
<aline.frank <at> wsl.ch> writes:
>
> Hello
> I'm working with mixed effects models using lmer() and have some
> problems to get all variance components of the model's random
> effects. I can get the variance of the random effect out of the
> summary and use it for further calculations, but not the variance
> component of the residual term. Could somebody help me with that
> problem? Thanks a lot! Below an example.
There's an r-sig-mixed-models at r-project.org mailing list that
is specifically for mixed-models (especially lme4) questions.
> Aline
> ## EXAMPLE
> #----------
>
require(lme4)
## Simulate data for the example
set.seed(6)
x1 <- runif(n=100, min=10, max=100) ## a continuos variable
x2 <- runif(n=100, min=10, max=100) ## a continuos variable
treat <- rep(letters[1:4], times=25) ## a fixed factor with 4 levels
treat.effect <- 20*rep(1:4, times=25)
group.label <- rep(LETTERS[1:5], each=20) ## the random effect
group.effect <- 10*rep(1:5, each=20) ## there are 5 groups
## Response variable:
y <- 2*x1 + (-5)*x2 + treat.effect + group.effect + rnorm(100)
## Dataframe
d.ex <- data.frame(y, x1, x2, "Group"=group.label, treat)
## Apply model
mod1 <- lmer(y~x1+x2+treat+x1:treat+ (1|Group), data=d.ex)
output <- summary(mod1); output
## # ok, there is the variance component of the random effect group and the
residual term
>
## Now I'd like to get the variance components of the random effect "Group"
and of the residual term
"Residual" in order
## to do further calculations with these numbers
output$varcor[1] ## reveals the variance of the random effect "Group"
output$varcor[2] ## does not reveal the residual term! what other command do
I need to use then?
VarCorr(mod1) ## prints standard dev by default in latest lme4
## print both var and std dev
print(VarCorr(mod1),comp=c("Variance","Std.Dev."))
unlist(VarCorr(mod1)) ## group variance (as raw number)
sigma(mod1) ## std. dev. of residual variance
More information about the R-help
mailing list