[R-sig-ME] extracting info from a mer object
Douglas Bates
bates at stat.wisc.edu
Wed Mar 30 18:32:49 CEST 2011
On Wed, Mar 30, 2011 at 9:49 AM, Robert Kinley <KINLEY_ROBERT at lilly.com> wrote:
> hello List
>
> having fitted a model with lmer() I'd like to better understand how it
> works
> by looking at components of the fitting process ... for example, to see
> the
> 'Lambda' matrix, Doug Bates's "LME4:Mixed effects modelling with R"
> suggests the command
>
>> env(mymodel)$Lambda
>
> but in my case I get:-
>
> > require(lme4)
> > ... fit model etc ...
> >env(QV.lmer)$Lambda
> >Error: could not find function "env"
>
>
> I daresay I'm making some fundamental error ... any idea what it is ?
Assuming that I will be consistent :-)
The approach of storing the pieces of the model in an environment gave
way to using certain S4 classes. It turns out that in lme4 the Lambda
matrix has to be calculated from a permutation, P, a diagonal matrix S
and a unit triangular matrix T available as the result of expand.
However, "the times they are a'changing" again and in the development
version called lme4a the Lambda matrix is back as an explicit matrix
in the re slot of the model.
I think the best approach would be to have an extractor function to
provide that matrix from the model so that I can hide the
implementation details, which keep changing. I'll discuss this with
Martin and Ben.
For now you should be able to get the Lambda matrix from a model fit by lme4 as
> library(lme4)
> (fm1 <- lmer(diameter ~ 1 + (1|plate) + (1|sample), Penicillin))
Linear mixed model fit by REML
Formula: diameter ~ 1 + (1 | plate) + (1 | sample)
Data: Penicillin
AIC BIC logLik deviance REMLdev
338.9 350.7 -165.4 332.3 330.9
Random effects:
Groups Name Variance Std.Dev.
plate (Intercept) 0.71691 0.84670
sample (Intercept) 3.73092 1.93156
Residual 0.30242 0.54992
Number of obs: 144, groups: plate, 24; sample, 6
Fixed effects:
Estimate Std. Error t value
(Intercept) 22.9722 0.8085 28.41
> image(Lambda <- with(expand(fm1), S %*% T))
More information about the R-sig-mixed-models
mailing list