Linear Models and Regularisation

Quantile regression

Instead of predicting the average, quantile regression predicts the optimistic, typical and pessimistic cases — a band of outcomes instead of one number.

Read these first

On this page 5
  1. Why it exists
  2. How it works
  3. A real example you have seen
  4. Remember this
  5. What to learn next

One lesson, three depths. Pick the one that fits you today — you can switch any time.

Beginner — No maths. Plain English.

Quantile regression predicts a chosen slice of the outcomes — the typical case, the best case, the worst case — instead of the average.

Ask a seasoned auto driver how long the airport run takes. He will not quote an average. He says: "usually forty minutes; leave an hour if you have a flight."

Two numbers, answering two different questions: usually and worst case. The second one is what saves you.

Ordinary regression answers only with averages. Quantile regression answers the driver's way.

Why it exists

Averages hide the shape of trouble. Deliveries mostly arrive on time, but the late ones are very late. Hospital stays are mostly short, but some stretch to weeks. Incomes bunch low with a long rich tail.

For planning, the tail is the question. "Average delivery: 34 minutes" does not tell the customer when to start worrying, and it does not tell the kitchen how much buffer to promise. "Nine times out of ten, under 52 minutes" does.

A quantile is a slice point in the outcomes: the 0.5 quantile is the median — half the outcomes fall below it. The 0.9 quantile is the near-worst case — nine in ten fall below. Quantile regression fits one line per slice, each depending on the features. Fit 0.1, 0.5 and 0.9 and you get a band: optimistic, typical, pessimistic — all as functions of distance, hour, traffic.

How it works

minutes │                         . ← the late tail
        │                  .   ────────  q=0.9  "leave this much time"
        │            .  ────────
        │      ──────────.────          q=0.5  "the usual"
        │  ────.───.──                  q=0.1  "if all goes well"
        └────────────────────────── km
           short trips: narrow band    long trips: wide band

Note the band widens with distance — long trips are not slower on average by a fixed amount; they are more unpredictable. One average line cannot say that. Three quantile lines say it automatically.

A real example you have seen

Delivery apps promising "arriving in 30–45 min" are quoting quantiles, not averages. Weather apps say "rain likely between 2 and 5 pm". Airlines pad their schedules. Electricity boards plan for peak demand. All of them live on the pessimistic slice, not the middle.

Remember this

  • Averages hide the tail; quantiles are the tail, made predictable.
  • Fit several quantiles and you get a band — best case to worst case, per input.
  • The band can widen with features, capturing growing unpredictability.

What to learn next

Developer — Code and libraries.

Setup

bash
pip install scikit-learn

Outputs verified with scikit-learn 1.7.2 on CPU.

Delivery times, three slices

quantile_delivery.py
import numpy as np
from sklearn.linear_model import LinearRegression, QuantileRegressor

rng = np.random.default_rng(5)
km = rng.uniform(1, 15, 120)                        # delivery distance
minutes = 10 + 3 * km + rng.gamma(2, 2 + 0.6 * km)  # delays stretch out on long trips

X = km.reshape(-1, 1)
for q in (0.1, 0.5, 0.9):
    m = QuantileRegressor(quantile=q, alpha=0).fit(X, minutes)
    print(f"q={q}: minutes = {m.intercept_:5.1f} + {m.coef_[0]:.2f} * km")

mean = LinearRegression().fit(X, minutes)
print(f"mean : minutes = {mean.intercept_:5.1f} + {mean.coef_[0]:.2f} * km")
Output
q=0.1: minutes =  10.6 + 3.43 * km
q=0.5: minutes =  13.0 + 4.13 * km
q=0.9: minutes =  14.3 + 6.00 * km
mean : minutes =  12.6 + 4.49 * km

The walkthrough

Read the slopes, top to bottom: 3.43, 4.13, 6.00. Each lucky kilometre costs 3.4 minutes; each unlucky kilometre costs 6. The band between the 0.1 and 0.9 lines widens by 2.6 minutes per km — that is growing unpredictability, printed as a number. The simulation built this in (rng.gamma(2, 2 + 0.6 * km): delay spread grows with distance), and the quantile fits recovered it.

The mean line (slope 4.49) sits above the median (4.13). Delays are skewed — many small, few huge — and the mean gets dragged toward the huge ones. When someone asks "how long does it usually take", the median is the honest answer; the mean answers a subtly different question about totals.

What to ship. A promise like "your order in 34–52 minutes" is the q=0.1 and q=0.9 predictions at the customer's distance. Roughly 80% of deliveries should land inside — and checking that coverage on held-out data is how you validate a quantile model.

alpha=0 disables the L1 penalty (QuantileRegressor regularises by default, like lasso). The solver relies on scipy's linear programming; on large datasets try solver="highs" variants or use GradientBoostingRegressor(loss="quantile", alpha=0.9) for a nonlinear quantile model.

Common mistakes

Fitting quantiles you never validate. A q=0.9 model should see ~90% of held-out outcomes below its prediction. Count it ((y_test <= pred).mean()). Miscoverage means a misspecified model, and nothing in training warns you.

Crossing quantiles. Fit independently, the q=0.5 line can dip below the q=0.1 line in sparse regions — logically impossible statements like "median 20, best-case 24". Check for crossings; fixes include fitting on shared features with monotonicity constraints or re-sorting predictions.

Reading the median model as "robust mean". The median is genuinely resistant to outliers — a real bonus (see robust regression) — but it estimates a different quantity. If totals matter (revenue = mean × volume), you need the mean, outliers and all.

Expecting a probability distribution. Three quantiles are three points, not a full distribution. For every-quantile output, fit a dense grid or use distributional methods (researcher block).

Try it yourself

Compute coverage: generate 40 fresh points with the same recipe, and check what fraction fall below the q=0.9 line. Then swap the noise to symmetric — rng.normal(0, 3, 120) — and confirm the mean and median lines land on top of each other.

What to learn next

Researcher — Mathematics and papers.

Pinball loss

The $\tau$-quantile regression estimator (Koenker and Bassett, 1978, Regression quantiles) minimises:

$$ \hat\beta_\tau = \arg\min_{\beta} \sum_{i=1}^{n} \rho_\tau\big(y_i - x_i^\top \beta\big), \qquad \rho_\tau(u) = u\,(\tau - \mathbb{1}[u < 0]) $$

Where:

  • $\tau \in (0, 1)$ — the target quantile; $u$ — the residual.
  • $\rho_\tau$ — the pinball (check) loss: slope $\tau$ for positive residuals, $\tau - 1$ for negative.
  • $\mathbb{1}[\cdot]$ — indicator function.

Why this works: for a random variable $Y$, $\mathbb{E}[\rho_\tau(Y - c)]$ is minimised at the $\tau$-quantile of $Y$ — the asymmetry prices under- and over-prediction at $\tau : (1-\tau)$, and the optimum balances them. $\tau = 0.5$ gives absolute loss and the median; least squares' symmetric quadratic gives the mean. Quantile regression needs no distributional assumption — no Gaussianity, no constant variance — which is exactly why it thrives on heteroscedastic, skewed data.

Computation and inference

The objective is piecewise-linear and convex: a linear program, solved by simplex-type methods (Barrodale–Roberts), interior-point for large $n$ (Portnoy and Koenker, 1997), or smoothed approximations for GPU-era scale. scikit-learn wraps scipy.optimize.linprog (HiGHS).

Asymptotics: $\sqrt{n}(\hat\beta_\tau - \beta_\tau) \Rightarrow \mathcal{N}(0, \tau(1-\tau) D^{-1} \Omega D^{-1})$ with $D = \mathbb{E}[f_{y|x}(q_\tau(x))\, x x^\top]$ involving the conditional density at the quantile — hence standard errors typically come from bootstrap. The full theory lives in Koenker (2005), Quantile Regression.

The crossing problem and modern variants

Independent fits across $\tau$ can violate monotonicity $q_{\tau_1}(x) \leq q_{\tau_2}(x)$ for $\tau_1 < \tau_2$. Remedies: simultaneous estimation with non-crossing constraints (Bondell et al., 2010), or post-hoc rearrangement (Chernozhukov, Fernández-Val, Galichon, 2010), which provably improves fit.

Beyond linear: quantile forests (Meinshausen, 2006), gradient boosting with pinball loss (LightGBM objective="quantile"), and simultaneous quantile networks. Pinball loss also underpins probabilistic-forecast evaluation (CRPS is its integral over $\tau$) — the M5 forecasting competition's uncertainty track scored pinball directly. Conformalized quantile regression (Romano, Patterson, Candès, 2019) wraps any quantile model to deliver finite-sample coverage guarantees, and is the current default recipe for honest prediction intervals in ML pipelines.

What to learn next