# Using lmer for prediction

Hello I have two problems that sound like natural candidates for multilevel/mixed models, which I have never used. The simpler, and one that I hope to try as an introduction, is as follows: The data looks like many rows of the form

x y innergroup outergroup

where x is a numeric covariate upon which I want to regress y (another numeric variable), each y belongs to an innergroup, and each innergroup is nested in an outergroup (i.e, all the y in a given innergroup belong to the same outergroup). Unfortunately, innergroup has a lot of levels (many thousands), and each level has relatively few observations of y, so I thought this sort of model might be appropriate. My questions are

1. How do I write this sort of multilevel formula?

2. Once lmer fits the model, how does one go about predicting from it? I have fit some simpler toy examples, but have not found a predict() function. Most people seem more interested in inference than prediction with this sort of technique. I have several million rows, so the computations might be an issue, but I can always cut it down as appropriate.

I won't need to do the second for some time, but I might as well begin thinking about it and playing around with it. I have similar data as before, but without x, and y is now a binomial variable of the form $(n,n-k)$. y also exhibits a lot of overdispersion, even within innergroups. Most of the $n$ are no more than 2 or 3 (or less), so to derive estimates of the success rates of each $y_i$ I have been using the beta-binomial shrinkage estimator $(\alpha+k_i)/(\alpha+\beta+n_i)$, where $\alpha$ and $\beta$ are estimated by MLE for each innergroup separately. This is has been somewhat adequate, but data sparsity still plagues me, so I would like to use all the data available. From one perspective, this problem is easier since there is no covariate, but from the other the binomial nature makes it more difficult. Does anyone have any high (or low!) level guidance?

-
Please verify that the parentheses I inserted into the shrinkage formula are where you intended them to be. – whuber Nov 4 '10 at 4:10
Shouldn't the 2nd part of your question (with a binary variable) be a separate question? – chl Nov 4 '10 at 11:45

Expressing factors relationships using R formulas follows from Wilkinson's notation, where '*' denotes crossing and '/' nesting, but there are some particularities in the way formula for mixed-effects models, or more generally random effects, are handled. For example, two crossed random effects might be represented as (1|x1)+(1|x2). I have interpreted your description as a case of nesting, much like classes are nested in schools (nested in states, etc.), so a basic formula with lmer would look like (unless otherwise stated, a gaussian family is used by default):

y ~ x + (1|A:B) + (1|A)

where A and B correspond to your inner and outer factors, respectively. B is nested within A, and both are treated as random factors. In the older nlme package, this would correspond to something like lme(y ~ x, random=~ 1 | A/B). If A was to be considered as a fixed factor, the formula should read y ~ x + A + (1|A:B).

But it is worth checking more precisely D. Bates' specifications for the lme4 package, e.g. in his forthcoming textbook, lme4: Mixed-effects Modeling with R, or the many handouts available on the same webpage. In particular, there is an example for such nesting relations in Fitting Linear Mixed-Effects Models, the lme4 Package in R. John Maindonald's tutorial also provides a nice overview: The Anatomy of a Mixed Model Analysis, with R’s lme4 Package. Finally, section 3 of the R vignette on lme4 imlementation includes an example of the analysis of a nested structure.

There is no predict() function in lme4, and you have to compute yourself predicted individual values using the estimated fixed (see ?fixef) and random (see ?ranef) effects, but see also this thread on the lack of predict function in lme4. You can also generate a sample from the posterior distribution using the mcmcsamp() function. Sometimes, it might clash, though. See the sig-me mailing list for more updated information.

-

The ez package contains the ezPredict() function, which obtains predictions from lmer models where prediction is based on the fixed effects only. It's really just a wrapper around the approach detailed in the glmm wiki.

-

I would use the "logit.mixed" function in Zelig, which is a wrapper for lime4 and makes it very convenient to do prediction and simulation.

-

Stephen Raudenbush has a book chapter in the Handbook of Multilevel Analysis on "Many Small Groups". If you are only interested in the effects of x on y and have no interest in higher level effects, his suggestion is simply to estimate a fixed effects model (i.e. a dummy variable for all possible higher level groupings).

I don't know how applicable that is towards prediction, but I would imagine some of what he writes is applicable to what you are trying to accomplish.

-