How do I interpret the variance of random effect in a generalized linear mixed model

In a logistic Generalized Linear Mixed Model (family = binomial), I don't know how to interpret the random effects variance:

Random effects:
Groups   Name        Variance Std.Dev.
HOSPITAL (Intercept) 0.4295   0.6554
Number of obs: 2275, groups: HOSPITAL, 14

How do I interpret this numerical result?

I have a sample of renal trasplanted patients in a multicenter study. I was testing if the probability of a patient being treated with a specific antihypertensive treatment is the same among centers. The proportion of patients treated varies greatly between centers, but may be due to differences in basal characteristics of the patients. So I estimated a generalized linear mixed model (logistic), adjusting for the principal features of the patiens. This are the results:

Generalized linear mixed model fit by maximum likelihood ['glmerMod']
Family: binomial ( logit )
Formula: HTATTO ~ AGE + SEX + BMI + INMUNOTTO + log(SCR) + log(PROTEINUR) + (1 | CENTER)
Data: DATOS

AIC      BIC   logLik deviance
1815.888 1867.456 -898.944 1797.888

Random effects:
Groups   Name        Variance Std.Dev.
CENTER (Intercept) 0.4295   0.6554
Number of obs: 2275, groups: HOSPITAL, 14

Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept)               -1.804469   0.216661  -8.329  < 2e-16 ***
AGE                       -0.007282   0.004773  -1.526  0.12712
SEXFemale                 -0.127849   0.134732  -0.949  0.34267
BMI                        0.015358   0.014521   1.058  0.29021
INMUNOTTOB                 0.031134   0.142988   0.218  0.82763
INMUNOTTOC                -0.152468   0.317454  -0.480  0.63102
log(SCR)                   0.001744   0.195482   0.009  0.99288
log(PROTEINUR)             0.253084   0.088111   2.872  0.00407 **

The quantitative variables are centered. I know that the among-hospital standard deviation of the intercept is 0.6554, in log-odds scale. Because the intercept is -1.804469, in log-odds scale, then probability of being treated with the antihypertensive of a man, of average age, with average value in all variables and inmuno treatment A, for an "average" center, is 14.1 %. And now begins the interpretation: under the assumption that the random effects follow a normal distribution, we would expect approximately 95% of centers to have a value within 2 standard deviations of the mean of zero, so the probability of being treated for the average man will vary between centers with coverage interval of:

exp(-1.804469-2*0.6554)/(1+exp(-1.804469-2*0.6554))

exp(-1.804469+2*0.6554)/(1+exp(-1.804469+2*0.6554))

Is this correct?

Also, how can I test in glmer if the variability between centers is statistically significant? I used to work with MIXNO, an excellent software of Donald Hedeker, and there I have an standard error of the estimate variance, that I don't have in glmer. How can I have the probability of being treated for the "average" man in each center, with a confidene interval?

Thanks

migrated from stackoverflow.comDec 25 '14 at 18:57

This question came from our site for professional and enthusiast programmers.

It's probably most helpful if you show us more information about your model, but: the baseline value of the log-odds of whatever your response is (e.g. mortality) varies across hospitals. The baseline value (the per-hospital intercept term) is the log-odds of mortality (or whatever) in the baseline category (e.g. "untreated"), at a zero value of any continuous predictors. That variation is assumed to be Normally distributed, on the log-odds scale. The among-hospital standard deviation of the intercept is 0.6554; the variance (just the standard deviation squared -- not a measure of the uncertainty of the standard deviation) is $0.6554^2=0.4295$.

(If you clarify your question/add more detail about your model I can try to say more.)

update: your interpretation of the variation seems correct. More precisely,

cc <- fixef(fitted_model) ## intercept
ss <- sqrt(unlist(VarCorr(fitted_model))) ## random effects SD
plogis(qnorm(c(0.025,0.975),mean=cc,sd=ss))

should give you the 95% interval (not really quite confidence intervals, but very similar) for the probabilities of a baseline (male/average age/etc.) individual getting treated across hospitals.

For testing the significance of the random effect, you have a variety of choices (see http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html for more information). (Note that the standard error of a RE variance is usually not a reliable way to test significance, since the sampling distribution is often skewed/non-Normal.) The simplest approach is to do a likelihood ratio test, e.g.

pchisq(2*(logLik(fitted_model)-logLik(fitted_model_without_RE)),
df=1,lower.tail=FALSE)/2

The final division by 2 corrects for the fact that the likelihood ratio test is conservative when the null value (i.e. RE variance=0) is on the boundary of the feasible space (i.e. the RE variance cannot be <0).

• I edited the question – user2310909 Dec 26 '14 at 21:05
• Thanks, it's a pleasure to read to someone who knows well about the subject. – user2310909 Dec 27 '14 at 11:58