Comparison of Mixed Models Packages

[billdenney]

Thank you so much for all your work on mixed models throughout the R ecosystem!

I found these pages as a great summary of how mixed models are available within R. One thing that I have not found is a resource comparing the features of the various R packages for fitting mixed models. I wonder if I'm missing it somewhere that you're aware of or if it may be a good fit here?

As a specific example, I need to fit a linear mixed model with different factor levels and heteroskedastic residual error by crossing the factor levels. (To ensure that I'm being clear here, I have treatment effects for treatments A and B measured daily on days 1 to 5 and treatment effects for treatments C and D measured daily on days 1 to 7; days are treated categorically; and I need to estimate unstructured covariance with treatment by day and a random effect on the intercept.)

I know that I can do unstructured residual error with nlme::lme() and that I can have the treament*day interaction for random effects with lme4::lmer(), but I don't see an easy way to do both. Looking at various resources that you have helpfully compiled over time, I see that I can beat lme4::lmer() into submission (😄) for the unstructured residual covariance matrix (though it's not immediately clear how to use the result) that you wrote at https://bbolker.github.io/mixedmodels-misc/notes/varmats.html.

So, the question/ask is: Do you know of a table that identifies the features of mixed effects modeling packages for R?

[bbolker]

No, although I've had a student work on something like this in the past (having trouble locating it at the moment). Maybe I should set up a google sheet so this can be a living/crowdsourced document? There is a comparison table in the glmmTMB JSS paper, but it focuses on zero-inflation.

[billdenney]

It would be very helpful for me, and I assume others, too. Croudsourcing it seems like a good solution, if it doesn't already exist.

[jebyrnes]

I would contribute! And I bet the #rstats twitterverse would be all over this.

[bbolker]

OK, next question. (1) format? (github + pull requests would be fun, but a google sheet would be much more accessible and probably easier to edit - GH diffs on a CSV file might not be all that useful ...) (2) Once that's settled, I might put up a preliminary list of categories for review. There's a very out-of-date table from Bolker et al. 2009 MEE and a zero-inflation-focused one in Brooks et al 2017:

The possibilities are endless: LMM vs GLMM, REML or not, approx integration (PQL/Laplace/AGQ/Monte Carlo), deterministic or stochastic, Bayesian vs frequentist, choice of response distributions, choice of link functions, zero (and zero-one) inflation, covariance structures, crossed REs, dispersion models, multi-membership, ... (maybe we should just look at brms and list all of the things it can do ...)

[bbolker]

By the way @billdenney , the dispersion-modeling feature of glmmTMB (dispformula=...) might do what you want. It's less flexible than nlme's heteroscedasticity-handling with the weights argument/varFunc classes (which allow you to fit models for dispersion as functions of the predicted mean), but can handle any linear model with continuous or categorical predictors.

[billdenney]

I do think of brms as a pretty inclusive list. My thought is to start simple-ish with the list we can readily come up with (and, I'm being generous to include myself in "we"). And then, to encourage others to add to the list where a new category makes sense adding all they know across all of the packages. And, I would assume that reaching out to the package authors would also be productive.

But, having been party to similar efforts in the past, I think that the most important parts are to make it simple and low-effort to contribute. So, I think that pull requests will likely get complicated quickly because the contributions won't be rows but will be spread throughout.

I like the Google sheets idea (though I've not administered one with open contributions before).

[jebyrnes]

Heck, this is an amazing start! Google doc-ing this and letting people add more packages, features, and something more extensive than a check, maybe, would be a great start! I'm sure the community will take it from there in wild ways.

[billdenney]

I agree that adding more packages, features, and details would be useful. My initial thought of how to sketch it out would be:

• Take some of the packages that we know better (again, including myself in "we" is very generous to myself)
  • Link to CRAN
  • Link to the website for the package (e.g. GitHub pages)
• Make a list of the features by group
  • Estimation method: Frequentist or Bayesian, and allowing both of those, and maybe something else?
  • Linear, nonlinear (and these would be separate since some packages allow only linear, some allow both, and I guess there may be a package that only focuses on nonlinear)
  • Supported response variable distributions
  • Multivariate response variable support
    • (To me, this seems like it would be a subset of heteroskedastic residual variable support, but it may be important enough to include on its own.)
  • Residual variable distribution support
    • Extensible response variable distribution support (i.e. can the user easily write their own distribution function for the response variable)

For all of these, either indicate "no" or indicate how to use it (e.g. "See the 'foo' argument in the function 'bar()'" or "See ?foo" or "See the page https://example.com/details"), and blank would indicate that the answer is unknown.

Contributors should be able to add new packages, new feature groups (but hopefully that would be rare), and new features. The feature list would likely get pretty long as there will likely be a lot of features that are important but only supported by a single package.

And (getting way ahead of myself), if we really take this all the way, then it eventually becomes a JSS article and a "Mixed/Hierarchical Models" CRAN task view (edit: and I see that Ben has a ~6-year-old version of the CRAN task view here: https://bbolker.github.io/mixedmodels-misc/MixedModels.html).

[bbolker]

OK, one more (last?) question: should the primary format (for now) be

• a table-like format (i.e. like the examples above: Google sheet, rows represent packages, columns represent features/feature sets. Cells could contain either yes/no check marks (possibly with footnotes), or lists of values from a semi-controlled vocabulary
• a tag-list like (YAMLish or JSONish or something) format, e.g.

package: lme4
   links:
       CRAN: ...
       GitHub: ...
   model_types: LMM, GLMM
   estimation:
        framework: frequentist
        approach: deterministic
        GLMM_approx: Laplace, AGQ
   response_distributions:
         family: Poisson, binomial, NB2, Gamma
         links: identity, cloglog, logit, probit

(I'm just a little bit scared/paralyzed by the thought of developing this hierarchy/vocabulary completely, and I'm not crazy about YAML)

• something more free-form?

[billdenney]

To me, a table is the easiest way to receive contributions from others, so I pretty strongly lean that way.

And, it is readily extensible for new structure of input.

[billdenney]

Related to glmmTMB/glmmTMB#522, it would be useful if a category were "Automatically handles rank-deficient fixed effects?"

[bbolker]

https://docs.google.com/spreadsheets/d/19itelYaVW0U0gtNtRfqh76ZGt1awlamNcJwT71u_5Uk/edit?usp=sharing

[billdenney]

I made some tweaks, most notably adding a "Column Descriptions" worksheet so that what the column contains is clear to people entering data.

[jebyrnes]

This is great! Let me know when you want me to unleash the Twitterverse on it!