Probabilistic forecasting: quantiles, CRPS and conformal intervals

Before this lesson: backtesting and point metrics (lesson 03), ARIMA prediction intervals (lesson 06), forecasting with lag features and gradient boosting (lesson 08). From the probability-and-statistics track: the cumulative distribution function (CDF), quantiles, and the normal distribution. Calibration of classifiers from the ml-fundamentals track helps, but this lesson defines everything it uses. Time: about 110 minutes, plus the notebook.

---

Why this matters #

• Decisions depend on quantiles, not on the mean. A shop that loses 9 per missing unit and 1 per unsold unit should stock the 90 % quantile of demand (section 3.2). A grid operator keeps reserve power for the 99 % quantile of load. A point forecast cannot answer these questions.
• Model intervals are often too narrow. The ARIMA and ETS intervals of lessons 04 and 06 assume that the model is right, that its parameters are exact and that the errors are Gaussian. In real backtests a "95 %" interval often covers clearly less than 95 % of the outcomes. Conformal methods correct this.
• Benchmarks and competitions score distributions. The M5 uncertainty competition scored nine quantiles with a scaled pinball loss; the GEFCom2014 energy competition used the pinball loss; weather ensembles (ECMWF, and ML models such as GenCast) are scored with the CRPS; time-series foundation models such as Chronos and TimesFM output quantiles and are compared by a weighted quantile loss (lesson 09).
• Monitoring. A value outside a calibrated 99 % interval is a natural anomaly alarm (lesson 11).

1. Predict first #

2. The picture #

Four forms of a probabilistic forecast #

Write $y_{t+h}$ for the unknown future value and $\hat y_{t+h\mid t}$ for a point forecast made at the forecast origin $t$ (the track notation). A probabilistic forecast describes the whole predictive distribution of $y_{t+h}$ given the data up to time $t$. It comes in four forms:

Form	What the model outputs	Example	Good for
distribution	a CDF $F_{t+h\mid t}$ or its parameters	$\mathcal N(\hat y, \hat\sigma_h^2)$ from ARIMA (lesson 06)	any summary of one horizon
quantiles	$\hat q_\tau$ for a few levels $\tau$	$\hat q_{0.1}, \hat q_{0.5}, \hat q_{0.9}$ from a quantile model	decisions at fixed levels
interval	$[\ell, u]$ with nominal coverage $1 - \alpha$	a 90 % interval	dashboards, alarms
sample paths	$M$ simulated futures $y^{(i)}_{t+1}, \dots, y^{(i)}_{t+H}$	DeepAR samples (lesson 09)	sums and maxima over many horizons

The first three describe one horizon at a time (the marginal distribution of each $y_{t+h}$). Only sample paths keep the dependence between horizons. You need them for questions about several horizons at once, such as "what is the 90 % quantile of the total demand next week?" (section 3.7).

The pinball loss turns a quantile into a training target #

The $\tau$-quantile of a distribution is the value $q$ with probability $\tau$ below it. To train a model to output a quantile, we need a loss whose best constant is that quantile. The squared error gives the mean, and the absolute error gives the median (lesson 03). The pinball loss (or quantile loss) tilts the absolute error: with the error $u = y - q$,

For $\tau = 0.9$, a forecast that is too low costs $0.9$ per unit and a forecast that is too high costs only $0.1$ per unit. So the best $q$ moves up, until only 10 % of the data lie above it.

Interactive figure: loading…

Calibration and sharpness #

How do we check a probabilistic forecast? We cannot check one distribution with one outcome. But over many forecasts, events with a forecast probability of 80 % should happen about 80 % of the time. This is calibration. A simple check is the coverage of central intervals: the share of outcomes inside the 50 %, 80 % and 95 % intervals. A finer check is the PIT (probability integral transform): for each forecast, compute $u_t = F_t(y_t)$, the forecast probability below the outcome. If the forecasts are calibrated, the $u_t$ are uniform on $[0, 1]$ and their histogram is flat.

Calibration is not enough. A forecast that ignores all current information and always predicts the long-run distribution of the series (the climatological forecast) is calibrated, but its intervals are wide. We also want sharpness: narrow intervals. The rule of Gneiting, Balabdaoui and Raftery (2007) is maximise sharpness subject to calibration.

Interactive figure: loading…

3. Definitions and main results #

3.1 Notation #

• $F$ is the predictive CDF of $y_{t+h}$ given the data up to $t$: $F(x) = P(y_{t+h} \le x)$ under the forecast. We drop the indices when only one forecast is meant.
• The quantile function is $F^{-1}(\tau) = \inf\{x : F(x) \ge \tau\}$ for $\tau \in (0, 1)$. A quantile forecast at level $\tau$ is written $\hat q_\tau$.
• A central $(1 - \alpha)$ interval is $[F^{-1}(\alpha/2), F^{-1}(1 - \alpha/2)]$. Its nominal coverage is $1 - \alpha$; its empirical coverage is the share of outcomes inside it in a backtest.

3.2 Quantiles and the pinball loss #

Definition. For a level $\tau \in (0, 1)$, the pinball loss of a quantile forecast $q$ when the outcome is $y$ is

For $\tau = 0.5$ it is half the absolute error. It is never negative and it is 0 only when $q = y$.

Theorem (the pinball loss is minimised by the quantile). Let $Y$ have a continuous CDF $F$ with finite mean. Then $L(q) = \mathbb E\,\rho_\tau(Y - q)$ is convex in $q$, and it is minimised at every $q$ with $F(q) = \tau$, that is at the $\tau$-quantile $F^{-1}(\tau)$.

Why this is true, in one sentence. Moving $q$ up by a small $\delta$ costs $(1 - \tau)\delta$ for each outcome below $q$ (probability $F(q)$) and saves $\tau\delta$ for each outcome above $q$ (probability $1 - F(q)$). The net change is $\delta\,[F(q)(1 - \tau) - (1 - F(q))\tau] = \delta\,[F(q) - \tau]$.

Proof.

1. Split the expectation at $q$: $L(q) = \tau\int_q^\infty (y - q)\,dF(y) + (1 - \tau)\int_{-\infty}^q (q - y)\,dF(y)$.
2. Differentiate with respect to $q$. The integrands are 0 at $y = q$, so only the derivative of the integrands counts: $L'(q) = -\tau\,(1 - F(q)) + (1 - \tau)\,F(q) = F(q) - \tau$.
3. $F$ is non-decreasing, so $L'$ is non-decreasing and $L$ is convex.
4. $L'(q) < 0$ when $F(q) < \tau$ and $L'(q) > 0$ when $F(q) > \tau$. So $L$ decreases up to the $\tau$-quantile and increases after it: the minimum is where $F(q) = \tau$. $\blacksquare$

The empirical version. For a sample $y_1, \dots, y_n$, the mean loss $L_n(q) = \frac1n\sum_i \rho_\tau(y_i - q)$ is piecewise linear with kinks at the data. Its slope is $F_n(q) - \tau$, where $F_n(q)$ is the share of data $\le q$. If $n\tau$ is not an integer, the unique minimiser is the $\lceil n\tau\rceil$-th smallest value (np.quantile(y, tau, method="inverted_cdf")). If $n\tau$ is an integer, every $q$ between the $n\tau$-th and the $(n\tau + 1)$-th smallest value is a minimiser.

Quantile regression. Replace the constant $q$ by a model $\hat q_\tau(x)$ and minimise the mean pinball loss over the training data. This is quantile regression (Koenker and Bassett, 1978). The theorem, applied for each $x$, says that the best model outputs the conditional $\tau$-quantile of $y$ given $x$. It works with any model:

• linear: sklearn.linear_model.QuantileRegressor(quantile=0.9) (a linear program);
• gradient boosting: HistGradientBoostingRegressor(loss="quantile", quantile=0.9) on the lag features of lesson 08; LightGBM has objective="quantile";
• neural networks: a multi-quantile head outputs $\hat q_{\tau_1}, \dots, \hat q_{\tau_K}$ together and is trained with $\sum_k \rho_{\tau_k}$ (MQ-RNN, the Temporal Fusion Transformer, Chronos-Bolt; lesson 09).

Quantile crossing. Models trained separately for each $\tau$ can give $\hat q_{0.9}(x) < \hat q_{0.8}(x)$ for some $x$. This is impossible for a real distribution. Simple fixes: sort the predicted quantiles at each $x$ (this never increases the total pinball loss; Chernozhukov, Fernández-Val and Galichon, 2010), or predict $\hat q_{\tau_1}$ and positive increments $\hat q_{\tau_{k+1}} - \hat q_{\tau_k} = \operatorname{softplus}(\cdot)$.

The stocking decision (newsvendor problem). Order $q$ units. Each unsold unit costs $c_o$ and each missing unit costs $c_u$. The cost is $c_u(y - q)$ if $y > q$ and $c_o(q - y)$ if $y < q$. This is $(c_u + c_o)\,\rho_\tau(y - q)$ with $\tau = c_u/(c_u + c_o)$. By the theorem, the best order is the quantile at the critical ratio $\tau = c_u/(c_u + c_o)$. With $c_u = 9$ and $c_o = 1$: $\tau = 0.9$. This answers question 1 of section 1: order the 90 % quantile of demand, not the mean or the median.

3.3 Calibration, sharpness and the PIT #

Definition. A sequence of forecasts $F_t$ for outcomes $y_t$ is probabilistically calibrated if the PIT values $u_t = F_t(y_t)$ are uniform on $[0, 1]$. Equivalently, for every level $\tau$, the share of outcomes below $F_t^{-1}(\tau)$ is $\tau$; so every central $(1 - \alpha)$ interval covers $1 - \alpha$ of the outcomes.

Theorem (the PIT of the true distribution is uniform). If $Y$ has a continuous CDF $F$, then $F(Y)$ is uniform on $[0, 1]$.

Proof. For $u \in (0, 1)$, let $x_u = F^{-1}(u)$; continuity gives $F(x_u) = u$. Then $P(F(Y) \le u) = P(Y \le x_u) = F(x_u) = u$, which is the CDF of the uniform distribution. $\blacksquare$

So a forecaster that knows the true conditional distribution produces a flat PIT histogram. How to read a histogram that is not flat:

PIT histogram	Coverage of central intervals	Diagnosis
flat	as nominal	calibrated
U-shaped	below nominal	too narrow (under-dispersed, over-confident)
hump in the middle	above nominal	too wide (over-dispersed)
sloped (mass near 0)	lower tail too often	forecasts too high (bias)
sloped (mass near 1)	upper tail too often	forecasts too low (bias)

Calibration is necessary but not sufficient. The climatological forecast (always the long-run distribution) is calibrated, because the outcomes really come from that distribution on average. But it uses no information, so it is not useful. Sharpness is a property of the forecasts alone: the width of their intervals (or the spread of $F_t$). The paradigm of Gneiting, Balabdaoui and Raftery (2007) is: maximise sharpness subject to calibration. This answers question 2 of section 1: both forecasts are calibrated, but forecast B is sharper (if it uses real information) and therefore better. A proper scoring rule such as the CRPS (next section) sees the difference.

Discrete outcomes and a flat PIT that hides problems

For counts, $F_t$ jumps, and $F_t(y_t)$ is not uniform even for the true distribution. Use the randomised PIT $u_t = F_t(y_t - 1) + V_t\,[F_t(y_t) - F_t(y_t - 1)]$ with $V_t$ uniform on $[0, 1]$, or the non-randomised version of Czado, Gneiting and Held (2009). Also note that a flat histogram over all times can hide errors that cancel: too narrow in winter and too wide in summer gives a flat total histogram (Hamill, 2001). Check the PIT in subgroups (season, horizon, level of the series) as well.

3.4 Proper scoring rules and the CRPS #

A scoring rule $S(F, y)$ gives a number (lower is better here) to a forecast $F$ when $y$ happens. It is proper if, when $y$ really comes from $G$, the expected score $\mathbb E_{Y \sim G}\, S(F, Y)$ is smallest at $F = G$; strictly proper if $F = G$ is the only minimiser. A proper score cannot be gamed: the best strategy is to report what you really believe. Coverage alone is not a score of this kind: an interval of infinite width always has perfect coverage.

Definition. The continuous ranked probability score is

It compares the forecast CDF with the "CDF" of the outcome, a step from 0 to 1 at $y$. It is in the units of $y$, and it is strictly proper for distributions with a finite mean (Gneiting and Raftery, 2007).

Interactive figure: loading…

Three equivalent forms. For $F$ with a finite mean:

(1) Integral (threshold) form: the definition above.

(2) Kernel form (Gneiting and Raftery, 2007): with $X, X'$ independent draws from $F$,

The first term rewards forecasts close to $y$; the second term rewards spread, so a forecast does not gain by collapsing to a point. For an ensemble or a set of sample paths $x_1, \dots, x_M$, replace the expectations by means over the members.

(3) Quantile form: the CRPS is twice the pinball loss, averaged over all levels:

Proof of the quantile form

Assume $F$ is continuous and strictly increasing, so $F^{-1}(\tau) > x \iff \tau > F(x)$.

1. For a fixed $q$, $\rho_\tau(y - q) = \int_y^q \big(\mathbf 1\{x \ge y\} - \tau\big)\,dx$. Check: if $q > y$ the integrand is $1 - \tau$ and the integral is $(1 - \tau)(q - y)$; if $q < y$ the integral is $-\int_q^y (0 - \tau)\,dx = \tau(y - q)$. Both agree with the definition.
2. Put $q = F^{-1}(\tau)$ and integrate over $\tau$: $2\int_0^1 \rho_\tau(y - F^{-1}(\tau))\,d\tau = 2\int_0^1 \int_y^{F^{-1}(\tau)} (\mathbf 1\{x \ge y\} - \tau)\,dx\,d\tau$.
3. Swap the order of integration (Fubini). For $x > y$, the point $x$ lies between $y$ and $F^{-1}(\tau)$ exactly when $\tau > F(x)$; it contributes $\int_{F(x)}^1 (1 - \tau)\,d\tau = \tfrac12 (1 - F(x))^2$.
4. For $x < y$, the point $x$ lies between $F^{-1}(\tau)$ and $y$ exactly when $\tau < F(x)$; the integral runs from right to left, so it contributes $-\int_0^{F(x)} (0 - \tau)\,d\tau = \tfrac12 F(x)^2$.
5. Multiply by 2 and add: $\int_{x < y} F(x)^2\,dx + \int_{x > y} (1 - F(x))^2\,dx$. This is the integral of $(F(x) - \mathbf 1\{x \ge y\})^2$: the CRPS. $\blacksquare$

Two special cases.

• Point forecast. If $F$ is a step at $\hat y$, then $(F(x) - \mathbf 1\{x \ge y\})^2 = 1$ between $\hat y$ and $y$ and 0 elsewhere, so $\operatorname{CRPS} = |y - \hat y|$. The CRPS generalises the absolute error, and you can compare a probabilistic forecast with a point forecast on the same scale.
• Gaussian forecast $F = \mathcal N(\mu, \sigma^2)$, with $z = (y - \mu)/\sigma$, $\Phi$ the standard normal CDF and $\varphi$ its density (Gneiting et al., 2005):

Even a perfect mean ($z = 0$) gives $\sigma(2\varphi(0) - 1/\sqrt\pi) \approx 0.234\,\sigma$: the price of the uncertainty.

In practice: a finite set of quantiles. Many models output only $K$ quantiles. The quantile form suggests the approximation $\operatorname{CRPS} \approx \frac{2}{K}\sum_k \rho_{\tau_k}(y - \hat q_{\tau_k})$ for evenly spaced levels. Summed over a test set and divided by $\sum_t |y_t|$, this is the weighted quantile loss (WQL) used, for example, with $\tau = 0.1, 0.2, \dots, 0.9$ to evaluate Chronos (Ansari et al., 2024) and in the GIFT-Eval benchmark (Aksu et al., 2024). Epidemic forecast hubs use the closely related weighted interval score (Bracher et al., 2021).

3.5 Split conformal prediction #

The intervals so far come from a model. If the model is wrong, the coverage is wrong. Conformal prediction builds intervals around any point forecaster from its errors on held-out data, with a coverage guarantee that does not need a correct model (Vovk, Gammerman and Shafer, 2005; Angelopoulos and Bates, 2021).

Split conformal for regression. You have a fitted point forecaster $\hat y$ and a calibration set of $n$ new pairs that were not used to fit it.

1. Compute the nonconformity scores $s_i = |y_i - \hat y_i|$, $i = 1, \dots, n$.
2. Set $k = \lceil (n + 1)(1 - \alpha)\rceil$. If $k \le n$, let $\hat q$ be the $k$-th smallest score; if $k > n$, set $\hat q = \infty$.
3. For a new case, report the interval $[\hat y - \hat q,\ \hat y + \hat q]$.

Exchangeability. Random variables $Z_1, \dots, Z_{n+1}$ are exchangeable if their joint distribution does not change when you permute them. I.i.d. data are exchangeable; a time series with autocorrelation or a trend is not.

Theorem (split conformal coverage). If the calibration scores and the test score $s_1, \dots, s_n, s_{n+1}$ are exchangeable, then

If the scores have no ties (for example, continuous data), the coverage is also below $1 - \alpha + \frac{1}{n+1}$.

Proof.

1. $y_{n+1}$ is in the interval exactly when $s_{n+1} \le \hat q$.
2. Assume no ties. Let $R$ be the rank of $s_{n+1}$ among all $n + 1$ scores (1 = smallest). By exchangeability, every rank from 1 to $n + 1$ has the same probability $\frac{1}{n+1}$.
3. $s_{n+1} \le \hat q$ (the $k$-th smallest calibration score) means that fewer than $k$ calibration scores are below $s_{n+1}$, that is $R \le k$.
4. So the coverage is $P(R \le k) = \frac{k}{n+1}$. Since $k \ge (n + 1)(1 - \alpha)$, this is $\ge 1 - \alpha$; since $k < (n + 1)(1 - \alpha) + 1$, it is $< 1 - \alpha + \frac{1}{n+1}$.
5. With ties, $s_{n+1} \le \hat q$ also holds in some cases with $R > k$, so the coverage can only grow. $\blacksquare$

Answer to question 3 of section 1. With $n = 19$ and $\alpha = 0.1$: $k = \lceil 20 \cdot 0.9\rceil = 18$. Use the 18th smallest error. The coverage is exactly $18/20 = 90\,\%$ for continuous exchangeable data. The 17th would give only $17/20 = 85\,\%$.

Interactive figure: loading…

What the guarantee says, and what it does not.

• It is marginal: an average over calibration sets and test cases. A single calibration set gives a coverage that varies around $1 - \alpha$ (the variation shrinks like $1/\sqrt n$).
• It is not conditional: a 90 % interval can cover 99 % in summer and 70 % in winter and still be 90 % on average. With a constant width, the interval is too wide where the noise is small and too narrow where it is large.
• Conformalized quantile regression (CQR; Romano, Patterson and Candès, 2019) gives intervals whose width changes with $x$. Fit quantile models $\hat q_{\alpha/2}(x)$ and $\hat q_{1-\alpha/2}(x)$; use the scores $s_i = \max\{\hat q_{\alpha/2}(x_i) - y_i,\ y_i - \hat q_{1-\alpha/2}(x_i)\}$; report $[\hat q_{\alpha/2}(x) - \hat q,\ \hat q_{1-\alpha/2}(x) + \hat q]$. The same theorem gives coverage $\ge 1 - \alpha$.
• The calibration set must not be used to fit the model. Training residuals are too small, and the intervals then under-cover.

Why time series break exchangeability. Calibration errors come from the past and the test error from the future. (1) Errors are autocorrelated: a run of large errors in the calibration window makes all of them similar, so the window is effectively smaller than $n$. (2) The series drifts: the level, the noise or the model quality change, as in the widget. (3) Errors at horizon $h$ are larger than at horizon 1, so scores must come from the same horizon. The theorem then gives no guarantee, though the method often still works reasonably when the changes are slow.

3.6 Conformal methods for time series #

Practical recipe. Run a rolling-origin backtest (lesson 03). For each horizon $h$, keep the absolute $h$-step errors of the last $n$ origins as the calibration scores, and recompute $\hat q_h$ at each new origin. The window follows slow changes; the guarantee is only approximate.

Adaptive conformal inference (ACI; Gibbs and Candès, 2021) adds feedback. It keeps a working level $\alpha_t$ and changes it after each outcome:

1. Build the interval $\hat C_t(\alpha_t)$ with the conformal recipe at level $\alpha_t$ (for example, the $\lceil (n + 1)(1 - \alpha_t)\rceil$-th smallest of the last $n$ scores). If $\alpha_t \le 0$, the interval is the whole line; if $\alpha_t \ge 1$, it is empty.
2. Observe $y_t$ and set $\mathrm{err}_t = 1$ if $y_t \notin \hat C_t(\alpha_t)$, else 0.
3. Update the level with a step size $\gamma > 0$:

After a miss, $\alpha_t$ drops by $\gamma(1 - \alpha)$ and the next interval is wider; after a hit, $\alpha_t$ rises by $\gamma\alpha$ and the next interval is a little narrower.

Theorem (long-run coverage of ACI; Gibbs and Candès, 2021). For any sequence of data (no assumption at all) and $\alpha_1 \in [0, 1]$,

Proof.

1. $\alpha_t$ stays in $[-\gamma, 1 + \gamma]$. One step changes $\alpha_t$ by at most $\gamma$. If $\alpha_t < 0$, the interval is the whole line, so $\mathrm{err}_t = 0$ and $\alpha_t$ increases. If $\alpha_t > 1$, the interval is empty, so $\mathrm{err}_t = 1$ and $\alpha_t$ decreases. So $\alpha_t$ can leave $[0, 1]$ by at most one step.
2. Telescoping. Summing the update over $t = 1, \dots, T$ gives $\alpha_{T+1} - \alpha_1 = \gamma\sum_{t=1}^T(\alpha - \mathrm{err}_t)$.
3. Divide by $\gamma T$: $\big|\alpha - \frac1T\sum_t \mathrm{err}_t\big| = \frac{|\alpha_{T+1} - \alpha_1|}{\gamma T}$.
4. By step 1, $|\alpha_{T+1} - \alpha_1| \le \max(\alpha_1 + \gamma,\ 1 + \gamma - \alpha_1)$. $\blacksquare$

Interactive figure: loading…

What ACI does not promise. The guarantee is about the long-run average of misses. It says nothing about coverage at a given time, in a given season, or conditional on anything. The interval can be infinite for a few steps (when $\alpha_t < 0$). A large $\gamma$ reacts fast but makes the width jumpy; a small $\gamma$ is smooth but slow. Later methods choose $\gamma$ online or combine several values.

Two related methods, in short.

• EnbPI (Xu and Xie, 2021) fits an ensemble of models on bootstrap samples of the training data. The out-of-bag residuals are the scores, and the newest residuals replace the oldest ones as data arrive. It has an approximate guarantee when the errors are stationary and the model is consistent.
• Conformal PID control (Angelopoulos, Candès and Tibshirani, 2023) updates the quantile $\hat q_t$ directly, $\hat q_{t+1} = \hat q_t + \eta(\mathrm{err}_t - \alpha)$ (the "P" part, a quantile tracker), and adds an integral term and an optional forecast of the next score ("scorecaster"), as in a PID controller. It has the same kind of long-run guarantee.

3.7 Many horizons and sample paths #

A forecast for each horizon is not enough for a question about several horizons. Let $e_1, e_2$ be the errors of two days, each $\mathcal N(0, 1)$ and independent. The 90 % quantile of each is $1.28$. The 90 % quantile of the two-day total error is $1.28\sqrt 2 = 1.81$, not $1.28 + 1.28 = 2.56$. The quantile of a sum is not the sum of the quantiles; they agree only when the errors are perfectly correlated. With positive correlation, the truth is between the two numbers.

So for weekly totals, maxima over a day, or "the first time the stock runs out", use sample paths: simulate many futures $y^{(i)}_{t+1}, \dots, y^{(i)}_{t+H}$, compute the quantity on each path, and take its empirical quantiles. Paths come from a model (simulate the future shocks of an ARIMA or ETS model, or resample past residuals: a bootstrap) or from a generative network such as DeepAR (lesson 09). Score them with the sample CRPS (kernel form) per horizon, or with a multivariate score such as the energy score.

4. Worked example #

(a) A quantile from the pinball loss. Data $\{2, 4, 5, 7, 12\}$, $\tau = 0.7$. The mean pinball loss at $q = 7$: the outcomes 2, 4, 5 are below $q$ and cost $(1 - 0.7)(7 - y)$; the outcome 12 is above and costs $0.7 \cdot 5$:

In the same way $L(5) = \tfrac15[0.9 + 0.3 + 0 + 1.4 + 4.9] = 1.5$ and $L(12) = \tfrac15[3 + 2.4 + 2.1 + 1.5 + 0] = 1.8$. The slope just left of 7 is $F_n = 3/5 = 0.6$ minus $0.7$, negative; just right of 7 it is $4/5 - 0.7 = 0.1$, positive. So $q = 7$ is the minimiser. It is the 4th smallest value, since $\lceil 5 \cdot 0.7\rceil = 4$.

(b) The CRPS three ways. Forecast $\mathcal N(0.3, 1.5^2)$, outcome $y = 1$. Then $z = 0.7/1.5 = 0.466667$, $\Phi(z) = 0.679631$, $\varphi(z) = 0.357783$ and $1/\sqrt\pi = 0.564190$. The three terms in the bracket are $z(2\Phi(z) - 1) = 0.466667 \cdot 0.359262 = 0.167655$, $2\varphi(z) = 0.715567$ and $-0.564190$:

The kernel form gives $\mathbb E|X - y| - \frac12\mathbb E|X - X'| = 1.324833 - 1.5/\sqrt\pi = 1.324833 - 0.846284 = 0.478549$, and the integral of $2\rho_\tau$ over $\tau$ gives $0.4785$ too (section 5 checks all three). For comparison, the absolute error of the median forecast is $|1 - 0.3| = 0.7$: the spread of the forecast lowers the score because the outcome was close to the centre.

(c) A split conformal interval. The absolute errors of 9 calibration months are $0.4, 1.1, 0.7, 2.3, 0.2, 1.6, 0.9, 3.0, 1.3$. For $\alpha = 0.2$:

1. $k = \lceil (9 + 1)(1 - 0.2)\rceil = \lceil 8 \rceil = 8$.
2. Sorted: $0.2, 0.4, 0.7, 0.9, 1.1, 1.3, 1.6, 2.3, 3.0$. The 8th smallest is $\hat q = 2.3$.
3. If the next point forecast is $\hat y = 50$, the 80 % interval is $[47.7,\ 52.3]$.
4. For exchangeable data the coverage is exactly $8/10 = 80\,\%$.

(d) Two steps of ACI. $\alpha = 0.1$, $\gamma = 0.05$, $\alpha_t = 0.1$.

• The outcome falls outside: $\mathrm{err}_t = 1$, so $\alpha_{t+1} = 0.1 + 0.05(0.1 - 1) = 0.055$. The next interval uses the $\lceil (n + 1) \cdot 0.945\rceil$-th score: it is wider.
• Had it fallen inside: $\alpha_{t+1} = 0.1 + 0.05 \cdot 0.1 = 0.105$, slightly narrower.

Nine hits move $\alpha_t$ up by $9 \cdot 0.005 = 0.045$, and one miss moves it down by $0.045$: at a miss rate of exactly 10 %, $\alpha_t$ stays where it is.

5. In code #

NumPy and SciPy check the worked example; the notebook adds quantile gradient boosting and ACI.

import numpy as np
from scipy.stats import norm

def pinball(y, q, tau):
    u = y - q
    return np.maximum(tau * u, (tau - 1) * u)

# (a) the empirical minimiser of the mean pinball loss is a sample quantile
y = np.array([2., 4., 5., 7., 12.])
grid = np.linspace(0, 14, 1401)
L = [pinball(y, q, 0.7).mean() for q in grid]
print(grid[np.argmin(L)], np.quantile(y, 0.7, method="inverted_cdf"))   # 7.0 7.0

# (b) CRPS of N(0.3, 1.5^2) at y = 1, three ways
mu, s, yo = 0.3, 1.5, 1.0
z = (yo - mu) / s
closed = s * (z * (2 * norm.cdf(z) - 1) + 2 * norm.pdf(z) - 1 / np.sqrt(np.pi))
x = np.linspace(-15, 15, 300001)
integral = np.trapezoid((norm.cdf(x, mu, s) - (x >= yo)) ** 2, x)
taus = (np.arange(10000) + 0.5) / 10000
quant = 2 * pinball(yo, norm.ppf(taus, mu, s), taus).mean()
print(round(closed, 4), round(integral, 4), round(quant, 4))   # 0.4785 0.4785 0.4785

# (c) split conformal half-width
res = np.array([0.4, 1.1, 0.7, 2.3, 0.2, 1.6, 0.9, 3.0, 1.3])
n, alpha = len(res), 0.2
k = int(np.ceil((n + 1) * (1 - alpha)))
print(k, np.sort(res)[k - 1])                                  # 8 2.3

# (d) its coverage over many exchangeable data sets is k/(n+1) = 0.8
rng = np.random.default_rng(0)
scores = np.abs(rng.normal(size=(200000, n + 1)))   # n calibration scores + 1 test score
qhat = np.sort(scores[:, :n], axis=1)[:, k - 1]
print(round(np.mean(scores[:, n] <= qhat), 3))      # 0.801, close to k/(n+1) = 0.8

np.trapezoid needs NumPy 2.0 or newer (older versions call it np.trapz).

6. Where it shows up #

• Retail and supply chains. Stock levels are quantiles at the critical ratio. The M5 uncertainty competition (Makridakis et al., 2022) asked for nine quantiles of Walmart sales and scored them with a weighted scaled pinball loss. Many of the top teams used gradient-boosted trees (LightGBM) on lag features (lesson 08) with quantile losses.
• Energy. Load, wind and solar forecasts are quantile forecasts; the GEFCom2014 competition scored 99 quantiles with the pinball loss. Reserves and trading bids are set from upper quantiles.
• Weather and climate (ML for science). Ensemble forecasts are evaluated with the CRPS and with rank or PIT histograms. GenCast (Price et al., 2024), a diffusion model for global weather, reports the ensemble CRPS against the ECMWF ensemble.
• Epidemiology. The US COVID-19 Forecast Hub collected 23 quantiles per forecast and scored them with the weighted interval score, a quantile approximation of the CRPS (Bracher et al., 2021).
• Foundation models. Chronos, Chronos-Bolt, TimesFM and Moirai give quantiles or samples (lesson 09). Check their coverage on your own data with a backtest, and wrap them with conformal intervals if it is wrong.
• Software. Conformal intervals for forecasts exist in common libraries (for example MAPIE, and the conformal intervals of Nixtla's statsforecast and mlforecast). Read what the calibration set is: the guarantee depends on it.

7. Common mistakes #

• Judging intervals by coverage only. A very wide interval always covers. Report coverage and width, or use a proper score (CRPS, pinball loss, interval score).
• Trusting model intervals without a backtest. ARIMA and ETS intervals assume a correct model; check their empirical coverage over many origins.
• Adding quantiles over horizons. The 90 % quantile of a weekly total is not the sum of the daily 90 % quantiles. Use sample paths.
• Calibrating on training residuals. In-sample errors are too small, so the interval is too narrow. The calibration data must be new to the model.
• Shuffling time series before conformal calibration. Random splits mix future and past (leakage) and hide the drift. Calibrate on the most recent past and test on the future.
• Using np.quantile(scores, 1 - alpha) instead of the $\lceil (n + 1)(1 - \alpha)\rceil$-th smallest score. The difference is small for large $n$, but the guarantee needs the $n + 1$.
• Ignoring quantile crossing. Sort the predicted quantiles, or use a monotone parametrisation.
• Confusing the two meanings of $\alpha$. In this lesson $\alpha$ is the miscoverage, not the smoothing parameter.

8. Check yourself #

1. Show that the pinball loss at $\tau = 0.5$ is half the absolute error, and say which constant forecast minimises its mean.
   Answer

   For $u \ge 0$, $\rho_{0.5}(u) = 0.5u$; for $u < 0$, $\rho_{0.5}(u) = -0.5u$. So $\rho_{0.5}(u) = |u|/2$. The mean is minimised by the median (the 0.5-quantile), as for the absolute error.
2. An unsold newspaper costs 0.2 and a missing one costs 0.6. Which quantile of demand should you stock?
   Answer

   $\tau = c_u/(c_u + c_o) = 0.6/0.8 = 0.75$: the 75 % quantile.
3. The forecast is $\mathcal N(10, 2^2)$ and the outcome is $y = 10$. Compute the CRPS.
   Answer

   $z = 0$, so $\operatorname{CRPS} = 2\,[0 + 2\varphi(0) - 1/\sqrt\pi] = 2\,[0.7979 - 0.5642] = 0.467$. Even a perfect mean has a positive CRPS, about $0.234\,\sigma$.
4. Why is a forecast that always gives the long-run distribution calibrated, and why is it still poor?
   Answer

   The outcomes really come from the long-run distribution on average, so the PIT values are uniform and every interval has the right coverage. But it uses no current information, so its intervals are as wide as possible. A proper score such as the CRPS rewards a forecast that is calibrated *and* sharp.
5. You have $n = 9$ calibration scores. What is the smallest $\alpha$ for which split conformal gives a finite interval?
   Answer

   We need $k = \lceil 10(1 - \alpha)\rceil \le 9$, so $10(1 - \alpha) \le 9$, that is $\alpha \ge 0.1$. For $\alpha < 0.1$ the interval is infinite.
6. In ACI with $\alpha = 0.1$ and $\gamma = 0.01$, the method misses 5 times in a row. How much has $\alpha_t$ changed?
   Answer

   Each miss changes $\alpha_t$ by $\gamma(\alpha - 1) = -0.009$, so five misses lower it by $0.045$: from 0.1 to 0.055. The intervals become wider.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

1. Define the pinball loss $\rho_\tau(u)$ for the error $u = y - q$.

   Show answer

2. Which constant $q$ minimises the expected pinball loss $\mathbb E\,\rho_\tau(Y - q)$, and why?

   Show answer

3. What is the best stock level when an unsold unit costs $c_o$ and a missing unit costs $c_u$?

   Show answer

4. What is the PIT value of a forecast $F_t$, and what should its histogram look like?

   Show answer

5. What does a U-shaped PIT histogram mean?

   Show answer

6. What is the difference between calibration and sharpness?

   Show answer

7. Give the definition and the kernel form of the CRPS.

   Show answer

8. How is the CRPS related to the pinball loss?

   Show answer

9. In split conformal with $n$ calibration scores and miscoverage $\alpha$, which score is the half-width?

   Show answer

10. Write the update rule of adaptive conformal inference.

    Show answer

Add these cards to my reviewOpen the review →

9. Interview angle #

• "How do you get prediction intervals from a gradient-boosting forecaster?" — Train two more models with the quantile loss (for example $\tau = 0.05$ and $0.95$), sort the outputs to avoid crossing, and check the coverage in a rolling-origin backtest. Or keep the point model and add conformal intervals from recent backtest errors, one set of scores per horizon; use CQR to combine both.
• "Why does the pinball loss give a quantile?" — Derivative $F(q) - \tau$: a step up costs $1 - \tau$ for the outcomes below and saves $\tau$ for those above; the balance point has a share $\tau$ below.
• "How would you evaluate a probabilistic forecast?" — Calibration (coverage of several intervals, PIT histogram, also per subgroup) and a proper score (CRPS or the mean pinball loss over quantiles, scaled across series); compare with a seasonal naive forecast with bootstrapped intervals.
• "What does conformal prediction guarantee, and what breaks it for time series?" — Marginal coverage $\ge 1 - \alpha$ under exchangeability; dependence and drift break it; ACI restores a long-run average guarantee for any sequence, but not coverage at each time.

10. Practice #

The notebook for this lesson runs in Colab or in the browser. It has five exercises, each with a check cell and a folded solution:

1. Write pinball(y, q, tau) and check that its empirical minimiser over a grid is the sample quantile.
2. Write crps_gaussian(mu, sigma, y) and compare it with the sample kernel form and with twice the mean pinball loss over a grid of levels.
3. Train HistGradientBoostingRegressor(loss="quantile") on lag features for $\tau = 0.1, 0.5, 0.9$ and measure the coverage of the 80 % interval on a test period.
4. Write split_conformal(cal_scores, alpha) and check its coverage over 2000 simulated exchangeable data sets (between $1 - \alpha$ and $1 - \alpha + 1/(n + 1)$).
5. Write adaptive_conformal(scores, alpha, gamma) and compare its coverage after a variance jump with static split conformal.

11. Watch and read #

A Tutorial on Conformal PredictionAnastasios Nikolas Angelopoulos (with Stephen Bates) · 38 min

A hands-on tour of conformal prediction that follows the Angelopoulos and Bates paper in the reading list.

Quantile LossMeerkat Statistics · 6 min

The quantile loss as a formula and as a graph, and a proof that it gives the quantile.

Conformal Prediction for Time Series: Uncertainty Quantification for Trustworthy SystemsPython Italia · 28 min

Conformal prediction intervals for time-series forecasts in Python, including data that are not stationary.

GFDW02 | Mr. Romain Pic | Mathematical Properties of Continuous Ranked Probability Score ForecastingINI Satellite Events · 35 min

Advanced and optional; mathematical properties of the CRPS as a training and evaluation score.

• Hyndman and Athanasopoulos, Forecasting: Principles and Practice, 3rd ed., §5.5 distributional forecasts and prediction intervals and §5.9 evaluating distributional forecast accuracy (quantile scores, Winkler score, CRPS).
• Gneiting and Raftery (2007), Strictly proper scoring rules, prediction, and estimation, JASA 102: the CRPS, its kernel form and proper scores (author PDF).
• Gneiting and Katzfuss (2014), Probabilistic forecasting, Annual Review of Statistics and Its Application 1: calibration, sharpness and the PIT.
• Angelopoulos and Bates, A Gentle Introduction to Conformal Prediction and Distribution-Free Uncertainty Quantification (arXiv 2107.07511): split conformal, CQR, and the proofs.
• Gibbs and Candès (2021), Adaptive conformal inference under distribution shift (NeurIPS): ACI and the long-run guarantee.
• Romano, Patterson and Candès (2019), Conformalized quantile regression.
• Xu and Xie (2021), Conformal prediction interval for dynamic time-series (ICML; extended arXiv version Conformal prediction for time series): EnbPI; Angelopoulos, Candès and Tibshirani (2023), Conformal PID control for time series prediction.
• Koenker and Bassett (1978), Regression quantiles, Econometrica 46: the origin of quantile regression.