Michael Brenndoerfer | Data & AI, Private Equity, Technology

Michael Brenndoerfer

Quantile Regression › Visualizing Quantile Regression

Michael BrenndoerferOctober 26, 202611 min read

Machine Learning from ScratchAll chapters

Part 5 of 15Linear and Generalized Linear Models

Chapter 47Quantile Regression › Visualizing Quantile Regression

Reading section 7 of 10 · Visualizing Quantile Regression · 11 minSection 7 of 10 · 11 min

Section 7 of 10 · 11 min read

Quantile Regression › Visualizing Quantile RegressionSection view · Sections in this chapter7 / 10 ↓

Sections in this chapter10

Reading tools Definitions, links & notes

Choose your expertise level to adjust how many terms are explained. Beginners see more tooltips, experts see fewer to maintain reading flow. Hover over underlined terms for instant definitions.

Article links

Make inline references clickable

Notes & comments

Notes & commentsLeave notes and join the discussionSelection toolbarShow comment actions when selecting text. Turn off for the browser menu only.Comment indicatorsShow small counts beside annotated passagesPassage highlightsGently shade annotated text selections

A table of coefficients is a weak instrument for understanding a quantile fit. The signature picture is the fan: the data, with the fitted quantile curves drawn over them. Where the fitted fan opens up, the selected quantile gaps increase. Parallel fitted levels have constant gaps, though finitely many levels cannot establish the shape of the entire conditional distribution.

PythonIn [26]
CollapseExpandfrom book_plot_style import PALETTE, polish_axes, use_book_style

use_book_style()

fig, ax = plt.subplots(figsize=(6.0, 4.0))
ax.scatter(
    x_all[train_rows],
    y_tr,
    s=9,
    color=PALETTE["muted"],
    alpha=0.30,
    label="training data",
)
ax.fill_between(
    grid_x,
    grid_pred[0.05],
    grid_pred[0.95],
    color=PALETTE["primary"],
    alpha=0.18,
    label="5 to 95 percent band",
)
ax.fill_between(
    grid_x,
    grid_pred[0.25],
    grid_pred[0.75],
    color=PALETTE["amber"],
    alpha=0.28,
    label="25 to 75 percent band",
)
ax.plot(
    grid_x,
    grid_pred[0.50],
    color=PALETTE["muted"],
    linewidth=1.1,
    label=r"median, $\tau=0.5$",
)
ax.set_title("Fitted quantile curves form a widening fan")
ax.set_xlabel("Predictor $x$")
ax.set_ylabel("Response $y$")
ax.legend(loc="upper left")
polish_axes(ax)
plt.show()

VisualizationOut [26]
CollapseExpand

Scatter of six hundred points rising from left to right with increasing vertical spread, overlaid with a fan of shaded bands that are narrow near x equals zero and wide near x equals ten, plus a central solid median line. — Fitted quantile curves on the training data. The shading marks the central 50 and 90 percent bands, and the solid line is the median fit. The bands widen with $x$ because the response spread widens with $x$, so quantile regression captures the heteroscedasticity directly, whereas the separate constant-variance plug-in reference used below has one band width. The fan shape makes the change in conditional spread visible without fitting a separate variance model.

The response generator uses symmetric Gaussian noise, so its true conditional mean and median coincide. Any difference between their fitted lines is a finite-sample result. In the executed fit, the 90 percent band width grows from about 1.69 at $x=0$ to 17.13 at $x=10$ ; the 50 percent width grows from about 0.98 to 6.52. Those estimates are separate from the illustrative coefficients in the hand calculation.

The width of the bands is linear in $x$ , and it is worth isolating that fact. Plotting width against $x$ makes the comparison with a constant-variance model concrete.

PythonIn [27]
CollapseExpandfitted_width_90 = grid_pred[0.95] - grid_pred[0.05]
fitted_width_50 = grid_pred[0.75] - grid_pred[0.25]
ols_plugin_width = 2.0 * ols_halfwidth

PythonIn [28]
CollapseExpandfrom book_plot_style import PALETTE, polish_axes, use_book_style

use_book_style()

fig, ax = plt.subplots(figsize=(6.0, 4.0))
ax.plot(
    grid_x,
    fitted_width_90,
    color=PALETTE["primary"],
    linewidth=1.1,
    label="quantile band width, 5 to 95",
)
ax.plot(
    grid_x,
    fitted_width_50,
    color=PALETTE["amber"],
    linewidth=1.1,
    linestyle=":",
    label="quantile band width, 25 to 75",
)
ax.axhline(
    ols_plugin_width,
    color=PALETTE["muted"],
    linewidth=1.8,
    linestyle="--",
    label="OLS band width (constant variance)",
)
ax.set_title("Fitted quantile widths and a constant-width reference")
ax.set_xlabel("Predictor $x$")
ax.set_ylabel("Band width in response units")
ax.legend(loc="upper left")
polish_axes(ax)
plt.show()

VisualizationOut [28]
CollapseExpand

Three lines plotted against x: a steeply rising line for the 90 percent quantile band width reaching about seventeen, a rising 50 percent width reaching about 6.52, and a flat dashed plug-in reference at about 10.49. — Width of the fitted bands as a function of the predictor. The fitted 5 to 95 percent width grows from about 1.69 to 17.13 units, and the 25 to 75 percent width grows from about 0.98 to 6.52, while a separate constant-variance normal plug-in reference stays fixed at about 10.49. The widening widths show the estimated conditional spread growing with $x$, in contrast with the flat reference line. These are computed fitted summaries, not finite-sample coverage guarantees.

The dashed line uses one residual standard deviation for every predictor value. For this heteroscedastic generator it is too wide relative to the true conditional band at small $x$ and too narrow at large $x$ . An overly wide band includes more ordinary responses; an overly narrow one excludes more. The fitted quantile widths vary with $x$ , but their finite-sample accuracy still needs checking. An OLS mean model could instead be paired with a separate conditional-variance model.

Quantile crossing, interval coverage, and limitationsLink Copied

Each $\tau$ was fit by a separate linear program, and none of those programs knew about the others. Nothing in our procedure enforces $\hat Q(\tau_1 \mid x) \le \hat Q(\tau_2 \mid x)$ for $\tau_1 < \tau_2$ . Even this 600-row fit can have adjacent-level crossings near the edge of the predictor range. The displayed 50 and 90 percent band endpoints are separate comparisons; their ordering does not establish order among every fitted level.

The algebra of a crossing is simple. Two fitted lines

\hat Q(\tau_1 \mid x) = \hat\beta_0(\tau_1) + \hat\beta_1(\tau_1) x, \qquad \hat Q(\tau_2 \mid x) = \hat\beta_0(\tau_2) + \hat\beta_1(\tau_2) x

intersect where the two right-hand sides are equal. Collecting terms in $x$ gives

\left[\hat\beta_0(\tau_2) - \hat\beta_0(\tau_1)\right] + \left[\hat\beta_1(\tau_2) - \hat\beta_1(\tau_1)\right] x = 0,

so the crossing point is

x^\star = \frac{\hat\beta_0(\tau_1) - \hat\beta_0(\tau_2)}{\hat\beta_1(\tau_2) - \hat\beta_1(\tau_1)},

where:

$x^\star$ : the value of the predictor at which the two fitted quantile lines intersect
$\tau_1 < \tau_2$ : the two quantile levels whose curves are compared
$\hat\beta_0(\tau_k)$ : the estimated intercept at level $\tau_k$
$\hat\beta_1(\tau_k)$ : the estimated slope at level $\tau_k$

provided the slopes differ. If the slopes are equal the lines are parallel and either never cross or coincide. For an example, suppose two independently fitted lines returned $\hat Q(0.4 \mid x) = 3.0 + 1.0x$ and $\hat Q(0.6 \mid x) = 3.4 + 0.9x$ . Then $x^\star = (3.0 - 3.4)/(0.9 - 1.0) = (-0.4)/(-0.1) = 4$ , and, consistently, substituting $x = 4$ into both lines gives $3.0 + 1.0 \times 4 = 3.4 + 0.9 \times 4 = 7.0$ . To the left of $x = 4$ the 40 percent curve sits below the 60 percent curve, as it must; to the right the ordering reverses, which is impossible for a valid conditional distribution since $Q(0.4 \mid x) \le Q(0.6 \mid x)$ must hold for every $x$ .

The practical remedy is rearrangement. At any fixed $x$ , take the raw fitted values across all estimated levels and sort them into increasing order. Formally, if we have fitted curves at levels $\tau_1 < \dots < \tau_K$ , let $\hat Q_{\text{raw}}(\tau_k \mid x)$ denote the raw fitted value at level $\tau_k$ . Define the rearranged value $\hat Q_{\text{re}}(\tau_k \mid x)$ as the $k$ -th smallest of the raw fitted values $\{\hat Q_{\text{raw}}(\tau_1 \mid x), \dots, \hat Q_{\text{raw}}(\tau_K \mid x)\}$ at that $x$ , so that $\hat Q_{\text{re}}(\tau_1 \mid x) \le \dots \le \hat Q_{\text{re}}(\tau_K \mid x)$ . In words, $\hat Q_{\text{re}}(\tau_k \mid x)$ is the value that a sorted version of the raw fits places in $k$ -th position.

where:

$\hat Q_{\text{raw}}(\tau_k \mid x)$ : the raw, independently fitted value at level $\tau_k$ and covariate $x$ , before any reordering
$\hat Q_{\text{re}}(\tau_k \mid x)$ : the rearranged value at level $\tau_k$ and covariate $x$ , obtained by sorting the raw fits across levels at that $x$
$K$ : the number of estimated quantile levels

Sorting happens separately at each fixed $x$ , so rearrangement changes only the assignment of raw values to levels, not the set of values themselves. Sorting is monotone by construction, so the rearranged curves cannot cross. Chernozhukov, Fernández-Val and Galichon (2010) establish finite-sample improvements in function-distance measures for rearrangement toward a monotone target. This does not promise a better estimate at every individual level or a finite-sample coverage guarantee. Here sorting is used to repair the order of our finite set of fitted values. Rearrangement is a patch, not a cure: the underlying problem is that independent estimation ignores the coupling between levels, and a joint estimator with explicit noncrossing constraints avoids crossings on the predictor domain covered by those constraints, at the cost of a larger coupled optimization problem.

PythonIn [29]
CollapseExpandsmall_taus = np.array([0.20, 0.40, 0.50, 0.60, 0.80])
crossing_rows, small_fits = [], {}
for n_small in (25, 50, 100, 200):
    B_small = np.column_stack(
        [fit_quantile_lp(X_tr[:n_small], y_tr[:n_small], t) for t in small_taus]
    )
    P_small = G @ B_small
    small_fits[n_small] = P_small
    crossing_rows.append(
        {
            "training rows": n_small,
            "crossing adjacent-level pairs": int(
                np.sum(np.diff(P_small, axis=1) < 0)
            ),
        }
    )
crossing_table = pd.DataFrame(crossing_rows)
crossing_table

Out[30]:

	training rows	crossing adjacent-level pairs
0	25	85
1	50	24
2	100	5
3	200	7

Each row counts the (grid point, adjacent level) pairs where the fitted curve at the lower level sits above the fitted curve at the higher level. These are diagnostics for four nested subsets, not a proof that crossings decrease whenever sample size increases. Finite-sample uncertainty and misspecification can both produce crossings. Rearranging the 25-row fit removes every crossing by construction:

PythonIn [31]
CollapseExpandraw_25 = small_fits[25]
sorted_25 = np.sort(raw_25, axis=1)
print(
    f"crossings before rearrangement: {int(np.sum(np.diff(raw_25, axis=1) < 0))}"
)
print(
    f"crossings after  rearrangement: {int(np.sum(np.diff(sorted_25, axis=1) < 0))}"
)
print(
    f"average absolute change from sorting: {np.mean(np.abs(sorted_25 - raw_25)):.3f} response units"
)

ConsoleOut [32]
CollapseExpandcrossings before rearrangement: 85
crossings after  rearrangement: 0
average absolute change from sorting: 0.017 response units

The second count is zero for a reason that requires no data at all: sorting guarantees it. The third number tells us the price of the repair. When the raw curves barely cross, rearrangement changes almost nothing; when they cross badly, it moves the fitted values and the resulting loss changes accordingly.

PythonIn [33]
CollapseExpandline_lo = 3.0 + 1.0 * grid_x
line_hi = 3.4 + 0.9 * grid_x
cross_x = (3.0 - 3.4) / (0.9 - 1.0)
cross_y = 3.0 + cross_x

PythonIn [34]
CollapseExpandfrom book_plot_style import PALETTE, polish_axes, use_book_style

use_book_style()

fig, ax = plt.subplots(figsize=(6.0, 4.0))

ax.plot(
    grid_x,
    line_lo,
    color=PALETTE["primary"],
    linewidth=1.1,
    label=r"$\hat Q(0.4 \mid x)$",
)
ax.plot(
    grid_x,
    line_hi,
    color=PALETTE["amber"],
    linewidth=1.1,
    linestyle="--",
    label=r"$\hat Q(0.6 \mid x)$",
)
ax.scatter(
    [cross_x],
    [cross_y],
    s=70,
    color=PALETTE["muted"],
    zorder=5,
    label=r"crossing at $x^\star = 4$",
)
ax.set_title("Independently fitted quantile lines can cross")
ax.set_xlabel("Predictor $x$")
ax.set_ylabel("Fitted quantile")
ax.legend(loc="upper left")
polish_axes(ax)
plt.show()

VisualizationOut [34]
CollapseExpand

Two straight lines, one for tau 0.4 and one for tau 0.6, intersecting at x equals 4, with the lower quantile rising above the upper one to the right of the crossing, marked with a dot. — Two independently fitted quantile lines that cross at $x = 4$, a configuration no genuine conditional distribution can produce. To the left of the crossing the 40 percent curve sits below the 60 percent curve, and to the right the ordering reverses, with the 40 percent curve above the 60 percent curve. Rearrangement at each $x$ sorts the fitted values across levels and removes such crossings by construction, though it changes the fitted values themselves.

Now to coverage, which is the other place nominal statements can mislead. A band labelled "90 percent" is a claim about new observations: if we repeatedly draw $x$ and $y$ from the population and check whether $y$ falls inside $\left[\hat Q(0.05 \mid x), \hat Q(0.95 \mid x)\right]$ , we would like the long-run frequency to be 0.90. Our held-out estimate was close to nominal, but that is a result to be measured, not assumed. The per-level calibration plot makes the check visual. For each level $\tau$ , plot the empirical fraction of held-out responses at or below the fitted $\tau$ -quantile against $\tau$ itself. A calibrated fit lies on the diagonal.

PythonIn [35]
CollapseExpandfrom book_plot_style import PALETTE, polish_axes, use_book_style

use_book_style()

fig, ax = plt.subplots(figsize=(6.0, 4.0))
ax.plot(
    [0.0, 1.0],
    [0.0, 1.0],
    color=PALETTE["muted"],
    linewidth=1.2,
    linestyle="--",
    label="perfect calibration",
)
ax.plot(
    taus,
    eval_table["test share below fit"],
    marker="o",
    color=PALETTE["primary"],
    linewidth=1.0,
    markersize=6,
    label="fitted quantile curves",
)
ax.set_title("Held-out calibration of the quantile curves")
ax.set_xlabel(r"Nominal level $\tau$")
ax.set_ylabel("Share of test responses below the fitted quantile")
ax.legend(loc="upper left")
polish_axes(ax)
plt.show()

VisualizationOut [35]
CollapseExpand

Scatter of seven points close to a dashed forty-five degree diagonal line, with nominal level on the horizontal axis and observed share below the fitted quantile on the vertical axis, ranging from zero to one. — Calibration of the fitted quantile curves on held-out data. Each point compares a nominal level with the empirical fraction of test responses falling at or below the corresponding fitted quantile. Points near the diagonal indicate approximate marginal calibration in this sample. Departure from the diagonal does not identify its cause; assess sampling uncertainty and check conditional fit and crossings separately. This is a marginal check and does not guarantee calibration conditional on every predictor value.

The plotted fractions are close to the diagonal in this held-out sample. The measured central 90 percent coverage is 0.870 and central 50 percent coverage is 0.480. Neither measurement supplies a finite-sample guarantee or establishes calibration conditional on every $x$ : overcoverage in one predictor region can cancel undercoverage elsewhere. Misspecification, limited tail information, crossing and distribution changes can all affect performance. Standard split-conformal calibration can provide marginal finite-sample coverage under exchangeability when its calibration data are reserved appropriately; it does not automatically give conditional coverage at every predictor value.

A quantile band describes the variability of individual responses. A confidence band for a quantile curve concerns uncertainty about that fitted function. Under regularity conditions estimation uncertainty can decrease with sample size while response variability remains. A restricted linear quantile fit may also be misspecified. A marginal calibration plot can expose some errors, but it need not reveal every conditional error or produce a particular S-shape.

Comments

Add to the conversation

A question, a correction, or an example from your own work.

No comments yet. Be the first to share your thoughts!

Quantile crossing, interval coverage, and limitationsLink Copied

Comments

Add to the conversation

Stay up to date