On Fri, 2007-11-09 at 18:55 +0000, Prof Brian Ripley wrote:
> I think Bert's point is important: I picked up a student on it in a case
> study presentation on this week because I could think of three
> interpretations, none strictly confidence intervals. I think 'tolerance
> interval' is fairly standard for prediction of a random quantity: see
> ?predict.lm.
>
I think prediction interval is what is usually used. Regardless, I'm not
sure how "predict.lm" will be of much help because I asked specifically
about BLUP's for random effects and the last time I checked lm did not
handle mixed effects models. Neither predict.lme and predict.lmer
provide intervals. Here is the code that I included in my original
e-mail. My simple question is, will this code correctly compute a
prediction interval for each subjects random effect? In particular, will
the code handle the bVar slot correctly? Some postings warned about
inappropriate access to slots. Here is the code that I asked about in my
original e-mail:
# OrthoFem has all the females from Orthodont from the nlme package
library(lme4)
fm1OrthF. <- lmer(distance~age+(age|Subject), data=OrthoFem)
lmer(distance~age+(age|Subject), data=OrthoFem)@bVar$Subject[2,2,]*
(attr(VarCorr(lmer(distance~age+(age|
Subject),data=OrthoFem)),"sc")^2)[1]
(attr(VarCorr(fm1OrthF.),"sc")^2)[1]
fm1.s <- coef(fm1OrthF.)$Subject
fm1.s.var <- [EMAIL PROTECTED](attr(VarCorr(fm1OrthF.),"sc")^2)[1]
fm1.s0.s <- sqrt(fm1.s.var[1,1,])
fm1.s0.a <- sqrt(fm1.s.var[2,2,])
fm1.s[,1]+outer(fm1.s0.s, c(-2,0,2))
fm1.s[,2]+outer(fm1.s0.a, c(-2,0,2))
> fm1.s
(Intercept) age
F10 14.48493 0.3758608
F09 17.26499 0.3529804
F06 16.77328 0.3986699
F01 16.95609 0.4041058
F05 18.36188 0.3855955
F07 17.28390 0.5193954
F02 16.05461 0.6336191
F08 19.40204 0.3562135
F03 16.35720 0.6727714
F04 19.02380 0.5258971
F11 19.13726 0.6498911
> fm1.s[,1]+outer(fm1.s0.s, c(-2,0,2))
[,1] [,2] [,3]
[1,] 12.21371 14.48493 16.75616
[2,] 14.99377 17.26499 19.53622
[3,] 14.50205 16.77328 19.04450
[4,] 14.68487 16.95609 19.22732
[5,] 16.09066 18.36188 20.63311
[6,] 15.01267 17.28390 19.55512
[7,] 13.78339 16.05461 18.32584
[8,] 17.13082 19.40204 21.67327
[9,] 14.08598 16.35720 18.62843
[10,] 16.75257 19.02380 21.29502
[11,] 16.86604 19.13726 21.40849
> fm1.s[,2]+outer(fm1.s0.a, c(-2,0,2))
[,1] [,2] [,3]
[1,] 0.1738325 0.3758608 0.5778890
[2,] 0.1509522 0.3529804 0.5550087
[3,] 0.1966417 0.3986699 0.6006982
[4,] 0.2020775 0.4041058 0.6061340
[5,] 0.1835672 0.3855955 0.5876237
[6,] 0.3173671 0.5193954 0.7214236
[7,] 0.4315909 0.6336191 0.8356474
[8,] 0.1541852 0.3562135 0.5582417
[9,] 0.4707432 0.6727714 0.8747997
[10,] 0.3238688 0.5258971 0.7279253
[11,] 0.4478629 0.6498911 0.8519194
This web page describes "bVar":
http://hosho.ees.hokudai.ac.jp/~kubo/Rdoc/library/lme4/html/lmer-class.html
bVar:
A list of the diagonal inner blocks (upper triangles only) of
the positive-definite matrices on the diagonal of the inverse of
ZtZ+Omega. With the appropriate scale factor (and conversion to
a symmetric matrix) these are the conditional
variance-covariance matrices of the random effects.
Rick B.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
