[R-sig-ME] standard errors of estimaters in glmer

Lucy Browning lucybrowning at gmail.com
Tue Feb 2 16:23:33 CET 2010


Hi,

I am trying to calculate standard errors of the means (rather than the
differences between means) of parameter estimates for a model similar
to this one from the lme4 library:

gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 |
herd),family=binomial,data= cbpp)

I am aware that the non-intercept estimates are given as differences
from the intercept and non-intercept standard errors are given as
standard error of the difference between means:

fixef(gm1)
#(Intercept)     period2     period3     period4
 #   -1.3985     -0.9923     -1.1287     -1.5804

se.fixef(gm1)
#(Intercept)     period2     period3     period4
 #    0.2279      0.3054      0.3260      0.4288


Moreover, they are given on the logit scale and so must be
backtransformed.  I have assumed that the standard error of the
difference between the means of two treatment levels (with different
sample sizes) would be the following:

se.diff <- sqrt(variance1/n1 + variance2/n2)

hence, the standard errors of the parameter estimates should be
calculable from the output as follows:

se.estimate <- sqrt( se.diff^2 - se.intercept^2)

where se.intercept is the standard error of the intercept taken from
the output, and se.diff is any one of the standard errors given for
the fixed effects.

So the standard error of the estimate for period2 should be

sqrt(0.3054^2-0.2279^2) ## 0.2033.

To backtransform this, we have to (respectively) add this value to,
and subtract it from, the fixed-effect estimate to obtain upper and
lower standard errors, and then apply an inverse logit tranformation
(in the arm library).

invlogit(-1.3985-0.9923+0.2033) ##  0.1009
invlogit(-1.3985-0.9923-0.2033) ## 0.06952

I am worried that this may not be the correct way of doing it, because
when we simulate the distribution that comes from the model output, we
get something different:

test1<-rnorm(100000,-1.3985,0.2279)
test2<-rnorm(100000,-0.9923,0.3054)
test3<-test1+test2
sd(test3) ##  0.3809

invlogit(-1.3985-0.9923+0.3809)  ## 0.1182
invlogit(-1.3985-0.9923-0.3809)  ## 0.05887

for the model I am using, these differences are larger than this.  Can
anyone advise on the right way to calculate the standard errors of the
means for backtransformed estimates from a binomial glmer?

Many thanks in advance,

Lucy Browning

--




More information about the R-sig-mixed-models mailing list