[R-sig-eco] question on using ghlt (multcomp) for mixed models with lme (nlme)

[Hem Nalini Morzaria Luna]

Hi:

I am comparing biomass response in a an experiment with  three factors
(N-two levels, W-three levels, and S-two levels) using lme. Unit is
just specified as a block effect. (Data example below)

This is my model

> w.lme <- lme(fixed=log(Ry.above) ~ W * S * N, data= Mono, random=~1 | Unit/W/N)

I'd like to do multiple comparisons.
In this case, interaction terms of W:Species is slightly significant
(p 0.07) so if I want to only test for W within each Species.

I tried the following code to get a contrast matrix for those factors

> library(contrast)
> contrast(w.lme, list(W=c('High','Medium','Low'), Species = 'Tc'), list(W=c('High','Medium','Low'), Group = 'Sv'), type='average')

But I get

>Error in gendata.default(fit = list(modelStruct = list(reStruct = list( :
  not enough factors

I was also trying to understand what the correct way to specify
multiple comparions for the whole model is:

One factor is straightforward
> comp.Species <- glht(w.lme, linfct=mcp(Species="Tukey"))

To specify all comparisons I know mcp will not take interactions but
the linfct argument can use a matrix of comparisons but I am unclear
what is the correct approach to set up that matrix.

I tried:

1) matrix from fixed effects

> K <- diag(length(fixef(w.lme)))
> glht(w.lme, linfct=K)

2) matrix using contrasts

> matrix.con = contrasts(Mono$S:Mono$N:Mono$W)
#matrix dimension is [11,12]
#so I added one more column, results in this matrix being equal to K(above)
> missing.column = matrix(0,12,1)
> missing.column[1,1] = 1
con.matrix = cbind(missing.column,matrix.con)
>glht(w.lmer, linfct = con.matrix)

1 & 2 give the same result from glht. Is it correct to assume the rows
in K are the same as in the contrast matrix?

3) Model matrix using model.matrix
> lme.matrix <- model.matrix(w.lme, data=Mono)
glht(w.lmer, linfct = lme.matrix)

In 3) lme.matrix dimension is [36,12] and the results of glht repeat
three times each. I am having a hard time understanding why this is
and what the rows and columns in the matrix represent. Colnames are an
incomplete set of contrasts.

What is the correct approach to specify all comparisons?


Any help is appreciated. Thank you!

Hem Nalini Morzaria Luna

Data (Mono)
      N      W Unit Species  Ry.above
1  High   High    1      Sv 0.6156346
2  High   High    2      Sv 0.6235052
3  High   High    3      Sv 0.6985453
4   Low   High    4      Sv 1.7452381
5   Low   High    5      Sv 0.8424228
6   Low   High    6      Sv 0.6643357
7  High Medium    7      Sv 0.8251745
8  High Medium    8      Sv 0.7424867
9  High Medium    9      Sv 0.9362663
10  Low Medium   10      Sv 1.0038536
11  Low Medium   11      Sv 0.7767575
12  Low Medium   12      Sv 1.5007492
13 High    Low   13      Sv 0.9536121
14 High    Low   14      Sv 0.9799515
15 High    Low   15      Sv 0.5596921
16  Low    Low   16      Sv 0.8706301
17  Low    Low   17      Sv 1.0032415
18  Low    Low   18      Sv 1.3620130
19 High   High    1      Tc 0.6514929
20 High   High    2      Tc 0.4309266
21 High   High    3      Tc 0.2800695
22  Low   High    4      Tc 0.7883481
23  Low   High    5      Tc 0.4541833
24  Low   High    6      Tc 0.6007958
25 High Medium    7      Tc 0.6987778
26 High Medium    8      Tc 0.2240063
27 High Medium    9      Tc 0.1477910
28  Low Medium   10      Tc 0.2017751
29  Low Medium   11      Tc 0.3621421
30  Low Medium   12      Tc 0.2590503
31 High    Low   13      Tc 0.4265925
32 High    Low   14      Tc 0.5525204
33 High    Low   15      Tc 0.4196785
34  Low    Low   16      Tc 0.4507284
35  Low    Low   17      Tc 1.1027065
36  Low    Low   18      Tc 0.3676259