[R-sig-ME] Plotting results of interaction from GLMM using model summary coefficients
ONKELINX, Thierry
Thierry.ONKELINX at inbo.be
Wed May 16 15:06:38 CEST 2012
Dear Roslyn,
It seems to me like you are confusing coefficients and predictions. How did you make the plots? Did you back-transform the data of the glmer with poisson? Note that the poisson family uses the lok-link by default.
A reproducible example would ben ice.
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 Roslyn Anderson
Verzonden: zaterdag 12 mei 2012 19:58
Aan: r-sig-mixed-models op r-project.org
Onderwerp: [R-sig-ME] Plotting results of interaction from GLMM using model summary coefficients
Dear List Members,
I have a query regarding the plotting of results from a GLMM model.
I am working with both linear mixed effects models (LMMs) and generalised linear mixed models (GLMMs), mainly using the nlme and lme4 packages in R (2.14.2).
The study involved looking at the effects of *altitude* and the *percentage of traditional sheep on a farm* on various measures of plant diversity.
Here I will concentrate on the response variable *plant species richness*(which ranges from 4 - 26). It is count data so I have used a Poisson distribution.
The experimental design had a nested structure, so I chose to use LMMs and GLMMs. Model syntax for plant species richness is:
*SR <- glmer (plant.spp.rich ~ log10(alt) * TRADEWE.cat + (1|site/farm/habitat), family = poisson, data=veg)*
Fixed effects are:
1) log10(alt) - a log10 transformation on altitude (m) (range: 1.56 - 2.78) and
2) TRADEWE.cat - the percentage of traditional sheep on a farm split into four categories (1 = 0%, 2 = 50%, 3 = 70%, 4 = 100%). Groups are of unequal size, therefore data is unbalanced.
Random effects are:
habitats (categories 1-12), which are nested within farms (categories 1-12), which are nested within sites (categories 1-4).
After rigorous model simplification and validation I am happy with the model summary output and the interaction between log10(alt) and
TRADEWE.cat1 is significantly different to the interaction between
log10(alt) and TRADEWE.cat4. Therefore I wanted to plot this result with plant species richness on the y-axis, log10(alt) on the x-axis and each of the four traditional sheep categories as regression lines on the plot.
See summary below from the GLMM model and then the output of a linear regression model containing the same variables:
> *library(lme4)*
> *SR<-glmer(plant.spp.rich ~ log10(alt) * TRADEWE.cat +
(1|site/farm/habitat), family=poisson, data=veg)*
> *summary(SR)*
Generalized linear mixed model fit by the Laplace approximation
Formula: plant.spp.rich ~ log10(alt) * TRADEWE.cat + (1 |
site/farm/habitat)
Data: veg
AIC BIC logLik deviance
146.3 178.9 -62.13 124.3
Random effects:
Groups Name Variance Std.Dev.
habitat:(farm:site) (Intercept) 0.028298 0.16822
farm:site (Intercept) 0.000000 0.00000
site (Intercept) 0.000000 0.00000
Number of obs: 144, groups: habitat:(farm:site), 48; farm:site, 12; site, 4
Fixed effects:
Estimate Std. Error z value
Pr(>|z|)
(Intercept) 2.0134 0.5203 3.870
0.000109 ***
log10(alt) 0.1717 0.2365 0.726
0.467844
TRADEWE.cat2 -1.3003 1.7218 -0.755 0.450117
TRADEWE.cat3 0.5645 1.5340 0.368 0.712872
TRADEWE.cat4 1.6082 0.6260 2.569 0.010200 *
log10(alt):TRADEWE.cat2 0.5322 0.6907 0.771 0.440985
log10(alt):TRADEWE.cat3 -0.2772 0.7047 -0.393 0.694059
log10(alt):TRADEWE.cat4 -0.6733 0.2814 -2.393 0.016720 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) lg10() TRADEWE.2 TRADEWE.3 TRADEWE.4 l10():TRADEWE.2
l10():TRADEWE.3
log10(alt) -0.991
TRADEWE.ct2 -0.302 0.299
TRADEWE.ct3 -0.339 0.336 0.102
TRADEWE.ct4 -0.831 0.824 0.251 0.282
l10():TRADEWE.2 0.339 -0.342 -0.997 -0.115 -0.282
l10():TRADEWE.3 0.333 -0.336 -0.100 -0.996 -0.276 0.115
l10():TRADEWE.4 0.833 -0.841 -0.252 -0.282 -0.991 0.288
0.282
> *SRa<-lm(plant.spp.rich ~ log10(alt) * TRADEWE.cat, data=veg)*
> *summary(SRa)*
Call:
lm(formula = plant.spp.rich ~ log10(alt) * TRADEWE.cat, data = veg)
Residuals:
Min 1Q Median 3Q Max
-7.955 -2.091 -0.221 1.943 12.404
Coefficients:
Estimate Std. Error t value
Pr(>|t|)
(Intercept) 7.092 4.421 1.604
0.11096
log10(alt) 1.837 2.017 0.911
0.36397
TRADEWE.cat2 -17.959 15.078 -1.191 0.23570
TRADEWE.cat3 6.066 12.795 0.474 0.63618
TRADEWE.cat4 19.437 5.365 3.623 0.00041 ***
log10(alt):TRADEWE.cat2 7.318 6.058 1.208 0.22913
log10(alt):TRADEWE.cat3 -3.025 5.875 -0.515 0.60747
log10(alt):TRADEWE.cat4 -8.101 2.409 -3.363 0.00100 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 3.54 on 136 degrees of freedom
Multiple R-squared: 0.1819, Adjusted R-squared: 0.1398
F-statistic: 4.319 on 7 and 136 DF, p-value: 0.0002399
My question is why, when using a linear regression model, my coefficients (intercept and slope) for each traditional sheep category fit nice regression lines within the data points on the plot but when using an LMM or a GLMM the coefficients are much lower (and therefore the regression lines are also much lower) than the actual data points on the plot? Is this related to the random effects?
I was just wondering if anyone could suggest an alternative method I might use to plot this result using the coefficients from the GLMM output?
Many, many thanks in advance for any help, advice or guidance.
Best Wishes,
Roz
[[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.
More information about the R-sig-mixed-models
mailing list