Using mixed effects models to analyse longitudinal gene expression can highlight differences between sample groups over time. The most widely used differential gene expression tools are unable to fit linear mixed effect models, and are less optimal for analysing longitudinal data. This package provides negative binomial and Gaussian mixed effects models to fit gene expression and other biological data across repeated samples. This is particularly useful for investigating changes in RNA-Sequencing gene expression between groups of individuals over time, as described in: Rivellese, F., Surace, A. E., Goldmann, K., Sciacca, E., Cubuk, C., Giorli, G., ... Lewis, M. J., & Pitzalis, C. (2022) Nature medicine <doi:10.1038/s41591-022-01789-0>.

Fit linear and generalized linear mixed models with various
    extensions, including zero-inflation. The models are fitted using maximum
    likelihood estimation via 'TMB' (Template Model Builder). Random effects are
    assumed to be Gaussian on the scale of the linear predictor and are integrated
    out using the Laplace approximation. Gradients are calculated using automatic
    differentiation.

Extremely efficient procedures for fitting the entire lasso or elastic-net regularization path for linear regression, logistic and multinomial regression models, Poisson regression, Cox model,  multiple-response Gaussian, and the grouped multinomial regression; see <doi:10.18637/jss.v033.i01> and <doi:10.18637/jss.v039.i05>. There are two new and important additions. The family argument can be a GLM family object, which opens the door to any programmed family (<doi:10.18637/jss.v106.i01>). This comes with a modest computational cost, so when the built-in families suffice, they should be used instead.  The other novelty is the relax option, which refits each of the active sets in the path unpenalized. The algorithm uses cyclical coordinate descent in a path-wise fashion, as described in the papers cited.

Penalized methods are useful for fitting over-parameterized models. This package includes functions for restructuring an ordinal
 response dataset for fitting continuation ratio models for datasets  where the number of covariates exceeds the sample size or when
 there is collinearity among the covariates. The 'glmnet' fitting algorithm is used to fit the continuation ratio model after data restructuring.