Regularization

Given training data $\{(x_i, y_i)\}_{i=1}^n$ and a hypothesis class parameterized by $\theta$, ordinary empirical risk minimization solves

If the hypothesis class is rich enough, this objective is minimized by parameters that fit training-set idiosyncrasies — label noise, sampling quirks — rather than the underlying signal. Regularization modifies the objective to

where $R(\theta)$ is a penalty that grows with some notion of model complexity (typically the size of the parameters), and $\lambda \geq 0$ is the regularization strength controlling how much complexity costs relative to fit. At $\lambda = 0$ we recover unconstrained fitting; as $\lambda \to \infty$ the penalty dominates and $\theta$ is forced toward the simplest model the penalty prefers (often $\theta = 0$). The two canonical penalties are the squared L2 norm and the L1 norm:

Ridge regression (Tikhonov regularization) and the lasso (least absolute shrinkage and selection operator) are the resulting estimators. Both shrink coefficients toward zero relative to the unregularized solution, but they shrink differently: ridge shrinks proportionally and smoothly, never zeroing a coefficient unless it was already zero, while lasso can drive coefficients to exactly zero, yielding sparse solutions that perform automatic feature selection. A third option, elastic net, mixes both penalties, $\lambda_1\|\theta\|_1 + \lambda_2\|\theta\|_2^2$, to get sparsity while handling groups of correlated features more stably than lasso alone. Regularization is not limited to explicit penalty terms. Early stopping — halting gradient descent before the training loss is fully minimized — implicitly limits how far parameters travel from their initialization, acting as a penalty on parameter magnitude in models trained by iterative gradient methods. Dropout in neural networks randomly zeroes activations during training, which can be shown to approximate an L2 penalty (and an ensembling effect) at inference time. Data augmentation enlarges the effective sample size, reducing variance without touching $R(\theta)$ at all. All of these share the same purpose as explicit penalties: constrain effective model complexity relative to what the data can support.

The mechanics reduce to three linked pieces: (1) formulating the penalized objective, (2) solving or approximately solving it, and (3) selecting $\lambda$ by held-out performance rather than training loss (since training loss is always minimized at $\lambda = 0$).

Ridge regression has a closed form. With design matrix $X \in \mathbb{R}^{n \times p}$, target $y$, and objective

setting the gradient to zero gives

The addition of $\lambda I$ is exactly what makes this invertible even when $X^\top X$ is singular or ill-conditioned (e.g., $p > n$, or collinear features) — the original motivation of Tikhonov regularization in ill-posed inverse problems.

Lasso has no closed form in general because $\|\theta\|_1$ is non-differentiable at zero, but it is convex, and efficient solvers exist: coordinate descent (updating one $\theta_j$ at a time via a soft-thresholding operator) and least-angle regression (LARS) are the classical choices. The soft-thresholding update for coordinate $j$, holding all others fixed, is

where $\rho_j$ is the partial residual correlation with feature $j$ and $z_j = \sum_i x_{ij}^2$. This is the algorithmic heart of why lasso zeroes coefficients: any coefficient whose correlation with the residual falls below $\lambda$ in magnitude is snapped exactly to zero. ```

# Coordinate descent for Lasso regression

# minimizes (1/2n) ||y - X theta||^2 + lambda ||theta||_1 function LASSO_COORDINATE_DESCENT(X, y, lambda, tol, max_iter): n, p = shape(X) theta = zeros(p)

# precompute column norms z = [ sum(X[:, j]^2) for j in 1..p ] for iter in 1..max_iter: theta_old = copy(theta) for j in 1..p:

# partial residual: remove feature j's current contribution r_j = y - X @ theta + X[:, j] * theta[j] rho_j = dot(X[:, j], r_j) / n

# soft-threshold update if rho_j > lambda: theta[j] = (rho_j - lambda) / (z[j] / n) else if rho_j < -lambda: theta[j] = (rho_j + lambda) / (z[j] / n) else: theta[j] = 0.0 if max(abs(theta - theta_old)) < tol: break return theta ``` ```

# Selecting lambda via k-fold cross-validation

# candidate_lambdas: list of penalty strengths to try function SELECT_LAMBDA(X, y, candidate_lambdas, k): folds = SPLIT_INTO_K_FOLDS(X, y, k) best_lambda, best_score = None, +infinity for lambda in candidate_lambdas: fold_errors = [] for i in 1..k: X_train, y_train = folds except fold i X_val, y_val = fold i theta = FIT_REGULARIZED_MODEL(X_train, y_train, lambda) pred = X_val @ theta fold_errors.append(MEAN_SQUARED_ERROR(y_val, pred)) avg_error = mean(fold_errors) if avg_error < best_score: best_score = avg_error best_lambda = lambda return best_lambda ``` The cross-validation loop is essential and non-negotiable: because $L(\theta) + \lambda R(\theta)$ evaluated on the training set is minimized by the smallest $\lambda$ that still fits well, $\lambda$ can never be chosen by looking at training loss. It must be chosen by performance on data the fitting procedure never touched. The geometric picture below is the standard way to see why L1 penalization creates sparsity while L2 does not. Each penalty defines a constrained region $R(\theta) \leq t$ (equivalent to the penalized form by Lagrangian duality), and the unconstrained least-squares solution is pulled toward the nearest point on that region's boundary that also improves the fit. The elliptical contours are level sets of the least-squares loss.

The lasso's optimum tends to land on a corner of the diamond, where one coordinate is exactly zero, because corners are the points of the boundary most likely to first touch a shrinking elliptical contour as it grows from the unconstrained optimum. The ridge's circular boundary has no corners, so the contact point generically has all coordinates nonzero — shrunk, but not eliminated.

Time complexity. Closed-form ridge regression requires computing $X^\top X$ ($O(np^2)$), adding $\lambda I$, and solving a linear system or inverting a $p \times p$ matrix ($O(p^3)$ via Cholesky/LU). Total: $O(np^2 + p^3)$, identical in order to ordinary least squares — regularization adds no asymptotic cost here, only numerical stability. Lasso via coordinate descent costs $O(np)$ per full sweep over all coordinates (each coordinate update is $O(n)$, there are $p$ coordinates), so $O(np \cdot T)$ for $T$ sweeps to convergence; $T$ is typically small (tens to low hundreds) for well-conditioned problems but grows for highly correlated features. LARS computes the entire lasso solution path in $O(np \min(n,p))$ in the typical case, comparable to one ridge solve.

Space complexity. $O(np)$ to store the design matrix, plus $O(p^2)$ for ridge's normal-equations matrix (or $O(np)$ if solved via QR/SVD without forming $X^\top X$ explicitly, which is preferred for numerical stability). Lasso via coordinate descent needs only $O(np)$ and $O(p)$ working memory, no $p \times p$ matrix.

Cross-validation cost multiplies whichever base cost by $k$ (folds) times the number of candidate $\lambda$ values, but this is embarrassingly parallel across folds and $\lambda$ values.

Correctness — why ridge is the unique minimizer. The ridge objective $\|y - X\theta\|_2^2 + \lambda\|\theta\|_2^2$ is a strictly convex quadratic in $\theta$ whenever $\lambda > 0$, because its Hessian $2(X^\top X + \lambda I)$ is positive definite (adding $\lambda I$ to a positive semi-definite matrix makes all eigenvalues strictly positive, $\lambda + \sigma_{\min}(X^\top X) > 0$). A strictly convex quadratic has a unique global minimizer found by setting the gradient to zero, which gives the closed form directly — no local minima to worry about, unlike in non-convex neural network losses.

Why regularization reduces variance — a bias-variance derivation. For ridge regression with orthonormal design ($X^\top X = I$), the ridge estimate reduces to elementwise shrinkage of the OLS estimate: $\hat\theta_{\text{ridge},j} = \hat\theta_{\text{OLS},j} / (1 + \lambda)$. Since $\text{Var}(\hat\theta_{\text{OLS},j}) = \sigma^2$ (some constant depending on noise variance), shrinking by a factor $c = 1/(1+\lambda) < 1$ scales variance by $c^2$:

while introducing bias $\mathbb{E}[\hat\theta_{\text{ridge},j}] - \theta_j^* = (c - 1)\theta_j^*$, which grows with $\lambda$. The mean squared error decomposes as $\text{Bias}^2 + \text{Variance}$, and for small $\lambda$ the variance reduction (which enters quadratically in $c$) initially outpaces the bias increase, so MSE strictly decreases from its value at $\lambda = 0$ — this is the formal justification, not just intuition, for why some regularization is (almost) always better than none when noise is present. Beyond an optimal $\lambda^*$, further increases in bias dominate and MSE rises again — hence the U-shaped validation curve central to the Bias-Variance Tradeoff.

Why lasso produces exact zeros — the subgradient argument. At a coordinate-wise optimum of the lasso objective, the subgradient stationarity condition for $\theta_j$ is

where $\partial|\theta_j| = \text{sign}(\theta_j)$ if $\theta_j \neq 0$ and $\partial|\theta_j| \in [-1, 1]$ if $\theta_j = 0$. This means $\theta_j = 0$ is optimal whenever $\left|\frac{1}{n}X_j^\top(y - X\theta)\right| \leq \lambda$ — a whole interval of residual-correlations maps to the single point zero, unlike ridge's smooth gradient $2\lambda\theta_j$, which is only zero when $\theta_j$ itself is zero. This flat region around zero in the L1 penalty's subgradient is the precise mathematical reason sparsity emerges; it has no analogue in the differentiable L2 penalty.

Honest tradeoffs. Regularization requires an extra hyperparameter ($\lambda$, or two for elastic net) that must be tuned, adding computational cost and a risk of choosing badly on a small validation set. Ridge never performs feature selection, which can hurt interpretability with many irrelevant features. Lasso can behave erratically with strongly correlated features (arbitrarily picking one of a correlated group and zeroing the rest), which elastic net mitigates but does not eliminate. All penalty-based methods assume the "right" notion of complexity is captured by parameter magnitude, which is not always true (e.g., a large-magnitude but smooth function may generalize better than a small-magnitude jagged one in some bases).

Consider a tiny regression problem with 2 standardized features (mean zero, unit variance) and $n = 5$ points, deliberately including a spurious, weakly-correlated feature to show lasso zeroing it out.

An unregularized least-squares fit on this data gives approximately $\hat\theta_{\text{OLS}} = (2.05, -0.35)$ — feature 2 picks up a non-trivial coefficient purely from finite-sample noise, since it was constructed to be nearly uncorrelated with $y$ in the population. Applying ridge regression with increasing $\lambda$ (using the closed form $\hat\theta = (X^\top X + \lambda I)^{-1}X^\top y$):

Every ridge coefficient shrinks monotonically toward zero but $\hat\theta_2$ never reaches it exactly — consistent with the smooth-gradient argument above. Applying lasso with the same $\lambda$ grid (via coordinate descent), using the soft-threshold behavior:

At $\lambda = 0.30$, the partial correlation of feature 2 with the residual falls below the threshold $\lambda$, so the soft-thresholding update $S_\lambda(\rho_2) = \text{sign}(\rho_2)\max(|\rho_2| - \lambda, 0)$ evaluates to exactly zero — lasso has performed feature selection, discarding the spurious predictor while ridge kept a small trace of it at every $\lambda$. Cross-validation on a held-out fold would typically select $\lambda \approx 0.3$–$0.5$ for lasso here, since the true generative process only used feature 1; a validation curve plotted against $\lambda$ would show validation error decreasing sharply as the noise coefficient is suppressed, then flattening or rising slightly as $\lambda$ grows large enough to also shrink the true signal $\theta_1$.

Ridge regression entered statistics as a fix for a very practical problem: when predictor variables are highly correlated, $X^\top X$ becomes nearly singular, and OLS coefficients explode with huge variance even though predictions might be reasonable. Hoerl and Kennard showed that a small ridge penalty could dramatically stabilize the estimates with minimal added bias — the same mathematical trick Tikhonov had used decades earlier for inverse problems in physics, now reframed for statistical estimation. Tibshirani's 1996 lasso paper was motivated by the desire for a method that would do variable selection and shrinkage simultaneously, unlike the then-standard practice of stepwise selection followed by unregularized fitting on the selected subset (a two-stage process known to produce unstable, overconfident models). The elastic net addressed a specific weakness of lasso — its poor behavior with groups of correlated predictors — by convexly combining both penalties. As deep learning rose to prominence, the same underlying principle reappeared in a different vocabulary: dropout, weight decay, label smoothing, and early stopping are all regularization strategies, showing that the idea of penalizing complexity to control variance long outlived its origins in linear inverse problems.

Regularized linear models (ridge, lasso, elastic net) are the default choice in genomics and other high-dimensional, low-sample-size settings — for instance, predicting disease outcomes from thousands of gene expression measurements with only hundreds of patients, where $p \gg n$ makes unregularized regression literally underdetermined. Recommender systems and large-scale linear models in advertising and finance rely on L1 penalties to keep sparse, interpretable, and computationally cheap models at web scale. In deep learning, weight decay (L2-style penalties applied to network weights) is standard in nearly every production training recipe, from image classifiers to large language models, typically combined with dropout and data augmentation as complementary regularizers. Bayesian methods generalize the same idea: a Gaussian prior on weights yields exactly the ridge penalty and a Laplace prior yields exactly the lasso penalty under MAP Estimation, which is why regularization strength is sometimes reframed as a prior belief about parameter magnitude rather than an arbitrary tuning knob. Regularization is also the backbone of automated Model Selection pipelines, where a path of models across $\lambda$ values (e.g., the LARS-computed lasso path) replaces manually choosing which features to include.

Ask the tutor: How does the choice between ridge, lasso, and elastic net change when features are highly correlated versus nearly independent?