I have collected data on gas fluxes from plots of soil subjected to 5 different treatments ("D2", "K2", "M", "N", and "O2"), which also possessed variable clay contents. The experiment was laid out in a randomized complete block design, with 4 replications. Within each plot, two separate measurements of flux were performed. 4 blocks x 5 treatments x 2 repeated measurements for 40 total observations
I have built a linear mixed model with nlme, and performed multiple comparisons with Tukey's test, as shown below. The data.frame cum_fluxes can be downloaded from: https://www.dropbox.com/s/e58k4vdnevw2fsm/cum_fluxes
require(multcomp)
require(nlme)
cum_fluxes <- dget("cum_fluxes")
n2o_flux.lme <- lme(cum_flux_n2o_ln ~ treatment * clay + block, random = ~1|subsampling, data = cum_fluxes)
aov.n2o <- anova(n2o_flux.lme)
comp.treat <- glht(n2o_flux.lme, linfct=mcp(treatment="Tukey"))
print(summary(comp.treat))
I have two questions. On the one hand, the resulting ANOVA table looks like this (significant effect of treatment with p = 0.0279):
numDF denDF F-value p-value
(Intercept) 1 20 14815.691 <.0001
treatment 4 7 5.281 0.0279
clay 1 7 1.326 0.2874
block 3 7 0.674 0.5951
treatment:clay 4 7 0.810 0.5568
However, there are apparently no significant differences found when considering pair-wise comparisons of the five treatments. Did I make a mistake somewhere, or am I interpreting the results incorrectly?
Second, I get the following warning message when I run the glht() command:
Warning message:
In mcp2matrix(model, linfct = linfct) :
covariate interactions found -- default contrast might be inappropriate
This makes some sense to me: since I allow for interaction between a covariate (quantitative regressor – i.e., clay content) and a factor (so, non-parallel lines), the difference between two factor levels depends on the value of the covariate. However, I can't figure out how and to what value I should specify the covariate when I compare the factor levels?
