Prediction of random effects in glmer()

Hi Ben,

I note with interest below:

    Note that these SDs **are not appropriate for testing hypotheses

    about individual levels of the random effect**: one of the tradeoffs of

    using REs in the first place is that you trade parsimony for the ability

    to do hypothesis testing on individual predictions.


Does this mean that they are also inappropriate for computing interval 
predictions for individual predictions?

Where might one read more about this observation?

Many thanks,

Andrew

-- 
Andrew Robinson
Director, CEBRA and Professor of Biosecurity,
School/s of BioSciences and Mathematics & Statistics
University of Melbourne, VIC 3010 Australia
Tel: (+61) 0403 138 955 Email:?apro at unimelb.edu.au
Website:?https://researchers.ms.unimelb.edu.au/~apro at unimelb/
I?acknowledge the Traditional Owners of the land I inhabit, and pay my 
respects to their Elders.
On Feb 16, 2021, 8:16 AM +1100, Ben Bolker <bbolker at gmail.com>, wrote:

The penalty term ||u||^2 appears first in the likelihood expression
(eq 11 in the paper),

int exp( - (sum(d_i(y_obs,u)) + ||u||^2)/2) * C du

where C is a normalization constant I don't feel like writing out.

The sum(d_i(.)) part is the deviance of the GLM part of the model.
The other part would come out to exp(-||u||^2/2) if we factored it out;
this expresses the idea that the u values (the 'spherized' random
effects, i.e. on a scale where they are iid N(0,1) variables) are
assumed to be drawn from a Normal distribution.

This stuff is probably better explained in the original JSS lmer
paper (where the same penalty term appears, for the same reason).

An *analog* of standard errors is indeed available for the
conditional modes (this is also explained in the JSS paper); we call
these "conditional variances" and "conditional standard deviations" - we
can't quite call them "standard errors" because the u values are not
estimates in the technical sense (this is the same song-and-dance as
when we call something a "BLUP" rather than an "estimate").

When you call ranef(.) with condVar=TRUE (which is the default), the
conditional covariance matrices are indeed computed and returned, albeit
in a rather inconvenient form: attr(ranef(my_model)$grpvar,"postVar") is
a 3D array of stacked covariance matrices, or a list of 3D arrays (see
?ranef.merMod). If you just want the conditional standard deviations
for each random effect value, as.data.frame(ranef(my_model)) will
probably be more convenient.

(Again, the details behind this computation are explained in the JSS
lmer paper ...)

Note that these SDs **are not appropriate for testing hypotheses
about individual levels of the random effect**: one of the tradeoffs of
using REs in the first place is that you trade parsimony for the ability
to do hypothesis testing on individual predictions.

For your second question, I'm not sure what you mean by the
"conditional mean and conditional variance of the random effect"; can
you explain further/give an example?

Ben



On 2/15/21 3:24 PM, Ravi Varadhan wrote:

Dear Ben,
Thanks for your response. I went back and looked at the draft JSS 
paper you sent me. It does describe how the random effects are 
predicted as conditional modes, using a penalized, iteratively 
weighted least squares. However, I still have some questions. Why is 
the penalty term ||u||^2 added? What does this mean? Does glmer then 
provide standard errors for the predicted random effects (I don't 
think it does)?

One more question: it would be nice to also have an option for 
conditional mean and conditional variance of the random effect, 
although conditional variance would underestimate the true variance 
of the prediction.

Thank you,
Ravi