rregressionlme4lmertestmodelsummary

How to extract the goodness-of-fit statistics from lmer() model for msummary from modelsummary package


I am using lmerTest::lmer() to perform linear regression with repeated measures data.

My model contains a fixed effect (factor with 5 levels) and a random effect (subject):

library(lmerTest) 
model_lm  <- lmer(likertscore ~ task.f + (1 | subject), data = df_long)  

I would like to include the total number of observations, the number of subjects, total R^2, and the R^2 of the fixed effects in the regression table which I generate with modelsummary().

enter image description here

I tried to extract these and build a gof_map as described by the author of the package but did not succeed. Below my model output from lmerTest::lmer() the performance measures obtained:

Linear mixed model fit by REML ['lmerModLmerTest']
Formula: likertscore ~ factor + (1 | subject)
   Data: df_long
REML criterion at convergence: 6674.915
Random effects:
 Groups   Name        Std.Dev.
 subject  (Intercept) 1.076   
 Residual             1.514   
Number of obs: 1715, groups:  subject, 245
Fixed Effects:
                      (Intercept)                         factor1                         factor2  
                           3.8262                             1.5988                             0.3388  
                      factor3                             factor4                         factor5  
                          -0.7224                            -0.1061                            -1.1102  

library("performance")
performance::model_performance(my_model)

# Indices of model performance

AIC     |     BIC | R2 (cond.) | R2 (marg.) |  ICC | RMSE | Sigma
-----------------------------------------------------------------
6692.91 | 6741.94 |       0.46 |       0.18 | 0.34 | 1.42 |  1.51

Solution

  • The problem is that one of your statistics is not available by default in glance or performance, which means that you will need to do a bit of legwork to customize the output.

    First, we load the libraries and estimate the model:

    library(modelsummary)
    library(lmerTest) 
    mod  <- lmer(mpg ~ hp + (1 | cyl), data = mtcars) 
    

    Then, we check what goodness-of-fit statistics are available out-of-the-box using the get_gof function from the modelsummary package:

    get_gof(mod)
    #>        aic      bic r2.conditional r2.marginal       icc     rmse    sigma nobs
    #> 1 181.8949 187.7578      0.6744743   0.1432201 0.6200592 2.957141 3.149127   32
    

    You'll notice that there is no N (subject) statistic there, so we need to add it manually. One way to do this in a replicable way is to leverage the glance_custom mechanism described in the modelsummary documentation. To do this, we need to know what the class of our model is:

    class(mod)[1]
    #> [1] "lmerModLmerTest"
    

    Then, we need to define a method for this class name. This method should be called glance_custom.CLASSNAME. In lmerModLmerTest models, the number of groups can be retrieved by getting the ngrps object in the summary. So we do this:

    glance_custom.lmerModLmerTest <- function(x, ...) {
      s <- summary(x)
      out <- data.frame(ngrps = s$ngrps)
      out
    }
    

    Finally, we use the gof_map argument to format the result how you want it:

    gm <- list(
      list(raw = "nobs", clean = "N", fmt = 0),
      list(raw = "ngrps", clean = "N (subjects)", fmt = 0),
      list(raw = "r2.conditional", clean = "R2 (conditional)", fmt = 0),
      list(raw = "r2.marginal", clean = "R2 (marginal)", fmt = 0),
      list(raw = "aic", clean = "AIC", fmt = 3)
    )
    
    modelsummary(mod, gof_map = gm)
    
    Model 1
    (Intercept) 24.708
    (3.132)
    hp -0.030
    (0.015)
    N 32
    N (subjects) 3
    R2 (conditional) 1
    R2 (marginal) 0
    AIC 181.895