MMRM in R

Fitting the MMRM in R

Mixed models for repeated measures (MMRM) are a popular choice for analyzing longitudinal continuous outcomes in randomized clinical trials and beyond; see Cnaan, Laird and Slasor (1997) for a tutorial and Mallinckrodt, Lane and Schnell (2008) for a review.

This vignette shows examples from the mmrm package.

The mmrm package implements MMRM based on the marginal linear model without random effects using Template Model Builder (TMB) which enables fast and robust model fitting. Users can specify a variety of covariance matrices, weight observations, fit models with restricted or standard maximum likelihood inference, perform hypothesis testing with Satterthwaite or Kenward-Roger adjustment, and extract least square means estimates by using emmeans.

Main Features

  • Flexible covariance specification:
    • Structures: unstructured, Toeplitz, AR1, compound symmetry, ante-dependence, and spatial exponential.
    • Groups: shared covariance structure for all subjects or group-specific covariance estimates.
    • Variances: homogeneous or heterogeneous across time points.
  • Inference:
    • Supports REML and ML.
    • Supports weights.
  • Hypothesis testing:
  • C++ backend:
    • Fast implementation using C++ and automatic differentiation to obtain precise gradient information for model fitting.
    • Model fitting algorithm details used in mmrm.
  • Package ecosystems integration:
    • Integration with tidymodels package ecosystem
      • Dedicated parsnip engine for linear regression
      • Integration with recipes
    • Integration with tern package ecosystems
      • The tern.mmrm package can be used to run the mmrm fit and generate tabulation and plots of least square means per visit and treatment arm, tabulation of model diagnostics, diagnostic graphs, and other standard model outputs.

Using the mmrm::mmrm function

See also the introductory vignette

The code below implements an MMRM fit in R with the mmrm::mmrm function.

library(mmrm)
fit <- mmrm(
  formula = FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID),
  data = fev_data
)

The code specifies an MMRM with the given covariates and an unstructured covariance matrix for the timepoints (also called visits in the clinical trial context, here given by AVISIT) within the subjects (here USUBJID). While by default this uses restricted maximum likelihood (REML), it is also possible to use ML, see ?mmrm.

Printing the object will show you output which should be familiar to anyone who has used any popular modeling functions such as stats::lm(), stats::glm(), glmmTMB::glmmTMB(), and lme4::nlmer(). From this print out we see the function call, the data used, the covariance structure with number of variance parameters, as well as the likelihood method, and model deviance achieved. Additionally the user is provided a printout of the estimated coefficients and the model convergence information:

fit
mmrm fit

Formula:     FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID)
Data:        fev_data (used 537 observations from 197 subjects with maximum 4 
timepoints)
Covariance:  unstructured (10 variance parameters)
Inference:   REML
Deviance:    3386.45

Coefficients: 
                  (Intercept) RACEBlack or African American 
                   30.7774065                     1.5305950 
                    RACEWhite                     SEXFemale 
                    5.6435679                     0.3260274 
                     ARMCDTRT                    AVISITVIS2 
                    3.7744139                     4.8396039 
                   AVISITVIS3                    AVISITVIS4 
                   10.3421671                    15.0537863 
          ARMCDTRT:AVISITVIS2           ARMCDTRT:AVISITVIS3 
                   -0.0420899                    -0.6938068 
          ARMCDTRT:AVISITVIS4 
                    0.6241229 

Model Inference Optimization:
Converged with code 0 and message: 

The summary() method then provides the coefficients table with Satterthwaite degrees of freedom as well as the covariance matrix estimate:

fit |>
  summary()
mmrm fit

Formula:     FEV1 ~ RACE + SEX + ARMCD * AVISIT + us(AVISIT | USUBJID)
Data:        fev_data (used 537 observations from 197 subjects with maximum 4 
timepoints)
Covariance:  unstructured (10 variance parameters)
Method:      Satterthwaite
Vcov Method: Asymptotic
Inference:   REML

Model selection criteria:
     AIC      BIC   logLik deviance 
  3406.4   3439.3  -1693.2   3386.4 

Coefficients: 
                               Estimate Std. Error        df t value Pr(>|t|)
(Intercept)                    30.77741    0.88657 218.79000  34.715  < 2e-16
RACEBlack or African American   1.53059    0.62446 168.67000   2.451 0.015263
RACEWhite                       5.64357    0.66559 157.14000   8.479 1.56e-14
SEXFemale                       0.32603    0.53194 166.14000   0.613 0.540776
ARMCDTRT                        3.77441    1.07416 145.55000   3.514 0.000589
AVISITVIS2                      4.83960    0.80173 143.87000   6.036 1.27e-08
AVISITVIS3                     10.34217    0.82269 155.56000  12.571  < 2e-16
AVISITVIS4                     15.05379    1.31288 138.46000  11.466  < 2e-16
ARMCDTRT:AVISITVIS2            -0.04209    1.12933 138.56000  -0.037 0.970324
ARMCDTRT:AVISITVIS3            -0.69381    1.18764 158.17000  -0.584 0.559924
ARMCDTRT:AVISITVIS4             0.62412    1.85096 129.71000   0.337 0.736520
                                 
(Intercept)                   ***
RACEBlack or African American *  
RACEWhite                     ***
SEXFemale                        
ARMCDTRT                      ***
AVISITVIS2                    ***
AVISITVIS3                    ***
AVISITVIS4                    ***
ARMCDTRT:AVISITVIS2              
ARMCDTRT:AVISITVIS3              
ARMCDTRT:AVISITVIS4              
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Covariance estimate:
        VIS1    VIS2    VIS3    VIS4
VIS1 40.5544 14.3960  4.9760 13.3779
VIS2 14.3960 26.5714  2.7836  7.4773
VIS3  4.9760  2.7836 14.8980  0.9036
VIS4 13.3779  7.4773  0.9036 95.5565

Extracting effect estimates using emmeans

In order to extract relevant marginal means (LSmeans) and contrasts we can use the emmeans package. This package includes methods that allow mmrm objects to be used with the emmeans package. emmeans computes estimated marginal means (also called least-square means) for the coefficients of the MMRM.

if (require(emmeans)) {
  emmeans(fit, ~ ARMCD | AVISIT)
}
Loading required package: emmeans
mmrm() registered as emmeans extension
AVISIT = VIS1:
 ARMCD emmean    SE  df lower.CL upper.CL
 PBO     33.3 0.755 148     31.8     34.8
 TRT     37.1 0.763 143     35.6     38.6

AVISIT = VIS2:
 ARMCD emmean    SE  df lower.CL upper.CL
 PBO     38.2 0.612 147     37.0     39.4
 TRT     41.9 0.602 143     40.7     43.1

AVISIT = VIS3:
 ARMCD emmean    SE  df lower.CL upper.CL
 PBO     43.7 0.462 130     42.8     44.6
 TRT     46.8 0.509 130     45.7     47.8

AVISIT = VIS4:
 ARMCD emmean    SE  df lower.CL upper.CL
 PBO     48.4 1.189 134     46.0     50.7
 TRT     52.8 1.188 133     50.4     55.1

Results are averaged over the levels of: RACE, SEX 
Confidence level used: 0.95 

Note that the degrees of freedom choice is inherited here from the initial mmrm fit.