[R-sig-ME] contrasts in lmer models

[fernando barbero]

Why don t you test those contrasts using the LRT aproach? A full model
against a reduced one with all the effects except the one you are trying to
test. The book "Extending the linear model with R" (by Julian Faraway)
surely will help you
Best regards
Fernando

-----Mensaje original-----
De: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] En nombre de Frei Esther
Enviado el: miércoles, 06 de junio de 2012 03:40 p.m.
Para: mailman, r-sig-mixed-models
Asunto: [R-sig-ME] contrasts in lmer models

Dear R experts,

I am fitting a linear mixed model with lmer. Some of the a-priori-contrasts
I want to test involve interaction terms. I ran into serious trouble when I
tried to test these contrasts:

1. For testing interactions in the fixed factors, I tried using
fit.contrast() from the gmodels-package. In principle this works, but the
results I get are very suspicious. They vary quite a lot if I execute
fit.contrasts several times. I have increased the number of samples from the
default  1000 up to 200,000. However, the results don't seem to converge. It
seems also strange that sometimes p-values are exactly 0. Does it mean that
the contrasts are highly significant (i.e. smaller than the number of
decimals displayed in the result)? Or is there something wrong in the way I
use the function? Maybe I simply have to run more samples with my complex
model?

> require("gmodels")
> ffin <- lmer(y ~ CoV + Alt + Orig + Soil + (1 | Bed) + (1 | Site) + (1 |
Pair),  data = TETri,
+             na.action = na.exclude, REML=TRUE)
> cmat <- rbind("First Contrast (600m vs.1200m)"= c(1,0,-1),"Second Contrast
(1200m vs.1800m)"=
 +                         c(1,-1,0))

> fit.contrast(ffin,"Alt",cmat,conf.int=0.95,sim.mer=TRUE,n.sim=5000)
 
Estimate Std. Error   p value  Lower.CI    Upper.CI
AltFirst Contrast (600m vs.1200m)   -4.879821   0.979278 0.0000000 -6.317718
-3.86426640
AltSecond Contrast (1200m vs.1800m) -2.867170   1.796602 0.3333333 -4.593912
0.03296012


2. Some of the contrasts that I want to test involve interaction terms. As I
understand these cannot directly be approached with fix.contrasts(). Instead
I tried to recode the contrast terms into a new combined variable using
factorA.factorB = paste(factorA,factorB,sep="."). This only works when I
also refit the model with only the new combined factor as fixed effect but
results in suspicious p-values of '0'. As I understand, this is not the same
model anymore and so it is not very useful.

> TETri$O.A.Soil  <- as.factor(paste(TETri$Orig, TETri$Alt, TETri$Soil, sep=
levels(TETri$O.A.Soil) 
> cmat_1 <- rbind("HA Orig12 (O12.T12.SOL vs. O18.T18.GRA)"=
c(0,1,0,0,0,0,0,0,-1,0,0,0),"LF Orig12 +
(O12.T12.SOL vs. O18.T12.SOL)"= c(0,1,0,0,0,0,0,-1,0,0,0,0))
> f1c <- lmer(y ~ CoV + O.A.Soil + (1 | Bed) + (1 | Site) + (1 | Pair),
data = TETri,
+            na.action = na.exclude, REML=TRUE)
> fit.contrast(f1c,"O.A.Soil",cmat_1,conf.int=0.95,sim.mer=TRUE,n.sim=1000)
 
Estimate Std. Error p value   Lower.CI   Upper.CI
O.A.SoilHA Orig12 (O12.T12.SOL vs. O18.T18.GRA) -3.958218   2.110968       0
-6.6922765 -0.7359227
O.A.SoilLF Orig12 (O12.T12.SOL vs. O18.T12.SOL)  2.714403   1.387434       0
0.5868685  5.3193961


3. Using the function estimable() from gmodels, I have trouble specifying
the contrast matrix correctly. The help entry does not really give a good
example how to indicate the factor to which I want to refer to. Could anyone
give me a hint (example) how to define the contrasts correctly, please?
> f0 <- lmer(y ~ CoV + Alt * Orig * Soil + (1 | Bed) + (1 | Site) + (1 |
Pair),  data = TETri,
+             na.action = na.exclude, REML=TRUE)
> cmat <- rbind('Alt600m vs.Alt1200m'= c(1,0,-1),'Alt1200m vs.Alt1800m'=
c(1,-1,0))
> estimable(f0,cmat)
Error in FUN(newX[, i], ...) : 'param' has no names and does not match
number of coefficients of model. Unable to construct coefficient vector


4. My last question is a general question not specifically related to R: If
one of the fixed terms in my mixed model is not significant overall, is it
possible that one of the contrasts between levels of that term is still
significant?

Thanks a lot for your help!
Esther

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
-----
Se certificó que el correo no contiene virus.
Comprobada por AVG - www.avg.es
Versión: 2012.0.2177 / Base de datos de virus: 2433/5053 - Fecha de la
versión: 06/06/2012