PaperlayerView exact revision on arXiv

Source-linked AI summary

Regularized Ordinal Regression and the ordinalNet R Package

Michael J. Wurm, Paul J. Rathouz, Bret M. Hanlon

arXiv:1706.05003v1stat.CO

TL;DR

Regularized ordinal regression is not supported by many popular penalized-regression packages, while existing ordinal methods are limited in model coverage. This paper proposes an elastic-net coordinate descent algorithm and the ELMO class, with parallel, nonparallel, and semi-parallel forms. Simulations show that different forms can each achieve better out-of-sample prediction in different scenarios, and an expression-data example demonstrates prediction and variable selection.

  • Problem

    Many popular penalized-regression packages do not fit ordinal models, and available ordinal methods are limited in regularization or model coverage.

    Source: §1 · pp. 1–2

  • Method

    The paper proposes an elastic-net coordinate descent algorithm for the ELMO class, which includes parallel and nonparallel forms and a semi-parallel form that can shrink toward parallel models.

    Source: §1 · pp. 2–3, §8 · p. 32

  • Results

    Different simulation scenarios favored parallel, nonparallel, and semi-parallel models in turn; the GSE18081 example demonstrated penalized ELMO prediction and variable selection.

    Source: §8 · p. 33

  • Takeaways & Limitations

    The ELMO framework supports regularized prediction and variable selection across ordinal and unordered categorical response forms.

    Source: §8 · p. 32

    For cumulative probability models, the coordinate descent algorithm may seek solutions outside the constrained parameter space.

    Source: §3.11 · p. 18

Abstract

from arXiv · show

Regularization techniques such as the lasso (Tibshirani 1996) and elastic net (Zou and Hastie 2005) can be used to improve regression model coefficient estimation and prediction accuracy, as well as to perform variable selection. Ordinal regression models are widely used in applications where the use of regularization could be beneficial; however, these models are not included in many popular software packages for regularized regression. We propose a coordinate descent algorithm to fit a broad class of ordinal regression models with an elastic net penalty. Furthermore, we demonstrate that each model in this class generalizes to a more flexible form, for instance to accommodate unordered categorical data. We introduce an elastic net penalty class that applies to both model forms. Additionally, this penalty can be used to shrink a non-ordinal model toward its ordinal counterpart. Finally, we introduce the R package ordinalNet, which implements the algorithm for this model class.

1 Introduction

Regularized ordinal regression options are limited: existing methods either lack regularization or support only specific models. The paper proposes an elastic-net coordinate descent algorithm and the ELMO model class, including flexible nonparallel forms and a penalty that can shrink them toward parallel forms.

1 Introduction

  • The proposed coordinate descent algorithm fits a rich class of ordinal and non-ordinal multinomial regression models with an elastic net penalty.The paper presents this as a general algorithm for the model class.Source: §1 · p. 2
  • Existing ordinal regression algorithms either do not allow regularization or apply only to specific models, limiting options when predictors are numerous.Some existing approaches also have constraints: GMIFS fits a single path without an elastic-net mixing parameter and may be computationally expensive.Source: §1 · p. 2
  • The ELMO class contains common ordinal models and gives each a parallel form and a more flexible nonparallel form.The nonparallel form can accommodate unordered categorical responses.Source: §1 · p. 2
  • The elastic net penalty applies to both forms and can shrink a nonparallel model toward its parallel counterpart.The paper describes this as useful when relaxing the parallelism assumption.Source: §1 · p. 3
  • The paper introduces the R package ordinalNet to implement the proposed algorithm.The paper outlines a package demonstration and details about its availability on CRAN.Source: §1 · p. 3

2 Elementwise link multinomial-ordinal (ELMO) class

The ELMO framework represents multinomial probabilities through a model-family function and an elementwise link, with parallel, nonparallel, and semi-parallel forms. The section also describes elastic-net penalization and the inverse-link Jacobian used by the fitting algorithm.

2.1 Notation

  • The notation section defines vector and matrix conventions, including D for the Jacobian operator.It also uses ∇ for the gradient and suppresses K in identity and all-ones vectors when clear from context.Source: §2.1 · pp. 3–4

2.2 An Introduction to the ELMO model class

  • ELMO models link class probabilities to linear predictors by composing a multinomial-ordinal family function with an elementwise standard link.Examples include cumulative probability, stopping ratio, continuation ratio, and adjacent category families, paired with links such as logit or probit.Source: §2.2 · p. 4

2.3 Family function / 2.4 Elementwise link function

  • The multinomial-ordinal family function must be invertible and satisfy a monotonicity property so parallel forms are ordinal models.The elementwise link must also be monotone and invertible; symmetric links such as logit and probit yield certain equivalences between forward and backward models.Source: §2.3 · p. 5, §2.4 · p. 6

2.5 Parallel and nonparallel forms

  • Parallel forms share coefficient columns, while nonparallel forms allow coefficients to vary and can model unordered categorical data.For nonparallel cumulative probability models, monotonic linear predictors are needed to ensure non-decreasing cumulative probabilities; constraints based on training data may not guarantee this out of sample.Source: §2.5 · pp. 6–7

2.6 Semi-parallel form

  • The semi-parallel model over-parameterizes the nonparallel form to include parallel and nonparallel coefficients, allowing penalization to retain either or both for a covariate.With an elastic-net penalty, its penalized likelihood has a unique solution in most cases.Source: §2.6 · p. 7

2.7 Elastic net penalty

  • The elastic net combines lasso and ridge penalties, with α weighting them and λ controlling the degree of penalization.Cross-validation is used to choose λ.Source: §2.7 · p. 8
  • The semi-parallel penalty can include parallel and nonparallel models as special cases by tuning the penalty on parallel terms.A large ρ yields the nonparallel model; a large λ with λρ = λ* yields a regularized parallel model.Source: Semi-parallel model · p. 9

2.8 Inverse link function Jacobian

  • The coordinate descent algorithm requires the inverse-link Jacobian, which is computed using the composite structure of ELMO links.For the model-family component, the inverse and its derivative are calculated for specific families in Appendix B.Source: §2.8 · pp. 9–10

3 Coordinate descent optimization algorithm

The algorithm fits ELMO multinomial-ordinal models with an elastic net penalty by coordinate descent, using an IRLS-style quadratic approximation. Its specification covers parametrization, penalty choices, solution paths, implementation safeguards, and a boundary issue for cumulative probability models.

3.1 ELMO parametrization with a single coefficient vector

  • A single coefficient vector β represents all three ELMO forms through η_i = X_i^T β, enabling the coordinate-wise optimization formulation.The alternative parametrization is introduced because it is more convenient for the IRLS quadratic approximation.Source: §3.1 · p. 11

3.2 Elastic net penalty

  • The elastic net objective allows coefficient-specific nonnegative penalty weights; intercepts are unpenalized, while semi-parallel models can weight parallel and nonparallel coefficients differently.The weights can also be adjusted to leave selected covariates unpenalized or penalize them to varying degrees.Source: §3.2 · p. 12

3.3 Score and information matrix

  • The algorithm derives score and Fisher information for multinomial regression from the observation likelihood and combines them across independent observations.The information expression uses the inverse link’s Jacobian and the covariance structure of the multinomial response.Source: §3.3 · pp. 12–13

3.4 Optimization outer loop (quadratic approximation)

  • The outer loop updates a second-order Taylor approximation to the log-likelihood, replacing the Hessian with expected information; the inner loop optimizes the penalized approximation by coordinate descent.The approximation is centered at the current coefficients, and the inner loop cycles through coefficient estimates until convergence.Source: §3.4 · p. 14, §3.5 · p. 14

3.5 Optimization inner loop (coordinate descent)

  • Because the elastic net penalty is nondifferentiable at zero, the penalized quadratic step has no closed-form solution and is instead optimized by cycling coefficient updates.For the cumulative probability family, the algorithm may seek an optimum outside the constrained parameter space and terminates if the solution path leaves it.Source: §3.5 · p. 14, §3.11 · p. 18

3.6 Computational efficiency and numerical stability

  • Computational updates use score and information terms rather than explicitly forming the quadratic approximation’s working response, while ELMO design matrices have exploitable sparse structure.Near-zero fitted probabilities can destabilize the information matrix; probabilities are capped only for that calculation, not for the score or likelihood.Source: §3.6 · pp. 15–16

3.7 Regularization parameter sequence

  • For α > 0, λ_max is the smallest coefficient threshold at which a penalized coefficient enters the solution path; a decreasing log-scale sequence can then be chosen down to λ_min = 0.01 × λ_max.The procedure first fits an intercept-only model by unpenalized maximum likelihood, then calculates the quadratic approximation at that solution.Source: §3.7 · p. 16

3.8 Starting values; 3.9 Stopping rule; 3.10 Algorithm summary; 3.11 Issues with the cumulative probability family

  • Warm starts reuse each solution as the next starting value, and an active set focuses updates on nonzero coefficients; starting values can use observed category frequencies for intercepts and zeros elsewhere.The section also suggests convergence thresholds based on relative objective change and summarizes the algorithm as fitting a decreasing λ sequence; cumulative probability constraints can still cause optimization to leave the feasible parameter space.Source: §3.8 · pp. 16–17, §3.9 · p. 17, §3.10 · p. 17

4 Simulation

The simulations compare out-of-sample prediction across parallel, nonparallel, and semi-parallel ELMO models under three data-generating scenarios. Semi-parallel fitting was nearly as good as or better than both alternatives in every scenario.

4 Simulation

  • When data were generated by a nonparallel model with a relatively large sample, nonparallel and semi-parallel fits were better than the parallel fit.Simulation 1 used a forward stopping ratio data-generating model.Source: Table 2 · p. 19, §4 · p. 20
  • When data were generated by a parallel mechanism with only one-third of covariates having nonzero effects, parallel and semi-parallel fits were better than the nonparallel fit.Source: §4 · p. 20
  • With a parallel data-generating mechanism but a highly nonparallel effect for the first covariate, the semi-parallel fit was better than both alternatives.The simulation results summarize 100 replicates using average log-likelihood on 10,000 generated observations per replicate.Source: §4 · p. 20, Table 3 · p. 20
  • Across all three scenarios, the semi-parallel model fit was nearly as good as or better than both the parallel and nonparallel fits.The authors characterize this as evidence of versatility, not a reason to always prefer the semi-parallel model.Source: §4 · p. 20

5 Method comparison

The method comparison evaluates prediction on a small liver-disease dataset with many gene predictors, alongside variable selection. GMIFS performed best on both reported performance criteria, while standard lasso logistic regression was comparable.

5 Method comparison

  • The analysis used expression levels for 46 genes from 56 subjects across ordered healthy, cirrhosis, and hepatocellular carcinoma categories.The dataset is GSE18081; the comparison used five-fold cross-validation.Source: §5 · pp. 20–21
  • The comparison included parallel, nonparallel, and semi-parallel cumulative-logit lasso models, plus standard and group-lasso multinomial models.Source: §5 · p. 21
  • GMIFS performed best on both out-of-sample performance criteria, while standard lasso logistic regression performed comparably; the other methods did not perform well.The comparison used out-of-sample log-likelihood and misclassification rate on the liver-disease dataset.Source: §5 · pp. 21–22, Table 4 · p. 22
  • For nonparallel cumulative-logit lasso, out-of-sample non-monotone cumulative probabilities can lead cross-validation to select the null or a highly penalized model, neither good for prediction.Source: §5 · p. 22
  • The AIC stepwise method performed poorly on out-of-sample log-likelihood, likely from overfitting and occasionally assigning very small probabilities to observed classes.Its misclassification rate was more reasonable, though other methods performed better.Source: §5 · p. 22
  • Table 5 compares coefficient estimates or nonzero-coefficient counts across methods to show gene selection; group lasso selects all three coefficients together or none.Source: §5 · p. 22, Table 5 · p. 23

6 The ordinalNet R package

ordinalNet provides fitting and prediction tools for penalized parallel, nonparallel, and semi-parallel models, with cross-validation utilities. Its arguments configure model families, links, penalties, fitting controls, and data handling.

6 The ordinalNet R package

  • The ordinalNet package provides ordinalNet for fitting parallel, nonparallel, and semi-parallel elastic-net models, with summary, coefficient, and prediction methods for fitted objects.Source: §6 · p. 24
  • ordinalNetTune reports cross-validated out-of-sample log-likelihood and misclassification rates over lambda values; ordinalNetCV also tunes lambda within each fold.Source: §6 · p. 24

Arguments

  • Arguments configure response and covariate inputs, model family and link, parallel and nonparallel terms, penalty settings, convergence controls, and whether training data are retained.Supported families include cumulative, stopping ratio, continuation ratio, and adjacent category; links include logit, probit, complementary log-log, and cauchit.Source: §6 · p. 24, Arguments · pp. 25, 27

7 Demonstration in R

The demonstrations use ordinalNet on GSE18081 to fit penalized ELMO models, tune them by cross-validation, and compare model forms. Results vary by family: parallel and semi-parallel models perform similarly for cumulative probabilities, while nonparallel models are more competitive for reverse stopping ratios.

Example 1

  • The GSE18081 demonstration fits a parallel cumulative probability model with a logit link and uses the lasso over twenty λ values.The coefficient and summary methods display estimates and fit statistics, with the best-AIC fit selected by default.Source: §7 · p. 27, Example 1 · p. 27

Example 2

  • In the semi-parallel fit, coefficient columns are nearly identical for the first several predictors but differ for the first predictor.The cumulative probability model is susceptible to non-monotone out-of-sample probabilities; this issue did not occur for the semi-parallel fit on this dataset.Source: Example 2 · p. 28

Example 3

  • The nonparallel cumulative probability solution path terminates after the third λ value when the optimum leaves the constrained parameter space.The coordinate descent procedure is not designed to handle this issue, and better solutions may exist for later λ values.Source: Example 3 · pp. 28–29

Example 4

  • ordinalNetTune uses K-fold cross-validation to evaluate out-of-sample performance across a λ sequence and can select λ by average out-of-sample likelihood.The demonstration uses five folds and a sequence of twenty λ values.Source: Example 4 · p. 29

Example 5

  • ordinalNetCV estimates out-of-sample performance while tuning λ within each cross-validation fold.The comparison uses shared fold partitions for parallel, semi-parallel, and nonparallel cumulative probability models.Source: Example 5 · p. 30
  • Parallel and semi-parallel cumulative probability models have similar performance, while the nonparallel model performs much worse.Its out-of-sample log-likelihood becomes undefined with non-monotone cumulative probabilities, and cross-validation often selects the null model.Source: Example 5 · p. 30

Example 6

  • For reverse stopping ratios, the nonparallel model is much more competitive, though parallel and semi-parallel models retain a slight edge.The authors suggest that larger samples might make the more flexible models more beneficial.Source: Example 6 · p. 31

8 Discussion

The discussion presents ELMO models and the coordinate descent algorithm as tools for penalized ordinal and multinomial regression, including prediction and variable selection. It also notes that the package omits standard errors and identifies code speed and statistical inference as future research directions.

Contributions

  • ELMO models include common categorical regression models, with parallel forms suited to ordinal data and more flexible nonparallel forms for unordered responses.The semi-parallel form uses the elastic net penalty to shrink the nonparallel model toward the parallel model.Source: §8 · p. 32
  • The proposed coordinate descent algorithm fits ELMO models with the elastic net and can also be applied to multinomial models outside ELMO.The work was motivated by extending variable-selection tools for ordinal logistic regression.Source: §8 · p. 32

Experiments

  • In simulations, parallel, nonparallel, and semi-parallel models each achieve better out-of-sample prediction in different scenarios.The GSE18081 demonstration illustrates prediction and variable selection with penalized ELMO models.Source: §8 · p. 33

Software and future work

  • The ordinalNet package implements coordinate descent for parallel, nonparallel, and semi-parallel ELMO models.The paper identifies code speedup via C++ and statistical inference as future research directions.Source: §8 · p. 33

Limitations and inference

  • ordinalNet does not provide standard errors because the authors consider them uninformative for strongly biased penalized estimates.They explain that penalization reduces estimator variance by introducing bias, which can be a major component of mean squared error.Source: §8 · p. 33

A Uniqueness of the semi-parallel model estimator

The appendix examines identifiability of the semi-parallel model’s parametrization under elastic net penalization. Although the unpenalized model has infinitely many equivalent parametrizations, the penalty generally selects a unique optimum, with a lasso exception.

Problem formulation

  • The unpenalized semi-parallel model has infinitely many parametrizations for a given model, but the elastic net penalty generally selects a unique optimum.Equivalent parametrizations have the same likelihood but different penalty terms.Source: A Uniqueness of the semi-parallel… · p. 37
  • For fixed coefficients β1,…,βK, the appendix minimizes the elastic net penalty over ζ and δ1,…,δK subject to δk + ζ = βk.This constrained optimization chooses among equivalent parametrizations.Source: A Uniqueness of the semi-parallel… · p. 37

Lasso penalty

  • With the lasso penalty, the optimum is unique when ρ is noninteger; integer-valued ρ can yield an interval of optima.If ρ ≥ K, ζ = 0 and all parallel coefficients are penalized to zero, making the fit equivalent to the nonparallel model with elastic net penalty.Source: Lasso · p. 38

Elastic net penalty

B The inverse link function and its Jacobian for specific MO families

The appendix describes how to obtain inverse links and Jacobians for specific multinomial-ordinal families. It covers cumulative probability, stopping ratio, continuation ratio, and adjacent category families, including how to fit backward forms.

Inverse links and Jacobians

  • Each multinomial-ordinal family is defined by its multivariate link function gMO(p).The inverse link hMO(δ) and its Jacobian DhMO(δ) are also needed for optimization.Source: B The inverse link function and its… · p. 39
  • The appendix provides a method to compute the link, inverse link, and Jacobian for four families: cumulative probability, stopping ratio, continuation ratio, and adjacent category.Some functions can be computed recursively; δ and p have a one-to-one correspondence.Source: B The inverse link function and its… · p. 39

Forward and backward forms

  • A backward form can be fitted by reversing the response-category order and fitting the corresponding forward model.The paper presents only one of the forward or backward forms for each family.Source: B The inverse link function and its… · p. 39

Forward cumulative probability family

The forward cumulative probability, stopping ratio, continuation ratio, and adjacent category families are defined through distinct conditional probabilities, with inverse links and Jacobians computed in closed form or recursively.

Forward cumulative probability family

  • The cumulative family defines δ_j as P(Y ≤ j), with h_MO and its Jacobian given in closed form.The supplied passage states the probability definition and closed-form availability, but does not provide the formulas.Source: Forward cumulative probability family · p. 40

Forward stopping ratio family

  • The stopping ratio family defines δ_j as P(Y = j | Y ≥ j), and computes its inverse link and Jacobian recursively.The recursion begins with the first inverse-link component equal to δ_1.Source: Forward stopping ratio family · p. 40

Forward continuation ratio family

  • The continuation ratio family defines δ_j as P(Y > j | Y ≥ j) and derives its inverse link and Jacobian from the stopping ratio functions at 1 − δ.The inverse link uses the stopping ratio inverse link; the Jacobian uses its negative.Source: Forward continuation ratio family · p. 40

Forward adjacent category family

  • The adjacent category family defines δ_j as P(Y = j + 1 | j ≤ Y ≤ j + 1), with inverse link and Jacobian computed recursively.Its computation uses Δ_j = δ_j/(1 − δ_j).Source: Forward adjacent category family · p. 41

C Quadratic approximation to the log-likelihood

The quadratic approximation to the log-likelihood is a second-order Taylor expansion evaluated at the current estimate, with the Hessian replaced by its expected value.

C Quadratic approximation to the log-likelihood

D Statistical inference

Penalized regression coefficient bias complicates inference: standard errors should be interpreted cautiously and not used for confidence intervals. The cited asymptotic inference work is specific to least squares, though simulations suggest possible extension to generalized linear models.

D Statistical inference

  • Penalized regression coefficient estimators have substantial bias, so their standard errors should be interpreted with care and not used to construct confidence intervals.The cited package deliberately does not provide standard error estimates for this reason.Source: D Statistical inference · pp. 41–42
  • Even when observations outnumber predictors, usual likelihood-based confidence intervals and hypothesis tests can lead to incorrect inference.Source: D Statistical inference · p. 42
  • The cited asymptotic test-statistic results are for least squares; simulations suggest the statistic could also be used for generalized linear models.The passage presents this as a possible path toward high-dimensional inference for penalized multinomial and ordinal regression.Source: D Statistical inference · p. 42

Related papers

Loading 1706.05003v1…