Lesson 13 of 13

Random Slopes and GLMMs

In Lesson 12 you gave each of Dr. Reyes's eight clinics its own baseline recovery level: one random intercept per clinic, all drawn from a shared bell curve. But every clinic was forced to share the same slope for therapy sessions, giving eight parallel lines. That is often too rigid. At some clinics an extra therapy session lifts recovery a lot; at others it barely moves the needle.

This lesson lets the slope vary by clinic too, so the lines fan out instead of running parallel. Then it goes one step further: when the outcome is a yes/no or a count rather than a 0-100 score, the generalized linear mixed model (GLMM) carries everything you know onto that new scale.

By the end of this lesson you will be able to:

  • See why one shared slope (parallel lines) is too rigid, and fit a random-slopes model with lmer
  • Read the three numbers a random slope adds: the spread of intercepts, the spread of slopes, and their correlation
  • Test whether the random slope earns its keep, and watch partial pooling steady the noisy small-clinic slopes
  • Diagnose and fix the two warnings every mixed-model user meets: singular fits and convergence failures
  • Fit a GLMM for a yes/no outcome with glmer and read its effect as an odds ratio

Prerequisites: Lesson 12 on random intercepts (variance components, the ICC, partial pooling), an ordinary linear model, and the logit link and odds for the GLMM half.

The data, rebuilt

One shared slope, eight clinics

We are back with Dr. Anna Reyes and her eight walk-in clinics. Every lesson runs in its own session, so we rebuild her network right here. Alongside each patient's therapy sessions and 0-100 recovery score, we now record two more things: the patient's baseline illness severity (standardised), and whether, at a one-year follow-up, they recovered fully (1) or not (0). We use severity when we reach troubleshooting, and recovered for the GLMM at the end.

RInteractive R
library(lme4) set.seed(3) clinic_names <- c("Ashby","Brook","Cedar","Dale","Elm","Fern","Gale","Hill") n_per <- c(18, 34, 15, 40, 12, 30, 16, 25) # patients seen at each clinic clinic <- factor(rep(clinic_names, n_per)) N <- length(clinic) # 190 patients in all u_int <- setNames(rnorm(8, 0, 7), clinic_names) # each clinic's baseline shift w_slope <- setNames(rnorm(8, 0, 0.7), clinic_names) # each clinic's own extra slope sessions <- pmin(14L, rpois(N, 6)) # therapy sessions attended severity <- round(scale(rnorm(N))[, 1], 2) # standardised illness severity recovery <- 50 + u_int[as.character(clinic)] + (1.4 + w_slope[as.character(clinic)]) * sessions + (-2) * severity + rnorm(N, 0, 6) recovery <- round(pmin(100, pmax(0, recovery)), 1) # keep it on the 0-100 scale v_log <- setNames(rnorm(8, 0, 1.0), clinic_names) # clinic shift on the log-odds eta <- -1.3 + v_log[as.character(clinic)] + 0.25 * sessions recovered <- rbinom(N, 1, plogis(eta)) # full recovery at follow-up? 1/0 recov <- data.frame(clinic, sessions, severity, recovery, recovered) head(recov) #> clinic sessions severity recovery recovered #> 1 Ashby 4 -0.31 50.7 0 #> 2 Ashby 6 -0.44 42.6 0 #> 3 Ashby 4 1.46 40.3 0 #> 4 Ashby 5 -0.73 39.3 0 #> 5 Ashby 9 -0.48 70.3 0 #> 6 Ashby 4 -0.04 38.3 0

  

The question that drives this lesson: does an extra therapy session help every clinic by the same amount? Lesson 12's random-intercept model assumed yes. Let us look.