[R-sig-ME] Assessing interaction effects for GLMM

Tue May 29 10:58:11 CEST 2012

Dear Luke,

A workaround is to first create a new variable for the interaction.

TimeTertile <- interaction(Time, Tertile, drop = TRUE)

And then use that variable in the model

m1 <- glmer(y ~ TimeTertile + (1|Individual) + (1|Zone) + Max + Min, family=binomial)

Please not that I use glmer() instead of lmer() which give some more clear code. At the moment, lmer() will call glmer() for you when family is not Gaussian.
Then you can run multicomp

Timetest <- glht(m1, linfct = mcp(TimeTertile = "Tukey"))

Best regards,

Thierry

ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium
+ 32 2 525 02 51
+ 32 54 43 61 85
Thierry.Onkelinx op inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey

-----Oorspronkelijk bericht-----
Van: r-sig-mixed-models-bounces op r-project.org [mailto:r-sig-mixed-models-bounces op r-project.org] Namens Luke Duncan
Verzonden: maandag 28 mei 2012 9:36
Aan: r-sig-mixed-models op r-project.org
Onderwerp: [R-sig-ME] Assessing interaction effects for GLMM

Dear R gurus

I am running a GLMM that looks at whether chimpanzees spend time in shade more than sun (response variable 'y': used cbind() on counts in the sun and
shade) based on the time of day (Time) and the availability of shade (Tertile). I've included some random factors too which are the chimpanzee in question (Individual) and where they are in a given area (Zone). There are also two continuous predictors (Minimum daily temperature: Min; Maximum daily temperature: Max). I have run my GLMM and I know that Time and Min are significant predictors of the patterns of shade use while Tertile and Max are not. In addition, a Time*Tertile interaction effect is a good predictor as well.

I now need to assess how the specific interaction effect conditions differ to one another. So, for example, how does shade use differ between 10h00 at low shade and 10h00 at high shade? I tried using the package multcomp, but that will only allow me to work out the contrasts for the first-order effects (Time, Tertile) but won't allow me to do so for the interaction effect (Time*Tertile). Any ideas?

My code:

> m1 <- lmer(y ~ Time*Tertile + (1|Individual) + (1|Zone) + Max +
Min,family=binomial,REML=F)
> Anova(m1,type=3,test="Wald")
Analysis of Deviance Table (Type III tests)

Response: y
               Chisq Df Pr(>Chisq)
(Intercept)   0.9511  1     0.3294
Time         60.7807  4  1.988e-12 ***
Tertile       0.3391  1     0.5603
Max           1.3198  1     0.2506
Min          77.7736  1  < 2.2e-16 ***
Time:Tertile 38.9038  4  7.292e-08 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> summary(m1)
Generalized linear mixed model fit by the Laplace approximation
Formula: y ~ Time * Tertile + (1 | Individual) + (1 | Zone) + Max + Min
  AIC  BIC logLik deviance
 1168 1224 -569.9     1140
Random effects:
 Groups     Name        Variance Std.Dev.
 Zone       (Intercept) 0.81949  0.90526
 Individual (Intercept) 0.36417  0.60347 Number of obs: 412, groups: Zone, 8; Individual, 7

 Fixed effects:
                     Estimate Std. Error z value Pr(>|z|)
(Intercept)           0.77498    0.79465   0.975 0.329439
Time11h00            -1.54259    0.24351  -6.335 2.38e-10 ***
Time12h00             0.01695    0.77829   0.022 0.982627
Time13h00            -4.26913    0.78217  -5.458 4.81e-08 ***
Time14h00            -1.34503    0.43831  -3.069 0.002150 **
TertileLow            0.32614    0.56003   0.582 0.560323
Max                   0.03751    0.03265   1.149 0.250630
Min                  -0.30912    0.03505  -8.819  < 2e-16 ***
Time11h00:TertileLow  1.03079    0.28579   3.607 0.000310 ***
Time12h00:TertileLow -2.26187    0.79930  -2.830 0.004658 **
Time13h00:TertileLow  2.38129    0.79214   3.006 0.002646 **
Time14h00:TertileLow  1.72263    0.49397   3.487 0.000488 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) Tm1100 Tm1200 Tm1300 Tm1400 TrtlLw Max    Min    T1100:
Time11h00   -0.026
Time12h00   -0.035  0.177
Time13h00   -0.004  0.223  0.068
Time14h00   -0.073  0.259  0.081  0.103
TertileLow  -0.450  0.153  0.043  0.051  0.097
Max         -0.711 -0.169 -0.004 -0.061 -0.023  0.019
Min          0.146  0.186  0.014  0.055  0.099 -0.036 -0.455
Tm11h00:TrL  0.059 -0.851 -0.153 -0.190 -0.222 -0.198  0.096 -0.155 Tm12h00:TrL  0.095 -0.160 -0.974 -0.062 -0.081 -0.067 -0.079  0.012  0.192 Tm13h00:TrL  0.026 -0.208 -0.067 -0.983 -0.099 -0.075  0.024 -0.026  0.229 Tm14h00:TrL  0.126 -0.215 -0.069 -0.088 -0.876 -0.185 -0.047  0.006  0.254
            T1200: T1300:
Time11h00
Time12h00
Time13h00
Time14h00
TertileLow
Max
Min
Tm11h00:TrL
Tm12h00:TrL
Tm13h00:TrL  0.081
Tm14h00:TrL  0.098  0.116

> Timetest <- glht(m1, linfct = mcp(Time = "Tukey"))
Warning message:
In mcp2matrix(model, linfct = linfct) :
  covariate interactions found -- default contrast might be inappropriate
> Timecld<-cld(Timetest)
> Timecld
11h00 12h00 13h00 14h00 10h00
  "a"  "ab"   "c"   "a"   "b"
> summary(Timetest)

         Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: glmer(formula = y ~ Time * Tertile + (1 | Individual) + (1 |
    Zone) + Max + Min, family = binomial, REML = F)

Linear Hypotheses:
                   Estimate Std. Error z value Pr(>|z|)
11h00 - 10h00 == 0 -1.54259    0.24351  -6.335  < 0.001 ***
12h00 - 10h00 == 0  0.01695    0.77829   0.022  1.00000
13h00 - 10h00 == 0 -4.26913    0.78217  -5.458  < 0.001 ***
14h00 - 10h00 == 0 -1.34503    0.43831  -3.069  0.01497 *
12h00 - 11h00 == 0  1.55954    0.77322   2.017  0.22656
13h00 - 11h00 == 0 -2.72654    0.76570  -3.561  0.00285 **
14h00 - 11h00 == 0  0.19756    0.44278   0.446  0.99024
13h00 - 12h00 == 0 -4.28608    1.06546  -4.023  < 0.001 ***
14h00 - 12h00 == 0 -1.36198    0.86161  -1.581  0.47037
14h00 - 13h00 == 0  2.92410    0.85641   3.414  0.00464 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 (Adjusted p values reported -- single-step method)

Luke Duncan

*Post-doctoral** Fellow*
*School of Animal, Plant and Environmental Sciences* *University of the Witwatersrand* *Johannesburg, South Africa*
**
*+27 11 717 6452*

        [[alternative HTML version deleted]]

* * * * * * * * * * * * * D I S C L A I M E R * * * * * * * * * * * * *
Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document.
The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document.