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 effects. 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 Luke Duncan *Post-doctoral** Fellow* *School of Animal, Plant and Environmental Sciences* *University of the Witwatersrand* *Johannesburg, South Africa* ** *+27 72 312 0330* *+27 11 717 6452* [[alternative HTML version deleted]]
HI Luke, It would be better to ask this question on R mixed models (r-sig-mixed-models at r-project.org).? Just for curiosity (as I am doing a similar kind of light/dark response in insects), I am interested in the response variable (y) using cbind().? If I understand it correctly, you are using kind of a quasibinomial type of reponse (cbind(counts_in_sun, 12-counts_in_sun)).? I might be wrong.? It would be also great to have a workable small dataset using "dput". A.K. ----- Original Message ----- From: Luke Duncan <luke.mangaliso.duncan at gmail.com> To: r-help at r-project.org Cc: Sent: Saturday, May 26, 2012 6:19 AM Subject: [R] Assessing interaction effects in GLMMs 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 effects. 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 Luke Duncan *Post-doctoral** Fellow* *School of Animal, Plant and Environmental Sciences* *University of the Witwatersrand* *Johannesburg, South Africa* ** *+27 72 312 0330* *+27 11 717 6452* ??? [[alternative HTML version deleted]] ______________________________________________ R-help at 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.
Luke Duncan <luke.mangaliso.duncan <at> gmail.com> writes:> > 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).How many zones are there? It could be a toss-up between treating Zone as random vs fixed ...> 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.OK, be careful (this is a general point about interpreting models in R, not specific to GLMMs); the main effects of Time and Tertile are assessed **at the baseline level of the other** when the interaction is present in the model. You have to be very careful interpreting the meaning of main effects in the presence of interactions. The Type III Wald test you use below probably only makes sense if you have set sum-to-zero contrasts (options(contrasts="sum")).> 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 > effects. Any ideas?Are you trying to test a hypothesis, or just to ask about the values? You can construct the values yourself, or (with the development version of lme4) use the predict() method, to compute what the probabilities will be in each of these situations. By the way, the contrast you have suggested above (10AM, low vs high shade) is not an interaction; to test it, I would suggest that you relevel() your time variable to 10 AM, then look at the estimated effect of Tertile for high vs low.> > 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 *** > ---> > 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 >[snip]