Dear all,

I know that the subject of a predict() function for lme4 / lmer -- or lack
thereof -- has been kicked around before (I googled "lmer predict").  After
reading all of the previous e-mails, I'd like to get feedback on the various
ways that I've thought up to actually do prediction from lmer().  Most of
the ideas are briefly sketched out below; any feedback or comments would be
appreciated.

David


Problem Statement:

After estimating my model from a given specification, I'd like to get both
(a) predictions from the model and the original covariates, and (b)
predictions from the same model applied to new covariates.  In pseudocode:

Step 1:  Original covariates "origdat", columns "x1" ... "x4", groups
indicated by "group"
Step 2:  Estimate model, get object "M1".  My models are relatively simple,
something like:  M1 <- lmer(y ~ x1 + x2 + (x3 + x4 | group), data=origdat)
Step 3:  New covariates "newdat", with same column names as "dat", but
different data (changed covariates represent a counterfactual scenario).


After thinking about it, it seems like my options are:

[Option #1]  Write out the BLUPs and do simulations manually.  From the M1
object I could extract the fixed effects, random effects,
variance-covariance matrix, and residual sds, and so could just simulate by
hand.  In other words, write something like this:

orig.pred <- rnorm(fixef(M1) %*% origdat + mvnorm(ranef(M1) %*%
origdat[group], VarCorr(M1)), residsd).

Then apply the same parameters to the new data:

new.pred <- rnorm(fixef(M1) %*% newdat + mvnorm(ranef(M1) %*% newdat[group],
VarCorr(M1)), residsd).

Pros:  Probably fast.
Cons:  Really messy since I'm doing a number of different models, unless I
write a really clean function, and perhaps beyond my programming skills.


[Option #2]  Think of some creative use of refit() and simulate() -- or
other functions -- for lmer.  The documentation for refit in lme4.pdf is
pretty tempting, since it writes that "it is common to use the same model
specification and covariate data on many simulated responses.  That is, the
only thing that changes between model fits is the response vector."
 However, in this case I want to keep the model specification, but change
both the covariate data and get the response vector.  I'm not sure if this
is possible using lmer().

Cons: obviously, I don't know how quite to do this yet.


[Option #3]  Just use nlme().

Pros: nlme() will estimate this relatively simple model, and then I can
simply use the predict() function to get the predictions using the fit to
the original data, and then predict(M1, newdat) to get predictions with the
new data.
Cons: nlme() is slower than lmer(), and my model already takes about 2 hours
to run.






-- 
David Hsu
PhD Candidate, Urban Planning
University of Washington
(e-mail) dhsu2 at uw.edu

	[[alternative HTML version deleted]]


