Hello all,
I am running a glmer using the lme4 package and the binomial family and am
getting somewhat unexpected results, which I'm hoping someone can help me
make sense of. My data look something like the following:
id group successes total fe1_center fe2_center
1713 A 0 11 -0.0911 -17.2868
1717 A 0 155 -0.0911 -17.2886
2272 B 49 49 -0.0911 -32.2868
2289 B 7 22 -0.2416 -32.2868
1487 B 0 20 0.0537 2.7132
8199 C 10 127 -0.2416 -59.2868
.....
Where my response variable is the proportional of successes. I have
centered the two fixed effects variables to alleviate some problems of
multicollinearity and am also interested in their interaction. The data are
clustered spatially within groups so I am using group as a random/grouping
variable. When running the glmm, the coefficients of the fixed effects and
their interaction seem reasonable (see below). However, when plotting the
predictions vs. the response the curve is consistently lower than i would
expect. E.g., the predicted proportion is lower than the mean proportion of
the data across the full range of data.
#running the model
resp = cbind(data$successes, (data$total - data$successes))
model = glmer( resp ~ fe1_c * fe2_c + (1|group) ,
data=data, family = binomial, REML=F)
summary(model)
Random effects:
Groups Name Variance Std.Dev.
group (Intercept) 4.6432 2.1548
Number of obs: 12271, groups: group, 392
Fixed effects:
Estimate Std. Error z
value Pr(>|z|)
(Intercept) -2.3776295 0.1112830 -21.37
<2e-16 ***
fe1_c -0.8771395 0.0362946 -24.17
<2e-16 ***
fe2_c 0.0109161 0.0001074 101.65
<2e-16 ***
fe1_c:fe2_c -0.0528655 0.0010090 -52.39 <2e-16
***
---
Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
Correlation of Fixed Effects:
(Intr) fe1_c fe2_c
fe1_c -0.012
fe2_c 0.000 -0.687
fe1_c:fe2_c -0.022 0.411 -0.071
I suspect my problem has something to do with how the "average" intercept
is estimated (-2.378). Since I have centered my predictor variables I would
expect the intercept to be equal to the grand mean (is this a correct
assumption?), but in fact it is quite a bit lower.
mean(data$successes/data$total) # equal to 0.2008
logistic (-2.3776295) # equal to 0.0849
Perhaps the model is weighting the unique group intercepts differently
leading to something other than a true average intercept? My group sizes
vary greatly (data comes from messy observations, not experiments) so could
this be affecting the estimate?
Any incite you could give me would be much appreciated. Thank you for the
help.
Zack
--
Zack Steel
Landscape Ecologist
University of California, Davis
zacksteel at gmail.com
[[alternative HTML version deleted]]