[R-sig-ME] Fwd: extracting coefficients from lmer2 model

[Douglas Bates]

I meant to cc: the list on this reply

---------- Forwarded message ----------
From: Douglas Bates <bates at stat.wisc.edu>
Date: Feb 1, 2007 1:44 PM
Subject: Re: [R-sig-ME] extracting coefficients from lmer2 model
To: Colin Beale <c.beale at macaulay.ac.uk>


On 2/1/07, Colin Beale <c.beale at macaulay.ac.uk> wrote:
> Hi all,
>
> Some of you may have seen my recent message to the r-help list
> concerning a modification to lmer code that caused R to crash. Happily,
> the same modifications to the lmer2 code do not cause this crash so as
> far as I am concerned this is a major improvement to the original. Now,
> however, I'm not sure how to get at the bits of the fitted object that I
> want to play with to post-process the results as the methods aren't yet
> implemented for this structure - I need, for example, the coefficients
> of the fixed and random effects and the full variance-covariance matrix
> of fixed _and_ random effects (together!). I'm sure there must be ways
> around this issue, since I can see that there are slots with names that
> look very suitable at least to hold the fixed and random effects - but I
> can't work out how to get them. And I'm not sure if the full
> variance-covariance matrix is going to be in there anyway (looks like it
> might be for lmer objects).

I recommend the vignette in the lme4 package that describes the lmer2
representation.

In answer to some of your specific questions:

fixef is straightforward for the mer2 representation

The ranef slot is currently populated when you call for fixef but the
values in that slot are the random effects with the "spherical" prior
density (independent, constant-variance).  These are named b* in the
vignette.   I need to add code to undo the permutation and to undo the
orthogonalization transformation defined by S and T.  It is not
incredibly difficult but I will need some quiet time to consider
exactly what sequence of operations needs to be performed and then to
write and debug the code.

The variance-covariance matrix for the random effects and the fixed
effects can be derived from the inverse of L but do you really want to
create that object?  It has the potential of being very large and very
slow to calculate.  Earlier today there was discussion of the need to
set control = list(niterEM = 0, gradient = FALSE) in a call to lmer
that fits a complicated model to a large data set.  This reason that
the ECME iterations and the gradient calculations are so slow (on the
order of hundreds or thousands of evaluations of the deviance) is
because those calculations require part of the inverse of L.  It's
even simpler than that, all you need to do is to calculate some groups
of columns in the inverse and store some values derived from them.
You don't ever need to store the inverse itself - just these values
that are derived from groups of columns in the inverse.

Calculating the inverse of L won't be a problem for simple models and
small data sets but it will be a problem for multiple, non-nested
grouping factors with many levels.  Furthermore the sparse matrix
methods take care to minimize the amount of storage in L.  That's what
the fill-reducing permutation is designed to do.  However, the inverse
of L can be a mess.

You should consider carefully exactly what you need from that
variance-covariance matrix.  For example, if you are going to evaluate
a quadratic form using this matrix then don't calculate the matrix.
Write the quadratic form using the factorization and solve the system
involving L or L-transpose then collect an inner product.

> For example:
>
> fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
> fm2 <- lmer2(Reaction ~ Days + (Days|Subject), sleepstudy)
>
> using:
>
>  slotNames(fm1)
>  slotNames(fm2)
>
> shows me that there are slots with the ranef and fixef in both classes
> of object - I can use fixef() on both objects to get the fixed effects,
> but ranef() only works on fm1 at the moment - how else can I get these
> objects directly? And there are obviously many more slots in the fm1
> object than the fm2 object - can anyone point out if any of these slots
> are the full variance-covariance matrix I need, containing both fixed
> and random effects and if so, how I can get at it (or if not, how I
> might go about creating it?).
>
> Thanks in advance,
>
> Colin
>
>
> Dr. Colin Beale
> Spatial Ecologist
> The Macaulay Institute
> Craigiebuckler
> Aberdeen
> AB15 8QH
> UK
>
> Tel: 01224 498245 ext. 2427
> Fax: 01224 311556
> Email: c.beale at macaulay.ac.uk
>
>
> --
> Please note that the views expressed in this e-mail are those of the
> sender and do not necessarily represent the views of the Macaulay
> Institute. This email and any attachments are confidential and are
> intended solely for the use of the recipient(s) to whom they are
> addressed. If you are not the intended recipient, you should not read,
> copy, disclose or rely on any information contained in this e-mail, and
> we would ask you to contact the sender immediately and delete the email
> from your system. Thank you.
> Macaulay Institute and Associated Companies, Macaulay Drive,
> Craigiebuckler, Aberdeen, AB15 8QH.
>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>