Summarize a Mixed Model

Intro

This is pretty printing of results of a mixed model. Correlations of random and fixed effects are turned off by default, as the focus of interest is usually the variance components and fixed effects themselves. 95% confidence intervals are provided by default.

Examples

library(lme4)
library(glmmTMB)
library(nlme)
library(brms)
library(mgcv)


lmer_model <-
  lmer(Reaction ~ Days + (1 + Days | Subject), data = sleepstudy)

lme_model  <-
  lme(Reaction ~ Days, random = ~ 1 + Days | Subject, data = sleepstudy)

tmb_model  <-
  glmmTMB(Reaction ~ Days + (1 + Days | Subject), data = sleepstudy)

# brms_model <-
#   brm(Reaction ~ Days + (1 + Days | Subject), 
#       data = sleepstudy,
#       cores = 4,
#       refresh = -1,
#       verbose = FALSE
#   )

# this is akin to (1 | Subject) + (0 + Days | Subject) in lme4
mgcv_model <-
  gam(
    Reaction ~  Days +
      s(Subject, bs = 're') +
      s(Days, Subject, bs = 're'),
    data   = lme4::sleepstudy,
    method = 'REML'
  )

library(mixedup)

summarize_model(lmer_model)
Computing profile confidence intervals ...

Variance Components:
    Group    Effect Variance    SD SD_2.5 SD_97.5 Var_prop
  Subject Intercept   612.10 24.74  14.38   37.72     0.47
  Subject      Days    35.07  5.92   3.80    8.75     0.03
 Residual             654.94 25.59  22.90   28.86     0.50

Fixed Effects:
      Term  Value   SE     t P_value Lower_2.5 Upper_97.5
 Intercept 251.41 6.82 36.84    0.00    238.03     264.78
      Days  10.47 1.55  6.77    0.00      7.44      13.50


summarize_model(lme_model)

Variance Components:
    Group    Effect Variance    SD SD_2.5 SD_97.5 Var_prop
  Subject Intercept   612.08 24.74  15.52   39.44     0.47
  Subject      Days    35.07  5.92   3.90    8.98     0.03
 Residual        NA   654.94 25.59  22.80   28.73     0.50

Fixed Effects:
      Term  Value   SE     Z P_value Lower_2.5 Upper_97.5
 Intercept 251.41 6.82 36.84    0.00    237.93     264.88
      Days  10.47 1.55  6.77    0.00      7.41      13.52


summarize_model(tmb_model)

Variance Components:
    Group    Effect Variance    SD SD_2.5 SD_97.5 Var_prop
  Subject Intercept   565.52 23.78  15.02   37.66     0.45
  Subject      Days    32.68  5.72   3.81    8.59     0.03
 Residual        NA   654.94 25.59     NA      NA     0.52

Fixed Effects:
      Term  Value   SE     Z P_value Lower_2.5 Upper_97.5
 Intercept 251.41 6.63 37.91    0.00    238.41     264.40
      Days  10.47 1.50  6.97    0.00      7.52      13.41


summarize_model(brms_model)

Variance Components:
    Group    Effect Variance    SD SD_2.5 SD_97.5 Var_prop
  Subject Intercept   792.95 28.16  15.85   46.34     0.53
  Subject      Days    42.20  6.50   4.33    9.28     0.03
 Residual        NA   669.47 25.87  22.46   29.57     0.44

Fixed Effects:
      Term  Value   SE Lower_2.5 Upper_97.5
 Intercept 250.41 7.73    234.64     265.29
      Days  10.48 1.70      6.96      13.54


summarize_model(mgcv_model)

Variance Components:
    Group    Effect Variance    SD SD_2.5 SD_97.5 Var_prop
  Subject Intercept   627.57 25.05  16.09   39.02     0.48
  Subject      Days    35.86  5.99   4.03    8.91     0.03
 Residual        NA   653.58 25.57  22.79   28.68     0.50

Fixed Effects:
      Term  Value   SE     t P_value Lower_2.5 Upper_97.5
 Intercept 251.41 6.89 36.51    0.00    237.80     265.01
      Days  10.47 1.56  6.71    0.00      7.39      13.55

Options

summarize_model(
  lmer_model,
  ci     = FALSE,
  cor_re = TRUE,
  cor_fe = TRUE,
  digits = 3
)

Variance Components:
    Group    Effect Variance     SD Var_prop
  Subject Intercept  612.100 24.741    0.470
  Subject      Days   35.072  5.922    0.027
 Residual            654.940 25.592    0.503

Fixed Effects:
      Term   Value    SE      t P_value Lower_2.5 Upper_97.5
 Intercept 251.405 6.825 36.838   0.000   238.029    264.781
      Days  10.467 1.546  6.771   0.000     7.438     13.497


summarise_model(lmer_model, ci = FALSE)

Variance Components:
    Group    Effect Variance    SD Var_prop
  Subject Intercept   612.10 24.74     0.47
  Subject      Days    35.07  5.92     0.03
 Residual             654.94 25.59     0.50

Fixed Effects:
      Term  Value   SE     t P_value Lower_2.5 Upper_97.5
 Intercept 251.41 6.82 36.84    0.00    238.03     264.78
      Days  10.47 1.55  6.77    0.00      7.44      13.50

Get Convergence

Often a model estimated by lme4 needs a bit of help for convergence. At present this function only works for lme4, but then, that’s mostly where you’re going to have the problem. This function will start with your current model and run it long enough to the point that if you still have problems, you probably really have data or model issues, and not just an estimation problem. See this for details.

# for some reason this example apparently only works interactively, not when knit
ss2 = sleepstudy
ss2$Days = ss2$Days * 10

lmer_not_converged <- lmer(
  Reaction ~ Days + (Days|Subject),
  data = ss2
)

# lmer_not_converged

lmer_converged <- converge_it(lmer_not_converged) # final result is a converged model

lmer_converged
Linear mixed model fit by REML ['lmerMod']
Formula: Reaction ~ Days + (Days | Subject)
   Data: ss2
REML criterion at convergence: 1748.233
Random effects:
 Groups   Name        Std.Dev. Corr
 Subject  (Intercept) 24.7401      
          Days         0.5922  0.07
 Residual             25.5918      
Number of obs: 180, groups:  Subject, 18
Fixed Effects:
(Intercept)         Days  
    251.405        1.047

Michael Clark

2022-10-10

Intro

Examples

Options

Get Convergence