[R] lme - Confidence interval for predicted values

[Douglas Bates]

"Urs Wiedemann" <wiedemann at fmp-berlin.de> writes:

> After having fit (using lme) a mixed effects model with a single
> random effect, I like to estimate the confidence interval for the
> predicted mean expectations.
> To my knowledge this is usually done (for fixed effects models) by
> calculating:
> 
> cibandwidth <- sqrt(diag(Xnew %*% solve(t(X) %*% X) t(Xnew))) *
> qt((1-level)/2, DF)
> 
> The CI is then the predicted value +/- cibandwidth. This is what
> predict.lm provides with ci.fit if I am not wrong.
> 
> So for lme there is a very nice post on the S-news list:
> http://www.biostat.wustl.edu/archives/html/s-news/2003-09/msg00021.html
> 
> Hopefully I quote the message correctly:
> >var(y.hat) = X2 %*% inverse(t(X)%*%inverse(Sig)%*%X) %*% t(X2)
> >
> >The hard part for lme is deciding what goes in X and what is Sig:
> >If you want a confidence interval on y for a particular subject
> >included in X, then that subject is part of X; otherwise,
> >it must be included in Sig.

I think the answer to this question will soon become easier.  Later
today (if I can escape administrative duties, otherwise tomorrow) I
will make available an alpha-level release candidate of a completely
rewritten lme4 package.  And I do mean completely rewritten.  This is
the fourth and, I swear, the last time I have designed and implemented
from scratch methods for parameter estimation in linear mixed models.

The new version of lme (that is, from the not-yet-released lme4_0.6-1
or later) allows the usual optional arguments for model-fitting
functions in R, including "x=TRUE" which will cause the slot x in the
fitted model object to contain a list of model matrices.  The last
element of this list is the X that you want.

> My question is now where to obtain "Sig" from lme. Probably this is
> obvious and I simply overlooked it.

The new version of the lme4 package has a method for the "vcov"
generic.  This method returns this Sig.