Dear R experts
I am running a negative-binomial GLM (glm.nb) to test the null hypotheses
that species 1 and 2 are equally abundant between site 1 and site2, and
between each other. So, I have a 2x2 factorial design with factors Site
(1,2) and Taxon (1,2).
Since the Site:Taxon interaction is significant, I need to do the equivalent
to a "post-hoc test" for ANOVA, however, the same tests (e.g. Tukey HSD) do
not seem to be applicable for the GLM. I tried specifying orthogonal
contrasts, but could not figure out what the interaction contrast (see
Site1:Taxon1 in below example) means.
Could you please advise me how to specify a meaningful interaction contrast
(i.e. contrast species within sites)? Alternatively, could you recommend a
way to do posthoc comparisons?
Thanks for your time and kind regards
Maya
> library(MASS)
> counts <- c(1, 4, 9, 2, 1, 4, 2, 4, 1, 3, 2, 2, 1, 3, 1, 1, 2, 1, 113, 83,
49, 46, 13, 52, 4, 10, 14, 10, 3, 19, 8, 21, 151, 186, 99, 11, 29, 24, 24,
62, 15, 98, 30, 21, 63, 29, 48, 11, 16, 35, 21, 17, 6, 2, 2, 3, 3, 4, 4, 2,
1, 2, 2, 3, 4, 8, 10, 3, 14, 3, 11, 23, 3, 51, 8, 8, 7, 1, 13, 8, 1, 0, 1,
0, 1, 1, 1, 0, 1, 1, 1, 5, 8, 1, 1, 20, 2, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0,
0, 3, 0, 1, 0, 5, 1, 1, 9, 0, 34, 4, 1, 17, 0, 7, 33, 86, 73, 67, 79, 109,
27, 37, 23, 12, 17, 41, 8, 38, 4, 23, 14, 49, 64, 39, 31, 156, 110, 97, 33,
170, 137, 72, 28, 54)
> Site <- factor(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2))
> Taxon <- factor(c(2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1))
> contrasts(Site) <- cbind(c(1,-1))
> contrasts(Taxon) <- cbind(c(1,-1))
> s1 <- glm.nb(counts ~ Site*Taxon)
> summary(s1)
Call:
glm.nb(formula = counts ~ Site * Taxon, init.theta = 0.612672617555492,
link = log)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.269021 -1.215727 -0.454719 0.003515 2.517288
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.6139 0.1097 23.831 < 2e-16 ***
Site1 -0.8664 0.1097 -7.899 2.81e-15 ***
Taxon1 -0.3814 0.1097 -3.477 0.000506 ***
Site1:Taxon1 -0.5977 0.1097 -5.449 5.06e-08 ***
---
Signif. codes: 0 â***â 0.001 â**â 0.01 â*â 0.05 â.â 0.1 â
â 1
(Dispersion parameter for Negative Binomial(0.6127) family taken to be 1)
Null deviance: 242.94 on 153 degrees of freedom
Residual deviance: 176.20 on 150 degrees of freedom
AIC: 1165.4
Number of Fisher Scoring iterations: 1
Theta: 0.6127
Std. Err.: 0.0690
2 x log-likelihood: -1155.4010
>
> TukeyHSD(s1)
Error in UseMethod("TukeyHSD") : no applicable method for "TukeyHSD"
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
