[R] Extracting level-1 variance from lmer()

Doran, Harold HDoran at air.org
Tue Feb 5 20:35:11 CET 2008


Dave:

This is an attribute. So, you can get it as follows:

library(lme4)
example(lmer) 
attr(VarCorr(fm1), "sc")



> -----Original Message-----
> From: r-help-bounces at r-project.org 
> [mailto:r-help-bounces at r-project.org] On Behalf Of David Afshartous
> Sent: Tuesday, February 05, 2008 2:37 PM
> To: r-help at r-project.org
> Subject: [R] Extracting level-1 variance from lmer()
> 
> 
> All,
> 
> How does one extract the level-1 variance from a model fit via lmer()?
> 
> In the code below the level-2 variance component may be 
> obtained via subscripting, but what about the level-1 
> variance, viz., the 3.215072 term?
> (actually this term  squared)  Didn't see anything in the 
> archives on this.
> 
> Cheers,
> David
> 
> 
> 
> > fm <- lmer( dv ~ time.num*drug + (1 | Patient.new),  data=dat.new )
> 
> > VarCorr(fm)
> $Patient.new
> 1 x 1 Matrix of class "dpoMatrix"
>             (Intercept)
> (Intercept)    8.519916
> 
> attr(,"sc")
>    scale
> 3.215072 
> 
> > VarCorr(fm)[[1]][1]
> [1] 8.519916
> > VarCorr(fm)[[2]][1]
> Error in VarCorr(fm)[[2]] : subscript out of bounds
> 
> 
> 
> 
> 
> 
> 
> ##########################################################
> set.seed(500)
> n.timepoints <- 4
> 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.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(drug, dv)
> 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) dat.new$Patient.new = 
> rep(1:20, each=8)
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 



More information about the R-help mailing list