9
$\begingroup$

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

$\endgroup$

migrated from stackoverflow.com Dec 25 '14 at 18:57

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

9
$\begingroup$

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)[1] ## 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).

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

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.