GLMMLasso: An Algorithm for High-Dimensional Generalized Linear Mixed Models Using L1-Penalization

• GLMMs extend generalized linear models by adding random effects to accommodate clustered or overdispersed data.Source: §1 · p. 1
• The method first screens variables to reduce dimensionality, then refits the selected model by maximum likelihood for accurate parameter estimates.The authors describe refitting as especially important because shrinkage can severely affect variance-component estimates.Source: §1 · pp. 1–2
• The paper identifies a gap in research on truly high-dimensional GLMMs, noting that computation is much easier for high-dimensional Gaussian linear mixed models.Source: §1 · p. 2
• GLMMLasso combines ℓ1-penalization, Laplace approximation, coordinate gradient descent, and an active-set strategy for problems with thousands of variables.Source: §1 · p. 2

2.1 Model formulation

• GLMMs model grouped outcomes using fixed and random effects, with conditional responses in the exponential family and Gaussian random effects.Source: §2.1 · p. 2

2.2 Likelihood function

• The GLMM likelihood generally contains integrals that cannot be evaluated analytically, so numerical approximations are required.Source: §2.2 · p. 3

2.3 The GLMMLasso estimator

• The GLMMLasso estimator adds an ℓ1 penalty to the likelihood to encourage sparse fixed-effect estimates when many true coefficients are zero.Source: §2.3 · p. 3
• The estimator uses the Laplace approximation to handle intractable likelihood integrals; the resulting objective is non-convex despite its convex penalty.The paper describes the approximation as accurate with low computational burden for GLMMs.Source: §2.3 · pp. 3–4

3.1 The exact GLMMLasso algorithm

• The exact GLMMLasso cycles through parameter components, minimizing the objective one parameter at a time with quadratic approximations and an inexact line search.The line search is chosen to ensure that the objective decreases.Source: §3 · p. 4, §3.1 · p. 4, §3.2 · p. 5

3.2 The (approximate) GLMMLasso algorithm

• The approximate GLMMLasso holds the random-effects mode fixed during optimization, markedly speeding computation with little accuracy loss and motivating a refit.The approximation avoids recomputing derivatives of the mode with respect to fixed effects, a substantial cost in high dimensions.Source: §3.2 · pp. 4–5

3.3 Convergence behaviour and details of the GLMMLasso algorithm

• The exact algorithm converges to a stationary point; with no penalization in a low-dimensional setting, it matches lme4’s glmer results.An active-set strategy can reduce computation when fixed effects are sparse, cycling through all coordinates only periodically; initialization uses a cross-validated Lasso GLM.Source: §3.3 · pp. 5–6
• The algorithm is implemented in the R package glmmixedlasso.The paper also describes descent-direction and line-search details for the exact and approximate procedures.Source: §3.3 · p. 6

4 The two-stage GLMMLassoLA estimator(s)

• The two-stage procedure uses penalized screening to select candidate predictors, then unpenalized maximum likelihood refitting to mitigate bias and approximation error.The motivation includes shrinkage bias in fixed-effect estimates and, in some GLMMs, bias in covariance-parameter estimates.Source: §4 · pp. 6–7

4.1 The GLMMLassoLA-MLE hybrid estimator

• The GLMMLassoLA-MLE hybrid selects variables whose initial fixed-effect coefficients are nonzero, then refits using the selected covariates.Its selection set is defined by the nonzero coefficients in the initial estimate.Source: §4.1 · p. 7

4.2 The thresholded GLMMLassoLA estimator

• The thresholded GLMMLassoLA estimator selects variables whose initial fixed-effect coefficients exceed a positive threshold before refitting.The threshold is an additional regularization parameter.Source: §4.2 · p. 7

4.3 Selection of the regularization parameters

• The authors suggest BIC for estimators (8) and (10), and AIC to select λ followed by BIC to select λthres for estimator (11).The criterion uses −2 log likelihood plus a penalty factor times estimated degrees of freedom; the factor is log(n) for BIC and 2 for AIC.Source: §4.3 · p. 7

5 Simulation Study

• GLMMLasso screening is attractive in high-dimensional data and when p > 20, is numerically stable when p < n, and can outperform GLM procedures that ignore grouping.The authors report that glmpath and glmnet may fail to converge when high-dimensional models are misspecified.Source: §5 · p. 8, §5.2 · p. 10

5.1 Preview for the logistic mixed model

• The logistic GLMMLassoLA covariance-parameter estimate is notably biased, while the hybrid and thresholded procedures improve random-effects parameter accuracy over the original estimator.Fixed-effect estimates from glmmlasso are less biased than glmpath estimates, and the hybrid comparison shows the same pattern; thresholded procedures perform inadequately relative to their hybrid counterparts.Source: §5.1 · p. 8, §5.2 · pp. 9–10
• The preview reports that more random effects make GLMMLasso screening preferable to a Lasso that ignores grouping, and the algorithm often needs about 15 iterations.Source: §5.1 · p. 8

5.3 Logistic mixed model with a growing number of noise covariates

• As noise covariates grow, glmer's out-of-sample negative log-likelihood rises polynomially, while hybrid GLMMLasso remains stable.The figure compares minus twice out-of-sample log-likelihood across a growing number of covariates; GLMMLassoLA estimators remain stable and are comparable to p-glmer in the low-dimensional framework.Source: Figure 1 · p. 10, §5.3 · pp. 10–11

5.2 High-dimensional logistic mixed model

• Hybrid GLMMLasso outperforms thresholded GLMMLasso in parameter-estimation accuracy with similar true positives, and GLMMLasso procedures outperform glmpath for screening.The glmpath comparison fits a model without random effects.Source: §5.2 · p. 10

5.4 Correlated Random Effects

• With correlated random effects, results are very similar to the uncorrelated case, and correlation-estimate bias seems less severe than variance-component bias.The correlated-effects illustration uses a correlation of ρ = 0.5 between the two random effects.Source: §5.4 · p. 11

Results

• GLMMLasso selects the same epilepsy model as backward BIC selection, retaining Intercept and Base; hybrid GLMMLasso addresses penalized-estimate bias.The table covers methods using 4,000 variables, including 3,994 noise covariates, alongside Model III based on six fixed-effect covariates.Source: §6 · pp. 11–12, Table 3 · p. 12
• Thresholded GLMMLasso includes additional noise variables and has the smallest BIC score among the models considered.Compared with hybrid GLMMLasso, those additional covariates reduce variability while leaving fixed-effect estimates unaltered.Source: §6 · p. 12

7 Concluding Remarks

• GLMMLasso combines Laplace approximation, coordinatewise optimization, and an active-set strategy to handle thousands of variables in high-dimensional GLMMs.The method is intended primarily to screen variables, followed by maximum-likelihood refitting for accurate parameter estimates.Source: §7 · p. 12
• Maximum-likelihood refitting is especially important because L1 shrinkage can severely bias variance-component estimates in GLMMs.The paper characterizes its contribution as addressing numerical challenges in nonlinear mixed-effects models with non-convex loss.Source: §7 · p. 12

Appendix A: PIRLS algorithm

• PIRLS finds the random-effects mode by minimizing the convex function S(u), using Newton–Raphson to find its global minimum.The section frames the mode as the maximizer of −S(u), equivalently the minimizer of S(u).Source: Appendix A · p. 15
• Using the current mode estimate as the next starting value often reduces each PIRLS mode update to fewer than three iterations.The PIRLS algorithm typically converges fast.Source: Appendix A · p. 15

Appendix B: Comparison of the exact and approximate GLMMLasso algorithm

• Often more than 50% less computing time, with fewer outer iterations, for the approximate algorithm across logistic and Poisson mixed models.The algorithms yield similar likelihood values, although the Poisson differences are larger.Source: Appendix B · p. 16
• The comparisons use 50 simulation runs per model and a sequence of 21 λ values per run.Reported measures include relative time, outer iterations, likelihood difference, fixed-effect difference, and active-set coincidence.Source: Appendix B · p. 16
• For logistic mixed models, exact and approximate fits have very similar parameter estimates and often coinciding active sets.Table 4 reports means and standard deviations for simulated logistic mixed-model comparisons.Source: Table 4 · p. 16, Appendix B · p. 16
• For Poisson mixed models, active sets and parameter estimates differ more because some approximate fits select extra variables with very small fixed-effect coefficients.The two-stage procedures can effectively address this extra-selection problem.Source: Appendix B · p. 16

Appendix C: Low-dimensional logistic mixed model

• The low-dimensional logistic comparisons include unpenalized ML and PQL fits, with backward selection dropping covariates whose p-values exceed 5%.The selected models retain only variables deemed significant by this iterative procedure.Source: Appendix C · pp. 16–17
• In logistic simulations, thresholded and iterative procedures have the smallest active sets, while two-stage approaches estimate fixed effects better than single-stage GLMMLasso.GLMMLasso covariance-parameter estimates are considerably biased; the other procedures' covariance estimates are similar.Source: Appendix C · p. 17

Appendix D: Simulation study for the Poisson mixed model

• The Poisson appendix considers two low-dimensional and three high-dimensional designs because its behavior differs from the binary case.The supplied passage introduces these simulation settings but does not report their outcomes.Source: Appendix D · p. 17

D.1 Summary for the Poisson mixed model

• Covariance-parameter shrinkage from L1 penalization is not an issue for Poisson mixed models, unlike logistic models and as in the Gaussian case.Source: §D.1 · p. 17
• Ignoring observation grouping with the Lasso fails to recover the true active set.Source: §D.1 · p. 17
• Thresholded GLMMLasso performs best for Poisson mixed models, whereas hybrid GLMMLasso is preferable for logistic mixed models.Source: §D.1 · p. 18
• In real-data applications, total iterations are often greater than 100, and convergence problems also occur in glmnet.Source: §D.1 · p. 18

D.2 Low-dimensional Setting

• In two low-dimensional random-intercept Poisson mixed-model designs, L1 used N = 20, nC = 10, n = 200, p = 10, θ² = 1, and s0 = 5; L2 differed at p = 50.Results over 100 simulation runs are shown in Table 7.Source: §D.2 · pp. 18–19
• Thresholded methods selected fewer variables than the hybrid procedures.Table 7 reports medians and MADs for the L1 and L2 Poisson mixed-model simulations.Source: §D.2 · p. 19, Table 7 · p. 19
• GLMMLasso covariance-parameter estimates showed no dramatic bias, while glmnet-based estimators had the worst fixed-effect estimation accuracy.Source: §D.2 · p. 19

D.3 High-dimensional Setting

• The simulations used three designs with p = 500 or 1000, N = 30 or 40, and s0 = 5.H1–H3 results are reported as medians and rescaled MADs over 100 runs.Source: §D.3 · p. 19
• For high-dimensional Poisson mixed models, glmnet screening failed and produced large standard errors; its high median true-positive count had a mean below 5.This behavior was more pronounced than in the logistic mixed model.Source: §D.3 · p. 19, Table 8 · p. 20
• The Poisson results indicate that GLMMLasso screening is important and standard Lasso ignoring the grouping structure cannot substitute for it.Source: §D.3 · p. 20