Proximal Gradient Updates for Penalized Least Squares

Why Separate the Loss and the Penalty?

In high-dimensional expression analysis, the number of measured genes can exceed the number of samples. Lasso regression remains a useful way to define a sparse fit in that setting, but its absolute-value penalty has a kink at zero. Proximal gradient handles that kink by taking a gradient step on the smooth squared-error loss and then applying a separate thresholding operation. The two parts of the update have distinct roles: the gradient responds to prediction error, while thresholding can set coefficients exactly to zero.

This tutorial derives that update for a fixed design matrix and shows how to calculate an iteration. It concerns the optimization step itself, not how to select the penalty strength or diagnose convergence of a complete numerical fit.

Objective and Assumptions

Let \(y\in\mathbb{R}^n\) be the response vector and \(X\in\mathbb{R}^{n\times p}\) the design matrix. The \(j\)th coefficient is \(\beta_j\), and \(\|\cdot\|_2\) denotes Euclidean norm. For a penalty strength \(\lambda\geq 0\), consider the lasso objective

The factor \(1/(2n)\) fixes the scaling used below: it determines the gradient and the numerical meaning of \(\lambda\). Other conventions are possible, but their gradients and penalty values must be scaled consistently. We treat \(X\) and \(y\) as fixed finite data and optimize over \(\beta\). No probability distribution for the errors is required for this optimization derivation.

For clarity, the displayed objective has no intercept. If a model includes an unpenalized intercept, it must be handled separately: its squared-error gradient is updated without soft-thresholding. Any centering, scaling, or other preprocessing of predictors must be decided before fitting, since it changes the coefficient scale and the effect of a common penalty.

Deriving the Lasso Update

For the squared-error part, differentiating with respect to \(\beta\) gives

Thus, from the current coefficient vector \(\beta^{(k)}\), the gradient step forms

The remaining operation is the proximal map of \(g(\beta)=\lambda\|\beta\|_1\). Because the penalty and squared distance separate across coordinates, each coordinate minimizes \(\lambda|u|+(u-z_j)^2/(2\eta)\). For \(u>0\), the derivative is \(\lambda+(u-z_j)/\eta\), so a positive minimizer satisfies \(u=z_j-\eta\lambda\), which is positive only when \(z_j>\eta\lambda\). For \(u<0\), the derivative is \(-\lambda+(u-z_j)/\eta\), so a negative minimizer satisfies \(u=z_j+\eta\lambda\), which is negative only when \(z_j<-\eta\lambda\). When \(|z_j|\leq\eta\lambda\), neither nonzero candidate lies in its required region, and the minimum is at zero. These cases give the soft-thresholding operator

Applying this coordinate by coordinate yields the full update:

The threshold is \(\eta\lambda\), not just \(\lambda\): changing the step size changes how far the gradient step moves and how strongly the proximal operation thresholds that step. When \(\lambda>0\), the lasso penalty is nondifferentiable at coefficients equal to zero, so a plain gradient step on the whole objective is not defined there. When \(\lambda=0\), the penalty vanishes, the objective is differentiable everywhere, and this update reduces to an ordinary gradient step on squared error.

Choosing a Step Size

A standard fixed step size is based on how rapidly the squared-error gradient can change. For any \(\beta,\gamma\),

Here \(\|X\|_2\) is the largest singular value of \(X\), so \(L\) is the largest eigenvalue of \(X^\mathsf{T}X/n\). It is a Lipschitz constant for the gradient: the inequality states how much the gradient can change per unit change in the coefficient vector. For \(L>0\), the usual safe fixed-step condition is \(0<\eta\leq 1/L\). A smaller step may be used, but it changes the amount of progress per iteration. If \(X\) is the zero matrix, then \(L=0\) and the squared-error gradient is identically zero; the usual \(1/L\) formula does not apply.

The step-size condition is grounded in a quadratic upper bound. Since the squared-error loss has an \(L\)-Lipschitz gradient, for any \(u\) and current point \(x\),

When \(\eta\leq 1/L\), replacing \(L/2\) in this upper bound by \(1/(2\eta)\) gives a possibly looser upper bound. The proximal update minimizes that bound plus \(g(u)\), because completing the square in its \(u\)-dependent terms gives \(g(u)+\|u-(x-\eta\nabla f(x))\|_2^2/(2\eta)\). The minimized surrogate is no larger at the update than at \(u=x\). Consequently \(F\) does not increase in this exact-arithmetic step. This justifies the step-size restriction for the update; it does not mean that every larger step must fail, only that this guarantee no longer follows from the stated bound.

More Predictors Than Samples

The update uses matrix-vector products, a spectral norm or eigenvalue for the step-size bound, and coordinatewise thresholding. It does not require inverting \(X^\mathsf{T}X\). Thus it is defined when \(p>n\), including cases where the design is rank deficient. Rank deficiency can mean multiple coefficient vectors minimize the lasso objective, but it does not prevent calculating a proximal-gradient step from a chosen starting vector.

Comparison, Diagnostics, and Limitations

These methods target the same penalized least-squares objective when implemented with consistent scaling, but their individual updates are not interchangeable. In particular, a proximal-gradient step thresholds the whole vector after a joint gradient move; a coordinate-descent update recalculates a single coordinate in relation to the others.

Several checks help identify a mis-specified step. Confirm that the residual in the gradient is \(y-X\beta^{(k)}\), that the factor \(1/n\) matches the chosen loss scaling, and that the threshold is \(\eta\lambda\). Check that the chosen step satisfies the bound based on the actual design matrix. If a predictor is rescaled, both the gradient geometry and the coefficient scale change; using an unchanged step-size calculation or interpreting coefficients as if the scale were unchanged can be misleading. These are checks on the update, not a substitute for assessing convergence of the full algorithm.

A common mistake is to say, “The lasso objective has no gradient, so gradient-based updates cannot be used.” The correction is that, for \(\lambda>0\), the whole objective is indeed nondifferentiable at zero, but its squared-error part has a well-defined gradient. Proximal gradient uses that smooth gradient and handles the penalty through its proximal operator. At \(\lambda=0\), there is no lasso penalty, so the objective is simply differentiable squared error.

Proximal gradient’s simplicity comes with limitations. A global step size is controlled by the largest design-matrix direction, which may make updates conservative in other directions. The method supplies iterates, not an automatic guarantee that a chosen finite number of iterations is sufficient. Also, with correlated or rank-deficient predictors, different coefficient vectors can have similar or identical fitted values; the mechanics of one update do not resolve that identifiability issue.