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
Extracting Standard Errors of Uncorrelated Random Effects?
2 messages · Derek Dunfield, Douglas Bates
On Wed, Dec 14, 2011 at 1:05 PM, Derek Dunfield <dunfield at mit.edu> wrote:
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"
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.
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
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models