[R] fit.contrast - contrast testing with lme

Gang Chen gangchen at mail.nih.gov
Wed Aug 8 16:59:39 CEST 2007


I need some help on contrast testing with lme. It seems fit.contrast  
in 'gmodels' package works fine for simple contrasting among levels  
of a factor such as

 > fit.contrast(fit.lme, "Trust", c(1,-1))
                      Estimate   Std. Error   t-value   Pr(>|t|)
Trust c=( 1 -1 ) -0.001442638 0.0005836959 -2.471557 0.01566063

However if I want to compare some interaction terms, it does not work  
as I wanted:

 > fit.contrast(fit.lme, "Trust:Sex", c(1,-1,0,0))
Error in `contrasts<-`(`*tmp*`, value = c(0.5, -0.5, 0, 0,  
-0.353553390593274,  :
         contrasts apply only to factors

or,

 > fit.contrast(fit.lme, "TrustU:Sex", c(1,-1))
Error in `contrasts<-`(`*tmp*`, value = c(0.5, -0.5)) :
         contrasts apply only to factors

Is there a way I can run such tests? Also how can run an F test for a  
hypothesis like the following with either fit.contrast or lme?

H0: FreqHi - FreqLo = 0, and FreqNo - FreqLo = 0

I also tried to run contrasting directly in lme,  but could not get  
it work:

> anova(fit.lme, L=c("TrustT:Sex:Freq"=1, "TrustU:Sex:Freq"=-1))
Error in anova.lme(fit.lme, L = c("TrustT:Sex:Freq" = 1,  
"TrustU:Sex:Freq" = -1)) :
         Effects TrustT:Sex:Freq, TrustU:Sex:Freq not matched

I'm confused with the error message, and could not figure out what it  
means.

Here is more information:

>> fit.lme <- lme(Beta ~ Trust*Sex*Freq, random = ~1|Subj, Model)
>> summary(fit.lme)

Linear mixed-effects model fit by REML
Data: Model
          AIC       BIC   logLik
    -825.4663 -791.4348 426.7331

Random effects:
Formula: ~1 | Subj
          (Intercept)    Residual
StdDev: 0.001144573 0.001167392

Fixed effects: Beta ~ Trust * Sex * Freq
                             Value    Std.Error DF    t-value p-value
(Intercept)         0.0001090007 0.0005780194 77  0.1885762  0.8509
TrustU              0.0014426378 0.0005836959 77  2.4715572  0.0157
SexM                0.0008230359 0.0005836959 77  1.4100423  0.1626
FreqLo              0.0001998191 0.0005836959 77  0.3423343  0.7330
FreqNo              0.0004900107 0.0005836959 77  0.8394965  0.4038
TrustU:SexM        -0.0012598266 0.0008254707 77 -1.5261918  0.1311
TrustU:FreqLo      -0.0012383346 0.0008254707 77 -1.5001558  0.1377
TrustU:FreqNo      -0.0009141543 0.0008254707 77 -1.1074341  0.2716
SexM:FreqLo        -0.0008469211 0.0008254707 77 -1.0259857  0.3081
SexM:FreqNo         0.0006361012 0.0008254707 77  0.7705922  0.4433
TrustU:SexM:FreqLo  0.0013272173 0.0011673918 77  1.1369082  0.2591
TrustU:SexM:FreqNo  0.0006241524 0.0011673918 77  0.5346555  0.5944
...


> anova(lme(Beta ~ Trust*Sex*Freq, random = ~1|Subj, Model))
>
                 numDF denDF  F-value p-value
(Intercept)        1    77 4.813938  0.0312
Trust              1    77 3.113293  0.0816
Sex                1    77 3.535774  0.0638
Freq               2    77 6.083832  0.0035
Trust:Sex          1    77 1.634858  0.2049
Trust:Freq         2    77 0.678558  0.5104
Sex:Freq           2    77 2.165059  0.1217
Trust:Sex:Freq     2    77 0.647042  0.5264

Your help is highly appreciated.

Thanks,
Gang



More information about the R-help mailing list