[R-sig-ME] Low intercept estimate in a binomial glmm

[Jarrod Hadfield]

Hi,

plogis(-2.3776295) is the mode not the mean.

An approximation for the mean is:

c2<-((16*sqrt(3))/(15*pi))^2

plogis(-2.3776295/sqrt(1+c2*4.6432))

and this should be closer to the observed mean.

Cheers,

Jarrod



Quoting Zack Steel <zacksteel at gmail.com> on Wed, 3 Apr 2013 12:58:46 -0700:

> 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]]
>
>



-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.