Ridge & Lasso Regression - Supervised Learning

Ridge & Lasso Regression

Implement ridge_fit (closed form) and lasso_fit (cyclic coordinate descent), both with an unpenalised intercept, matching scikit-learn's Ridge and Lasso.

Background

Both add a penalty on the weights to least squares to fight overfitting and collinearity. Ridge (L2) shrinks all weights smoothly and has a closed form. Lasso (L1) has no closed form, but it sets some weights to exactly zero, so it doubles as feature selection. "Why does L1 give sparsity and L2 not?" is a standard interview question.

The intercept is not penalised

Penalising b would make the fit depend on where the origin of y is. Instead, centre the data, fit without an intercept, then recover it:

x_mean = X.mean(axis=0)    (d,)       y_mean = y.mean()
Xc = X - x_mean            (n, d)     yc = y - y_mean    (n,)
... fit w on (Xc, yc) ...
b = y_mean - x_mean @ w

This is what scikit-learn does with fit_intercept=True.

1. Ridge: ridge_fit(X, y, lam)

Minimises ||yc - Xc w||^2 + lam * ||w||^2 (sklearn Ridge(alpha=lam); note: no 1/n factor). Setting the gradient to zero gives

w = (Xc^T Xc + lam * I_d)^(-1) Xc^T yc

Solve the linear system with np.linalg.solve instead of forming the inverse.

2. Lasso: lasso_fit(X, y, lam, n_sweeps=100)

Minimises (sklearn Lasso(alpha=lam) scaling, with a 1/(2n) factor):

(1 / (2n)) * ||yc - Xc w||^2 + lam * ||w||_1

by cyclic coordinate descent:

1. Start from w = 0.
2. Each sweep updates j = 0, 1, ..., d-1 in order, each update using the latest values of all other weights:

r     = yc - Xc @ w                         # current residual, (n,)
rho_j = (1/n) * Xc[:, j] @ (r + Xc[:, j] * w_j)   # fit of feature j to the residual without it
z_j   = (1/n) * Xc[:, j] @ Xc[:, j]
w_j   = soft_threshold(rho_j, lam) / z_j    # leave w_j = 0 if z_j == 0 (constant feature)

soft_threshold(a, lam) = sign(a) * max(|a| - lam, 0)

3. Run exactly n_sweeps full sweeps (no early stopping).

Inputs

• X: float array (n, d); y: float array (n,)
• lam: penalty strength >= 0; n_sweeps: int

Expected Return Format

Both functions return a tuple (w, b):

• w: NumPy weight array of shape (d,)
• b: intercept as a Python float

Edge Cases and Conventions

• With lam = 0, ridge is ordinary least squares.
• Whenever |rho_j| <= lam, lasso sets w_j to exactly 0. If lam >= max_j |Xc[:, j] @ yc| / n, every weight is zero and b = mean(y).
• The two lams are on different scales (ridge has no 1/n), exactly like sklearn's alpha.

Constraints

• NumPy only; no sklearn. Keep the residual up to date incrementally inside the lasso loop rather than recomputing Xc @ w.

Follow-up

• Sketch the L1 and L2 constraint regions and explain geometrically why L1 hits corners (zeros).
• Elastic net mixes both penalties: how does the coordinate update change?

Input: ridge_fit(np.array([[1.0], [2.0], [3.0], [4.0]]), np.array([3.0, 5.0, 7.0, 9.0]), 0.0)
Output: (array([2.]), 1.0)

Input: Same data, lam = 5
Output: (array([1.]), 3.5)

Input: X = np.random.randn(40, 5), y = 3*X[:, 0] - 2*X[:, 1] + 0.1*noise + 4, lasso_fit(X, y, 0.3)
Output: w ≈ [2.7375, -1.6758, 0.0, 0.0, 0.0], b ≈ 3.9411

import numpy as np

def ridge_fit(X, y, lam):
    """
    Ridge regression in closed form (intercept not penalised).
    Objective: ||yc - Xc w||^2 + lam * ||w||^2 on centred data.

    Args:
        X: Features, shape (n, d)
        y: Targets, shape (n,)
        lam: L2 penalty strength

    Returns:
        (w, b): weights of shape (d,) and float intercept
    """
    # Your code here
    pass


def lasso_fit(X, y, lam, n_sweeps=100):
    """
    Lasso regression by cyclic coordinate descent with soft-thresholding.
    Objective: (1/(2n)) * ||y - Xw - b||^2 + lam * ||w||_1

    Args:
        X: Features, shape (n, d)
        y: Targets, shape (n,)
        lam: L1 penalty strength
        n_sweeps: Number of full passes over the coordinates

    Returns:
        (w, b): weights of shape (d,) and float intercept
    """
    # Your code here
    pass