Skip to content

Extracting Standard Errors of Uncorrelated Random Effects?

2 messages · Derek Dunfield, Douglas Bates

#
Hello all:
I was hoping to extract SE for my random effects to get a better idea
of their precision for each grouping factor.

Andrew Gelman has developed a nice package called "Arm" which does
this seamlessly for models with correlated random effects such as
Reaction ~ 1 + Days + (1 + Days | Subject)

When I try to run this on an uncorrelated model such as
Reaction ~ 1 + Days + (1 | Subject) + (0 + Days | Subject)

I get the error "Error in .local(object, ...) : Code not written yet"

I believe this may be due to some unfinished code in the lme4 package
for ranef, though I'm not sure.

Does anyone have an idea how I could go about calculating these standard errors?

Many thanks in advance for you help!
:Derek

--
Derek Dunfield, PhD
Postdoctoral Fellow, MIT Intelligence Initiative
Sloan School of Management
MIT Center for Neuroeconomics, Prelec Lab
77 Massachusetts Ave E62-585
Cambridge MA 02139
#
On Wed, Dec 14, 2011 at 1:05 PM, Derek Dunfield <dunfield at mit.edu> wrote:
Unfortunately, it is complicated.  Strangely enough, the case of
correlated random effects is easier than the case of uncorrelated
effects.  The values that you want to get are the solutions to a
system like L x = cbind(e_i, e_j) where the indices i and j are the
random effects corresponding to the k'th level of the grouping factor
and L is the sparse Cholesky factor.  For correlated random effects
the i and j are going to be adjacent.  For uncorrelated random effects
they may end up being adjacent after the fill-reducing permutation but
that is not guaranteed.

In this particular case things work out fine.  If I fit that model as
fm1 then I can check that the random effects for each subject are
adjacent by examining

image(fm1 at L)

where we can see that each of the pairs of blocks starting with 1:2
are correlated within the block but orthogonal to the other blocks.
So, setting

L <- fm1 at L
Id <- Diagonal(x = rep.int(1, ncol(L))   # Identity matrix of the same size as L

 the conditional variance-covariance matrix of the orthogonal random
effects for subject 1 is something like

crossprod(solve(L, Id[,1:2], sys="L")) * lme4:::sigma(fm1)^2

Before using such an expression I would check it against the results
from the "Arm" package for the model with correlated random effects.
It is entirely possible that I have left something out of the
calculation.

If you wonder why something that looks simple like this is not part of
the lme4 package, it is because a lot of the simplicity comes from the
fact that both random-effects terms have the same grouping factor.
The general case is much messier.