SAS Proc Mixed: random element

[WouterMpublic]

Replicate SAS Proc Mixed in MMRM; how to adress the "random" element from the SAS code in MMRM

In repeated measure studies where we have a randomized complete block design, the blocking of the subjects is an important characteristic that can be taken into account in the statistical model. Can this also be done in MMRM, and if yes how should this look like?

For example the SAS code would look similar to this:

PROC MIXED data=my_data;
CLASS subject block treatment visit;
MODEL y = treatment|visit / DDFM = KENWARDROGER;
RANDOM block;
REPEATED visit / subject= subject type=ar(1);

In MMRM the code looks something like this:
mmrm (
formula = y ~ treatment*visit + ar1(visit | subject),
data = my_data,
method = "Kenward-Roger"
)

Where:
subject = unit on which we retrieve repeated measures
block = group of x subjects equal to amount of treatments in which we randomly decide which subject receives which treatment
visit = time interval for the repeated measure
treatment = x amount of treatments tested in the study

In the MMRM code the line "RANDOM block;" from SAS is missing. Can I still incorporate this in MMRM, and how should this be done?

[danielinteractive]

Thanks @WouterMpublic , I guess this would need random effects, so a structure for the covariance which is not nested within the subjects - and we don't allow this in the mmrm package. Mainly to keep things simple in terms of implementation and scope. If there is enough interest we might think more in this direction but it would be a very big change I think...

[WouterMpublic]

Hi @danielinteractive , thank you for the quick response. Clear that it is for now not possible, and quite a big change to make. I noticed in the nlme package with function lme that they do have a random option, next to the correlation. Maybe this can be of inspiration, but I might be completely wrong here as well ;). (I'm still struggling to get the nlme::lme working as I intend, but making progress!).

[danielinteractive]

In the FDA draft guidance for bioequivalence studies also the combination of repeated measures and random effects is needed:

[danielinteractive]

Relatedly, for the bioequivalence analysis specification I found the Julia package https://github.com/PharmCat/ReplicateBE.jl and a corresponding manuscript by @PharmCat. It seems very relevant if we want to go in this direction

[gitrman]

I thought I'd chime in to say that I also would very much appreciate if the mmrm package could be extended to handle random effects, as SAS and SPSS can. I frequently analyse group psychotherapy datasets where there are repeated observations per patient, and the patients are nested within groups. Ideally I'd model the within-patient covariance structure using one of the options available in mmrm (AR1, UN, etc) and fit a random intercept to deal with the nesting of patients within groups. Currently there is no good solution in R to run those models. mmrm has a great array of repeated covariance structures and good options for hypothesis testing (Satterthwaite and Kenward-Roger DF), but doesn't handle random terms. nlme::lme() can be used to fit models with both repeated and random effects but there are fewer repeated covariance structures available compared to mmrm, and hypothesis testing options are poor (no Satterthwaite or Kenward Roger DF, despite people requesting them for years). The lmerTest package can fit models with random effects and calculate KR DF, but doesn't handle repeated terms. Hopefully this helps illustrate that there is currently no solution in R that can fit flexible models with repeated and random statements and which allows for defensible hypothesis testing. If mmrm could be extended to handle repeated terms it would solve this problem.

The random functionality would also be helpful for analysing trials with many sites, where it might be desired to model site as a random term rather than fixed factor.

[danielinteractive]

Thanks @gitrman ! Much appreciated.

So far, almost all the practical examples involved random intercept terms for groups of subjects. I wonder if we could get a specialized solution for this use case relatively easily. But I would need to dedicate some time for looking into the details...

[PharmCat]

I thought I'd chime in to say that I also would very much appreciate if the mmrm package could be extended to handle random effects, as SAS and SPSS can. I frequently analyse group psychotherapy datasets where there are repeated observations per patient, and the patients are nested within groups. Ideally I'd model the within-patient covariance structure using one of the options available in mmrm (AR1, UN, etc) and fit a random intercept to deal with the nesting of patients within groups. Currently there is no good solution in R to run those models. mmrm has a great array of repeated covariance structures and good options for hypothesis testing (Satterthwaite and Kenward-Roger DF), but doesn't handle random terms. nlme::lme() can be used to fit models with both repeated and random effects but there are fewer repeated covariance structures available compared to mmrm, and hypothesis testing options are poor (no Satterthwaite or Kenward Roger DF, despite people requesting them for years). The lmerTest package can fit models with random effects and calculate KR DF, but doesn't handle repeated terms. Hopefully this helps illustrate that there is currently no solution in R that can fit flexible models with repeated and random statements and which allows for defensible hypothesis testing. If mmrm could be extended to handle repeated terms it would solve this problem.

The random functionality would also be helpful for analysing trials with many sites, where it might be desired to model site as a random term rather than fixed factor.

Hi! It is very computationally hard problem. In pure R any solution will be slow and maybe memory greed. You can try to use Metida.jl (https://github.com/PharmCat/Metida.jl) for Julia - this package is "evolution" of ReplicateBE.jl for most LMM problems (random and repeated effects with flexible covariance structure models) (or MixedModels.jl if you don't need structured repeated effects). For Julia you also can find CSV.jl, DataFrames.jl for any work with tables (GLM.jl HypothesisTest.jl and other for statistics). Most statistical tasks can be solved with Julia packages.

[danielinteractive]

Thanks @PharmCat, great to hear from you. I definitely want to try out Metida.jl, having become a Julia fan recently.
Do you have a documentation of the algorithm you implemented in this Julia package? I guess it is https://github.com/PharmCat/Metida.jl/blob/master/docs/src/details.md ?
For our R package mmrm, we internally use C++ via the TMB package, so it could still be possible to port something to this package potentially without getting too slow.

[PharmCat]

Thanks @PharmCat, great to hear from you. I definitely want to try out Metida.jl, having become a Julia fan recently. Do you have a documentation of the algorithm you implemented in this Julia package? I guess it is https://github.com/PharmCat/Metida.jl/blob/master/docs/src/details.md ? For our R package mmrm, we internally use C++ via the TMB package, so it could still be possible to port something to this package potentially without getting too slow.

Hi! I just use REML optimization with Newton (Hessian required) method, some minimal details you can find in docs (https://pharmcat.github.io/Metida.jl/dev/details/). So as Julia ecosystem have ForwardDiff.jl for forward differentiation - Metida.jl doesnt need additional methods to get gradient and hessian. I dont use the penalized residual sum of squares (PRSS) like lme4 or MixedModels.jl - instead I use covariance matrix inversion for REML computation with sweep operator. So with this approach (how SPSS work) I can fit models with multiple random and repeated effects. Moreover - one of the main purpose to provide user defined constructor for any covariance structure. For example (https://discourse.julialang.org/t/make-model-matrix-with-raw-data-non-real-complex-numbers-and-other-structures-in-statsmodels-and-possible-application-with-metida-spatial-covariance-with-shortest-path/109123) it is possible to use A* algorithm to estimate distance in the maze for spatial repeated effect covariance.

[gitrman]

Thought I'd mention that I came across a paper today that evaluated the type of model we've been discussing in this issue in the the context of cluster RCTs. The used an unstructured residual covariance matrix plus random intercept for cluster. I thought I'd mention this as it highlights that this kind of model has practical relevance and that it would be great if the mmrm package could be extended to handle it. In the linked paper the authors used SAS and nlme::lme() to fit the models, but noted that the lme implementation is suboptimal because the Kenward-Roger DF adjustment is not available.

[jprice80]

I just wanted to chime in as well. I also would like to congratulate the package authors on a fantastic package. I know it would be a very computationally difficult and somewhat challenging problem to incorporate complicated crossed random effects structures into the package. However as someone who performs analysis in both agricultural statistics and veterinary medicine at my university, having the ability to incorporate complex random effect structures in R would be extremely useful.

Right now the available packages either lack many covariance structures or proper degrees of freedom computations. It would be amazing to have a such as mmrm package to encompass complicated mixed model anova experimental designs in combination with proper degrees of freedom estimates especially in agricultural statistics where we compute things like split-split plots, strip plots, doubly repeated measures, latin squares and other similar designs regularly. See here for a brief page on some of the types of agricultural experiments we see and compute via proc mixed: https://stats.oarc.ucla.edu/sas/library/sas-librarysas-code-for-some-advanced-experimental-designs/