Re: [R] bVar slot of lmer objects and standard errors

From: Spencer Graves <spencer.graves_at_pdf.com>
Date: Sun 26 Feb 2006 - 11:27:46 EST

          I shall provide herewith an example of what I believe is "the conditional posterior variance of a random effect". I hope someone more knowledgeable will check this and provide a correction if it's not correct. I say this in part because my simple rationality check (below) didn't match as closely as I had hoped.

          I don't have easy access to the "hlmframe" of your example, so I used the example in the "lmer" documentation:

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

    Data: sleepstudy

       AIC      BIC    logLik MLdeviance REMLdeviance
  1753.628 1769.593 -871.8141   1751.986     1743.628
Random effects:
  Groups   Name        Variance Std.Dev. Corr
  Subject  (Intercept) 612.090  24.7405
           Days         35.072   5.9221  0.066
  Residual             654.941  25.5918

# of obs: 180, groups: Subject, 18
<snip>

          I think we want the "Residual Variance" of 654.941 (Std.Dev. of 25.5918). To get this, let's try "VarCorr":

(vC.fm1 <- VarCorr(fm1))

$Subject
2 x 2 Matrix of class "dpoMatrix"

             (Intercept) Days
(Intercept) 612.09032 9.60428

Days 9.60428 35.07165

attr(,"sc")
[1] 25.59182

          Our desired 25.5918 is 'attr(,"sc")' here. We get that as follows:

(s. <- attr(vC.fm1, "sc"))

[1] 25.59182

          Next, by studying 'str(fm1)' and earlier emails you cite, we get the desired conditional posterior covariance matrix as follows:

 > (condPostCov <- (s.^2)*fm1@bVar$Subject) , , 1

          [,1] [,2]
[1,] 145.7051 -21.444432
[2,] 0.0000 5.312275
<snip>
, , 18

          [,1] [,2]
[1,] 145.7051 -21.444432
[2,] 0.0000 5.312275

          This actually gives us only the upper triangular portion of the desired covariance matrices. The following will make them all symmetric:

condPostCov[2,1,] <- condPostCov[1,2,]

          As a check, I computed the covariance matrix of the estimated random effects. To get this, I reviewed str(ranef(fm1)), which led me to the following:

  var(ranef.fm1@.Data[[1]])

             (Intercept) Days
(Intercept) 466.38503 31.04874


Days 31.04874 29.75939

          These numbers are all much larger than "condPostCov". However, I believe this must be a random bounce -- unless it's a deficiency in my understanding (in which case, I hope someone will provide a correction).

          Ulrich: Would you mind checking this either with a published example or a Monte Carlo and reporting the results to us?

	  Viel Glück,
	  spencer graves
	  	

Ulrich Keller wrote:

> Hello,
> 
> I'm sorry to resurrect this thread that I started almost two months ago. 
> I've been pretty busy since I posted my question and the issue is not 
> that high on my priority list. Thanks to all those who replied, and I 
> hope I can tickle your interest again.
> 
> As a reminder, my question was how one can extract the conditional 
> posterior variance of a random effect from the bVar slot of an lmer 
> model. Thanks to your answers, I now understand that I have to use the 
> diagonal elements of the conditional matrices. However, I am not quite 
> sure what this means:
> 
> Douglas Bates wrote:
> 
>>I'd have to go back and check but I think that these are the upper
>>triangles of the symmetric matrix (as Spencer suggested) that are the
>>conditional variance-covariance matrices of the two-dimensional 
>>random effects for each school up to a scale factor.  That is, I think
>>each face needs to be multiplied by s^2 to get the actual
>>variance-covariance matrix.
> 
> 
> What is s^2? Where can I find it in the lmer object? I tried reading the 
> source, but gave up fairly quickly. Thanks in advance for your replies, 
> and this time I promise I'll be more responsive.
> 
> 
> Uli
> 
> My original post:
> 
>>Hello,
>>
>>I am looking for a way to obtain standard errors for 
emprirical Bayes estimates of a model fitted with lmer
(like the ones plotted on page 14 of the document
available at
http://www.eric.ed.gov/ERICDocs/data/ericdocs2/content_storage_01/0000000b/80/2b/b3/94.pdf).

Harold Doran mentioned
(http://tolstoy.newcastle.edu.au/~rking/R/help/05/08/10638.html)
that the posterior modes' variances can be found in the bVar slot of lmer objects. However, when I fit e.g. this model:

>>
>>lmertest1<-lmer(mathtot~1+(m_escs_c|schoolid),hlmframe)
>>
>>then lmertest1@bVar$schoolid is a three-dimensional 
array with dimensions (2,2,28). The factor schoolid

has 28 levels, and there are random effects for the intercept and m_escs_c, but what does the third dimension correspond to? In other words, what are the contents of bVar, and how can I use them to get standard errors?

>>
>>Thanks in advance for your answers and Merry Christmas,
>>
>>Uli Keller
> 
> 
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sun Feb 26 11:36:46 2006

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:42:45 EST