Support for mmrm models

[ddsjoberg]

Hello @larmarange ! I hope you've been well!!

Some colleagues developed the mmrm package, which build repeated measures models. Is this a model we can support as it is, or do I need to request they add/export additional information from the model?

THANK YOU!

library(mmrm)

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

class(fit)
#> [1] "mmrm"     "mmrm_fit" "mmrm_tmb"

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)
#> Method:      REML
#> Deviance:    3386.45
#> 
#> Coefficients: 
#>                   (Intercept) RACEBlack or African American 
#>                   30.77747548                    1.53049977 
#>                     RACEWhite                     SEXFemale 
#>                    5.64356535                    0.32606192 
#>                      ARMCDTRT                    AVISITVIS2 
#>                    3.77423004                    4.83958845 
#>                    AVISITVIS3                    AVISITVIS4 
#>                   10.34211288                   15.05389826 
#>           ARMCDTRT:AVISITVIS2           ARMCDTRT:AVISITVIS3 
#>                   -0.04192625                   -0.69368537 
#>           ARMCDTRT:AVISITVIS4 
#>                    0.62422703 
#> 
#> Model Inference Optimization:
#> Converged with code 0 and message: convergence: rel_reduction_of_f <= factr*epsmch

^{Created on 2023-07-26 with reprex v2.0.2}

[larmarange]

Let's have a look at it.

First try

library(mmrm)
library(broom.helpers)

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

tidy_plus_plus(mod)
#> ! `broom::tidy()` failed to tidy the model.
#> ✔ `tidy_parameters()` used instead.
#> ℹ Add `tidy_fun = broom.helpers::tidy_parameters` to quiet these messages.
#> ✖ Unable to identify the list of variables.
#> 
#> This is usually due to an error calling `stats::model.frame(x)`or `stats::model.matrix(x)`.
#> It could be the case if that type of model does not implement these methods.
#> Rarely, this error may occur if the model object was created within
#> a functional programming framework (e.g. using `lappy()`, `purrr::map()`, etc.).
#> # A tibble: 11 × 19
#>    term              variable var_label var_class var_type var_nlevels contrasts
#>    <chr>             <chr>    <chr>         <int> <chr>          <int> <chr>    
#>  1 (Intercept)       (Interc… (Interce…        NA unknown           NA <NA>     
#>  2 RACEBlack or Afr… RACEBla… RACEBlac…        NA unknown           NA <NA>     
#>  3 RACEWhite         RACEWhi… RACEWhite        NA unknown           NA <NA>     
#>  4 SEXFemale         SEXFema… SEXFemale        NA unknown           NA <NA>     
#>  5 ARMCDTRT          ARMCDTRT ARMCDTRT         NA unknown           NA <NA>     
#>  6 AVISITVIS2        AVISITV… AVISITVI…        NA unknown           NA <NA>     
#>  7 AVISITVIS3        AVISITV… AVISITVI…        NA unknown           NA <NA>     
#>  8 AVISITVIS4        AVISITV… AVISITVI…        NA unknown           NA <NA>     
#>  9 ARMCDTRT:AVISITV… ARMCDTR… ARMCDTRT…        NA unknown           NA <NA>     
#> 10 ARMCDTRT:AVISITV… ARMCDTR… ARMCDTRT…        NA unknown           NA <NA>     
#> 11 ARMCDTRT:AVISITV… ARMCDTR… ARMCDTRT…        NA unknown           NA <NA>     
#> # ℹ 12 more variables: contrasts_type <chr>, reference_row <lgl>, label <chr>,
#> #   estimate <dbl>, std.error <dbl>, conf.level <dbl>, conf.low <dbl>,
#> #   conf.high <dbl>, statistic <dbl>, df.error <dbl>, p.value <dbl>, n <dbl>

^{Created on 2023-07-27 with reprex v2.0.2}

It seems that there is no tidy() method for this model but tidy_parameters() does work.

I need to check how to identify variables.

[larmarange]

@ddsjoberg could you have a look at #229 ? It would requires testing to check if everything works as expected.

[larmarange]

If you want to provide them some feedback:

• they could develop a tidy() method to be sure that the tidier is working as expected
• a terms() method is missing
• a model.matrix() method is missing
• however, model.frame() is working as expected
• to be noted, weights(model) return a vector of the same length of the original data frame passed to the model, but not of the length of model.frame(model) (where NA values have been removed). With lm() or glm(), the vector obtained with weights() is of the same length as model.frame(mod)

[ddsjoberg]

@larmarange thank you soooo much!! I'll give them this feedback, and let you know!

If they are amenable to these updates, would this affect PR #229? Should I review it now, it wait to here back from them?

[larmarange]

If they implement some of these feedbacks, #229 will have to be updated as some of the custom code I wrote would be useless.

However, you can already test if the results are correct and we can implement some unit tests, that will still be valid.

Regards

[ddsjoberg]

awesome, i'll review and add unit tests! What timeline do you need to have this included in the next release?

[larmarange]

The next release will be the 7 of August. It is a minor release to be retained on CRAN (see #225 ).

We will prepare another release for mmrm support later in August or September. No rush

[ddsjoberg]

Great, thank you :)