thr3ads.net - R help - [R] contrasts for lmer model [May 2010]

If this information is useful, please help other people find it:
Share via:

Kay Cichini

2010-May-19 11:36 UTC

[R] contrasts for lmer model

hello,

i found it most convenient to use package contrast for planned comparisons
on mixed models.
for instance i have a model with 2 fixed factors, one with 4 levels (stage)
and one with 2 levels (gap) and a nested random factor (site) and i tested
gap within level A of factor stage:

library(contrast)
library(nlme)

m1<-lme(rich ~ stage*gap, random=~1|site,data=richness)

contrast(m1,
a = list(stage="A",gap="0"),
b = list(stage="A",gap="1"))


..this works fine. now i wont to do the same with lmer() because the
response is a count - but this doesn't work with contrast():

library(lme4)
detach(package:nlme)

m2<-lmer(rich ~ stage*gap+(1|site),family=poisson,data=richness)
contrast(m2,
a = list(stage="A",gap="0"),
b = list(stage="A",gap="1")) 

Error in UseMethod("contrast") : 
  no applicable method for 'contrast' applied to an object of class
"mer"

i tried with glht() from multcomp, but was not able to set it up right -
could anyone please help me with how this sort of contrast is produced with
glht()?

thanks for any help,
kay

my data:

richness<-data.frame(list(structure(list(stage = structure(c(1L, 1L, 1L, 1L,
1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L), .Label = c("A", "B",
"C", "D"
), class = "factor"), site = structure(c(2L, 1L, 1L, 1L, 1L, 
2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 6L, 6L, 6L, 6L, 
7L, 7L, 7L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 10L, 10L, 10L, 10L, 
11L, 11L, 11L, 11L, 12L, 12L, 12L, 12L, 13L, 13L, 13L, 14L, 14L, 
14L, 14L, 15L, 15L, 15L, 15L, 16L, 16L, 16L, 16L, 17L, 17L, 17L, 
17L, 18L, 18L, 19L, 19L, 19L, 19L, 20L, 20L, 20L, 20L, 21L, 21L, 
21L, 21L, 22L, 22L, 22L, 22L, 23L, 23L, 23L, 24L, 24L, 24L, 24L, 
25L, 25L, 25L, 25L, 26L, 26L, 26L, 26L, 27L, 27L, 27L, 27L, 28L, 
28L, 28L, 28L, 29L, 29L, 29L, 30L, 30L, 30L, 30L, 31L, 31L, 32L, 
32L, 32L, 32L, 33L, 33L, 33L, 33L, 34L, 34L, 34L, 34L, 35L, 35L, 
35L, 35L, 36L, 36L, 36L, 36L, 37L, 37L, 38L, 38L, 38L, 38L, 39L, 
39L, 39L, 39L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
4L, 4L, 4L, 4L, 5L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 
9L, 9L, 9L, 9L, 10L, 10L, 10L, 10L, 11L, 11L, 11L, 11L, 12L, 
12L, 12L, 12L, 13L, 13L, 13L, 14L, 14L, 14L, 14L, 15L, 15L, 15L, 
15L, 16L, 16L, 16L, 16L, 17L, 17L, 17L, 17L, 18L, 18L, 19L, 19L, 
19L, 19L, 20L, 20L, 20L, 20L, 21L, 21L, 21L, 21L, 22L, 22L, 22L, 
22L, 23L, 23L, 23L, 24L, 24L, 24L, 24L, 25L, 25L, 25L, 25L, 26L, 
26L, 26L, 26L, 27L, 27L, 27L, 27L, 28L, 28L, 28L, 28L, 29L, 29L, 
29L, 30L, 30L, 30L, 30L, 31L, 31L, 32L, 32L, 32L, 32L, 33L, 33L, 
33L, 33L, 34L, 34L, 34L, 34L, 35L, 35L, 35L, 35L, 36L, 36L, 36L, 
36L, 37L, 37L, 38L, 38L, 38L, 38L, 39L, 39L, 39L, 39L), .Label =
c("A11",
"A12", "A14", "A15", "A16",
"A17", "A18", "A19", "A20",
"A5",
"A7", "A8", "B1", "B12",
"B13", "B14", "B15", "B17",
"B18", "B2",
"B4", "B7", "B8", "B9", "C1",
"C10", "C11", "C15", "C17",
"C18",
"C19", "C2", "C20", "C3",
"C4", "C6", "D1", "D4", "D7"),
class = "factor"),
    gap = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L), .Label = c("0", "1"), class =
"factor"),
    rich = c(0L, 2L, 3L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 0L, 
    1L, 3L, 3L, 3L, 1L, 2L, 3L, 7L, 4L, 2L, 1L, 3L, 2L, 2L, 2L, 
    2L, 3L, 4L, 1L, 1L, 2L, 1L, 3L, 1L, 1L, 2L, 0L, 0L, 3L, 1L, 
    0L, 1L, 3L, 2L, 3L, 5L, 6L, 7L, 2L, 5L, 8L, 10L, 10L, 9L, 
    3L, 9L, 6L, 2L, 5L, 1L, 6L, 8L, 7L, 7L, 2L, 6L, 7L, 5L, 3L, 
    2L, 3L, 9L, 9L, 9L, 9L, 5L, 7L, 4L, 2L, 8L, 4L, 8L, 5L, 8L, 
    5L, 3L, 7L, 6L, 7L, 8L, 6L, 7L, 5L, 6L, 10L, 10L, 15L, 10L, 
    11L, 5L, 10L, 6L, 7L, 7L, 10L, 7L, 11L, 6L, 4L, 4L, 6L, 5L, 
    6L, 5L, 5L, 8L, 5L, 5L, 7L, 5L, 9L, 8L, 6L, 9L, 3L, 4L, 8L, 
    9L, 7L, 10L, 9L, 9L, 12L, 10L, 12L, 12L, 13L, 11L, 9L, 9L, 
    9L, 0L, 2L, 1L, 4L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 3L, 1L, 7L, 4L, 2L, 4L, 2L, 2L, 3L, 1L, 3L, 1L, 1L, 1L, 
    2L, 2L, 3L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 0L, 3L, 2L, 2L, 3L, 
    1L, 2L, 2L, 4L, 5L, 4L, 5L, 7L, 6L, 4L, 9L, 5L, 2L, 2L, 0L, 
    3L, 6L, 1L, 4L, 3L, 3L, 4L, 2L, 5L, 7L, 4L, 2L, 1L, 1L, 5L, 
    4L, 4L, 4L, 1L, 2L, 1L, 2L, 6L, 5L, 10L, 1L, 4L, 4L, 1L, 
    3L, 5L, 4L, 4L, 4L, 4L, 4L, 3L, 3L, 5L, 9L, 4L, 5L, 6L, 7L, 
    4L, 5L, 4L, 7L, 4L, 6L, 5L, 2L, 3L, 2L, 6L, 1L, 4L, 5L, 6L, 
    3L, 5L, 4L, 4L, 3L, 2L, 4L, 1L, 4L, 3L, 4L, 3L, 6L, 2L, 9L, 
    8L, 10L, 5L, 7L, 10L, 7L, 7L, 4L, 5L, 5L)), .Names = c("stage", 
"site", "gap", "rich"), class =
"data.frame", row.names = c(NA,
-286L))))



-----
------------------------
Kay Cichini
Postgraduate student
Institute of Botany
Univ. of Innsbruck
------------------------

-- 
View this message in context:
http://r.789695.n4.nabble.com/contrasts-for-lmer-model-tp2222682p2222682.html
Sent from the R help mailing list archive at Nabble.com.

Kay Cichini

2010-May-19 12:29 UTC

head link

[R] contrasts for lmer model

.. i found this helpful thread:

http://r.789695.n4.nabble.com/Please-Post-Planned-Contrasts-Example-in-lme-nlme-td1598717.html#a1599031

m0<-lm(rich ~ gap*stage,data=richness)

cc<-contrast(m0,
a = list(stage="A",gap="0"),
b = list(stage="A",gap="1"))

m2<-lmer(rich~gap*stage+(1|site),family=poisson,data=richness)

# plug in matrix from contrast(m0)
summary(glht(m2, linfct = cc$X)) 

ok, this works fine - for each stage seperatly..

but what if i needed all 4 gap effects within each of the stages?
i think that then it may be needed to test "gap" within each stage
simultanously so that the p-values will be adjusted accordingly, isn't it?

is this possible? - if so, how would this be done?
anyone who can help with this..

best,
kay



-----
------------------------
Kay Cichini
Postgraduate student
Institute of Botany
Univ. of Innsbruck
------------------------

-- 
View this message in context:
http://r.789695.n4.nabble.com/contrasts-for-lmer-model-tp2222682p2222743.html
Sent from the R help mailing list archive at Nabble.com.

Kay Cichini

2010-May-19 13:23 UTC

head link

[R] contrasts for lmer model

plugged it together myself -
obviously this should work:

m0<-lm(rich ~ gap*stage,data=richness)

cc<-contrast(m0,
list(stage = levels(stage), gap = "1"),
list(stage = levels(stage), gap = "0"))

# plug in matrix from contrast(m0)
summary(glht(m3, linfct = cc$X)) 


         Simultaneous Tests for General Linear Hypotheses

Fit: glmer(formula = rich ~ gap * stage + (1 | site), family = poisson, 
    REML = 1)

Linear Hypotheses:
       Estimate Std. Error z value Pr(>|z|)    
1 == 0  0.04820    0.15594   0.309   0.9965    
2 == 0 -0.44679    0.10204  -4.379 4.78e-05 ***
3 == 0 -0.52847    0.09133  -5.786 2.88e-08 ***
4 == 0 -0.44396    0.15549  -2.855   0.0171 *  
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 
(Adjusted p values reported -- single-step method)

-----
------------------------
Kay Cichini
Postgraduate student
Institute of Botany
Univ. of Innsbruck
------------------------

-- 
View this message in context:
http://r.789695.n4.nabble.com/contrasts-for-lmer-model-tp2222682p2222820.html
Sent from the R help mailing list archive at Nabble.com.

Possibly Parallel Threads

Search for more possibly parallel threads

R help - May 2010 - contrasts for lmer model

[R] contrasts for lmer model

[R] contrasts for lmer model

[R] contrasts for lmer model

Possibly Parallel Threads