I am new to Cross Validated so please forgive me if this question has been asked before. However, I did not find any post that answered my question, so here it is:
I am running a 3 level multilevel binary logistic regression (one binary outcome variable and one binary predictor variable) with 839 observations nested in 171 study participants nested in 29 groups. I am using the glmer() function of the lme4 package in R. When I am specifying the empty model and testing it against a “normal” logistic regression without a random intercept for groups and participants, the results clearly tell me that my data is clustered at the participant level and that I do need to use multilevel modeling.
Models: M0_simple: Outcome ~ 1 M0: Outcome ~ (1 | Group/Person) Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq) M0_simple 1 975.15 979.87 -486.57 973.15 M0 3 831.07 845.25 -412.54 825.07 148.07 2 < 2.2e-16 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Moreover, when I look at the estimates of the empty multilevel model, the random intercept variance at the participant level (level 2) is very high. And when I calculated the VPC for the participant level, the result of .975 is also extremely high.
Random effects: Groups Name Variance Std.Dev. Person:Group (Intercept) 127.4 11.29 Group (Intercept) 0.0 0.00 Number of obs: 839, groups: Person:Group, 171; Group, 29 Fixed effects: Estimate Std. Error z value Pr(>|z|) (Intercept) 8.8693 0.8254 10.74 <2e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
These are the results of the random intercept model once I put in my predictor variable:
M1 <- glmer(Outcome ~ Predictor + (1 | Group/Person), family = binomial("logit"), data = data_M) Random effects: Groups Name Variance Std.Dev. Person:Group (Intercept) 129.1 11.36 Group (Intercept) 0.0 0.00 Number of obs: 839, groups: Person:Group, 171; Group, 29 Fixed effects: Estimate Std. Error z value Pr(>|z|) (Intercept) 9.5390 0.9430 10.116 <2e-16 *** Predictor -0.8154 0.4899 -1.664 0.0961 . --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
How should I interpret these huge random intercept variances at the participant level? I realize that the data is obviously strongly clustered at the participant level. Is this variance so high because the level 1 variance is fixed to 3.29 for multilevel logistic regressions? And does this large variance also affect the fixed effects?
I tried to calculate the predicted probabilities for the random intercept level and ended up with 99.9 %. Moreover, odds-ratios for the intercept of 13891.05 do seem weird. Did I misspecify the model somehow or what might be the issue here? When I run the same model with SPSS 23 it gives out much more reasonable results:
GENLINMIXED /DATA_STRUCTURE SUBJECTS=Group*Person /FIELDS TARGET=Outcome TRIALS=NONE OFFSET=NONE /TARGET_OPTIONS DISTRIBUTION=BINOMIAL LINK=LOGIT /FIXED EFFECTS=Predictor USE_INTERCEPT=TRUE /RANDOM USE_INTERCEPT=TRUE SUBJECTS=Group COVARIANCE_TYPE=VARIANCE_COMPONENTS /RANDOM USE_INTERCEPT=TRUE SUBJECTS=Group*Person COVARIANCE_TYPE=VARIANCE_COMPONENTS /BUILD_OPTIONS TARGET_CATEGORY_ORDER=DESCENDING INPUTS_CATEGORY_ORDER=DESCENDING MAX_ITERATIONS=100 CONFIDENCE_LEVEL=95 DF_METHOD=RESIDUAL COVB=ROBUST PCONVERGE=0.000001(ABSOLUTE) SCORING=0 SINGULAR=0.000000000001 /EMMEANS_OPTIONS SCALE=ORIGINAL PADJUST=LSD. Random effects: Groups Name Variance Std.Error Person:Group (Intercept) 5.052 .831 Group (Intercept) 0.0 … Number of obs: 839, groups: Person:Group, 171; Group, 29 Fixed effects: Estimate Std. Error z value Pr(>|z|) (Intercept) 2.571 0.2413 10.652 .000 Predictor -0.445 0.2400 -1.853 .064 ---
While the p-values do seem in the same ball park, the results from SPSS make much more sense to interpret. When I calculate the predicted probabilities from the fixed effects of the intercept and the predictor I get 89.3 % and 92.9 % and the odds-ratio for the intercept of 13.076 seem much more likely than 13891.05.
So what it comes down to, I guess, are the following questions:
I have read that if I want to use likelihood ratio tests to determine the significance of a predictor, I have to use a statistical program that uses Maximum Likelihood (ML) and not Restricted Maximum Likelihood (REML). This is why I use glmer() in R.
However, once I have established that a certain model (with the predictor) is correct or not, how do I interpret the results? Can I simply look at the estimates and interpret the odds-ratio and calculate VPCs and predicted probabilities? Or is this susceptible to mistakes, since the variance at level 1 for multilevel logistic regressions is fixed at 3.29 and the “higher” random variances are scaled accordingly?
Am I even allowed to calculate predicted probabilities from the random intercept model and why are the results of SPSS and R so different?