[R] The mathematics inside lme()

Carl-Göran Pettersson cg.pettersson at evp.slu.se
Fri Oct 7 14:57:22 CEST 2005


Hello all!

Consider a dataset with a grouping structure, Group (factor)
Several treatments, Treat (factor)
Some sort of yield, Yield (numeric)
Something, possibly important, measured for each group; GroupCov (numeric)

To look for fixed effects from Treat on Yield, a first attempt could be:

m1 <- lm(Yield ~ Treat)

which gives, in a symmetric situation, the same estimated fixed effects as:

m2 <- lme(Yield ~ Treat,
          random =~1| Group)

but m2 is a much better model with safer significances.

Now I want to evaluate GroupCov as a covariate to Treat. I can then start
with either m1 or m2 as base, but what is most correct when GroupCov has
only one value for each Group?

m3 <- lm(Yield ~ Treat + GroupCov + Treat:GroupCov)

gives the same fixed effects as

m4 <- lme(Yield ~ Treat + GroupCov + Treat:GroupCov,
          random =~1| Group)

but this time the prob.values for GroupCov is much stronger in m3 than in
m4. Needless to say, anova(m3,m4) tells that m4 is a better *model* than
m3. But is it better for my purpose? Trying an old-fashioned style model
with only fixed effects? (Don´t shout at me, it is only a dirty test of
the system):

m5 <- lm(Yield ~ Group + Treat + GroupCov + Treat:GroupCov)

is accepted by lm() but here GroupCov is silently removed from the
analysis by lm(). I accept this removal, but I get even more suspicious
that something fishy is going on in m4. My gut-feeling is that m3 is the
right starting point but I have got a general recommendation always to
start with m4-type calls when evaluating numeric covariates. I need one
(or two) "doctor second opinion" on this.

How is the mathematics inside lme() working? Is some part of the variation
I want to catch as an effect from GroupCov already removed by the random
call, or why do I get better significances in the pure fixed call? Could
these sigificances be some sort of artefact?

Cheers
/CG

PS.
I sent basically this question, but in a more special case and with
another header, to the list two days ago. Nobody was interested, hopefully
this is more tasty ;-)
DS.


-- 
CG Pettersson, MSci, PhD Stud.
Swedish University of Agricultural Sciences (SLU)
Dep. of Crop Production Ekology. Box 7043.
SE-750 07 Uppsala, Sweden
cg.pettersson at evp.slu.se




More information about the R-help mailing list