[R-sig-ME] lme vs. lmer
Peter Dalgaard
p.dalgaard at biostat.ku.dk
Tue Sep 29 22:58:42 CEST 2009
Ben Bolker wrote:
> Douglas Bates wrote:
>> On Tue, Sep 29, 2009 at 1:02 PM, Ben Bolker <bolker at ufl.edu> wrote:
>>
>>> Christopher David Desjardins wrote:
>>>> I've started working through Pinheiro & Bates, 2000 and noticed the use
>>>> of lme from the nlme package. I am curious if lmer from lme4 has
>>>> superseded lme or if lme still holds its own? The reason I ask is that I
>>>> have taken a few classes where we've solely used lmer and just read
>>>> about lme today. If both functions are on equal footing, can the
>>>> p-values from lme be trusted?
>>>> Thanks!
>>>> Chris
>>> You should read the extended discussion of p-values, degrees of
>>> freedom, etc. that is on the R wiki (I think) and referenced from the R
>>> FAQ. At least in my opinion, (n)lme is still fine (and indeed necessary
>>> at this stage for fitting heteroscedastic and correlated models). The
>>> df/p-value estimates, however, are "use at your own risk" -- you'll have
>>> to read the literature and decide for yourself.
>>> I still think there's room for someone to implement (at least)
>>> Satterthwaite and (possibly) Kenward-Roger corrections, at least for the
>>> sake of comparison, but I'm not volunteering.
>> You may need to define them first. Many of the formulas in the mixed
>> models literature assume a hierarchical structure in the random
>> effects - certainly we used such a formula for calculating the
>> denominator degrees of freedom in the nlme package. But lme4 allows
>> for fully or partially crossed random effects so you can't think in
>> terms of "levels" of random effects.
>>
>> Referring to the "Satterthwaite and Kenward-Roger corrections" gives
>> the impression that these are well-known formulas and implementing
>> them would be a simple matter of writing a few lines of code. I don't
>> think it is. I would be very pleased to incorporate such code if it
>> could be written but, as I said, I don't even know if such things are
>> defined in the general case, let alone easy to calculate.
>>
>> I am not trying to be argumentative (although of late I seem to have
>> succeeded in being that). I'm just saying that I don't think this is
>> trivial. (It I wanted to be argumentative I would say that it is
>> difficult and, for the most part, irrelevant. :-)
>
> Fair enough. Actually, I'm not sure I meant implementing them in lmer
> -- even implementing them in nlme would be useful (and perhaps more
> straightforward, if not trivial). I also wouldn't impose the requirement
> that they have to be feasible for huge data sets -- I'm just curious if
> they can be implemented within lme in a relatively straightforward/
> boneheaded way.
> But again, this is very far down my to-do list (and at the edge
> of my abilities) and completely off yours, so unless someone else bites
> it won't happen.
I don't think (n)lme is all that easy either; you still need to sort out
the connection between its multilevel formulation and the projection
matrices in the K-R paper. In both nlme and lme4, an implementation is
almost certainly possible, although probably complicated and perhaps at
the expense of all computational efficiency.
One main problem is that even when they can be calculated, the
corrections rely on a normal distribution assumption which is more than
likely wrong in practice. This isn't any different from ordinary
t-tests: once you get into the single-digit df regime, you really don't
know what you are doing, and if there are more than 30 df, the normal
approximation works well enough without the correction.
Accordingly, I tend to see low df more as a warning flag than as
something that gives accurate p values, and I sometimes wonder whether
there is a way to raise such a flag more expediently.
--
O__ ---- Peter Dalgaard Øster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907
More information about the R-sig-mixed-models
mailing list