Efficient Solvers for SLOPE in R, Python, Julia, and C++

TL;DR

This work presents a fast, memory-efficient suite of SLOPE solvers across R, Python, Julia, and C++ that tackle generalized linear models with the sorted $\ell_1$ penalty. It hinges on a hybrid proximal gradient/coordinate descent algorithm that updates cluster-structured coefficients, extended to GLMs via IRLS, and it emphasizes full path fitting, cross-validation, and robust convergence via a duality-gap criterion. The authors show substantial speedups over existing implementations through a C++ core with language wrappers, efficient data handling (dense, sparse, out-of-core), and strong screening rules, supported by extensive benchmarks on simulated and real data, plus a real-world metabolomics application. The practical impact is a scalable, cross-language tool for SLOPE that enables fast model selection, interpretable clustering of coefficients, and seamless integration into biomedical data analyses and beyond.

Abstract

We present a suite of packages in R, Python, Julia, and C++ that efficiently solve the Sorted L-One Penalized Estimation (SLOPE) problem. The packages feature a highly efficient hybrid coordinate descent algorithm that fits generalized linear models (GLMs) and supports a variety of loss functions, including Gaussian, binomial, Poisson, and multinomial logistic regression. Our implementation is designed to be fast, memory-efficient, and flexible. The packages support a variety of data structures (dense, sparse, and out-of-memory matrices) and are designed to efficiently fit the full SLOPE path as well as handle cross-validation of SLOPE models, including the relaxed SLOPE. We present examples of how to use the packages and benchmarks that demonstrate the performance of the packages on both real and simulated data and show that our packages outperform existing implementations of SLOPE in terms of speed.

Figures (8)

• Figure 1: An illustration of the SLOPE thresholding operator for $\beta = [0.5, -0.5, 0.3, 0.7]^\intercal$, and thus $c = \{0.7, 0.5, 0.3\}$, where we consider an update for the second cluster, $C_2 = \{1, 2\}$. We have omitted the particular value of $\omega$ for brevity. The example is adapted from larsson2023. Across regions where the function is constant, the operator sets the result to be either exactly 0 or merges the cluster with another cluster, setting their magnitudes to be equal.
• Figure 2: SLOPE and lasso paths on the diabetes data set. Note that the $x$ axis is reversed. Numbers indicate the indices of the predictors.
• Figure 3: SLOPE with various level of relaxation $\gamma$, with $\gamma = 0$ (fully relaxed), $\gamma = 0.5$ (semi-relaxed), and $\gamma = 1$ (standard SLOPE).
• Figure 4: Mean-squared error (MSE) from cross-validation of $q$ and $\alpha$ for SLOPE fit to the diabetes data set. The dashed line marks the optimal value of $\alpha$ in the panel corresponding to the optimal value of $q$.
• Figure 5: Performance of the solvers for the single solution benchmark for three different settings of simulated data (\ref{['tab:simulated-data']}). Please see the text for information about the data and setup of the experiment.