Skip to content

fix: asymmetric quantile score cancels the calibration-size guard - #978

Open
shivamlalakiya wants to merge 7 commits into
scikit-learn-contrib:masterfrom
shivamlalakiya:fix/quantile-score-effective-calibration-samples
Open

fix: asymmetric quantile score cancels the calibration-size guard#978
shivamlalakiya wants to merge 7 commits into
scikit-learn-contrib:masterfrom
shivamlalakiya:fix/quantile-score-effective-calibration-samples

Conversation

@shivamlalakiya

@shivamlalakiya shivamlalakiya commented Aug 4, 2026

Copy link
Copy Markdown
Contributor

Description

Follow-up to #958, which merged with a review still open. This picks up the blocking finding from that review and the smaller ones that came with it. All of them are still live on master, and all concern code that is new in #958 and unreleased as of v1.4.1, so nothing here can break an existing configuration.

Review threads: guard, suggested fix, sym, HISTORY.md, test range.

Four independent commits, so any one of them can be dropped without touching the others. The HISTORY.md de-duplication that came out of the same review is split off into #979, since it shares nothing with this code.

1. QuantileRegressionScore cancels the calibration-size guard (14e0e072)

get_signed_conformity_scores returns one row per side, of shape (2, n_samples). BaseRegressionScore.get_effective_calibration_samples counts non-NaN entries and halves the result for asymmetric scores:

n: int = np.sum(~np.isnan(scores))
if not self.sym:
    n //= 2

so the two-row layout returns 2n // 2 == n. The halving that accounts for each side being calibrated at alpha / 2 is exactly cancelled, and _check_alpha_and_n_samples validates alpha against n rather than n // 2.

What actually happens below the threshold

Both sides of the asymmetric branch resolve the same reference level: alpha_low = beta = alpha / 2 with reverse=True gives alpha_ref = 1 - alpha / 2, and alpha_up = 1 - alpha + beta = 1 - alpha / 2. So both clip together, as soon as

ceil((1 - alpha/2) * (n + 1)) / n  >=  1

When they do, and for fixed quantile predictions q_lo, q_hi, the interval collapses analytically:

  • upper: signed = +1, level 1, so the quantile is max(y_cal - q_hi) and bound_up = q_hi + max(y_cal - q_hi) = max(y_cal)
  • lower: signed = -1, level 1, so the quantile is min(y_cal - q_lo) and bound_low = q_lo + min(y_cal - q_lo) = min(y_cal)

q_lo and q_hi cancel. The interval degenerates exactly to [min(y_cal), max(y_cal)] — the range of the calibration sample. For exchangeable continuous draws that gives, distribution-free and without simulation,

coverage = 1 - 2/(n + 1) = (n - 1)/(n + 1)

which is below the requested 1 - alpha exactly when

(n - 1)/(n + 1) < 1 - alpha    <=>    n < (2 - alpha)/alpha

That is the correctness column below, derived rather than tabulated. Note the degeneracy region is wider than the bug: for 19 <= n <= 38 at alpha = 0.1 the interval is still the calibration range, but (n-1)/(n+1) >= 0.9 there, so it is merely conservative. The bug is the sub-region where the guard passes and the closed form falls short:

alpha guard needs n >= interval == [min, max] for n <= correctness needs n >= silent-clip window
0.05 20 78 39 20..38
0.10 10 38 19 10..18
0.20 5 18 9 5..8
0.30 4 12 6 4..5

Verified two ways. Structurally: running the real _compute_regression_quantile through the asymmetric branch with fixed predictions, bound_low == min(y_cal) and bound_up == max(y_cal) holds for every n the level condition predicts and fails at the first n past it, on all four alpha. Numerically, 1,000,000 trials of standard normal y:

n=15   sim 0.8752   (n-1)/(n+1) = 14/16 = 0.8750   <- guard passes, coverage short
n=18   sim 0.8948   (n-1)/(n+1) = 17/19 = 0.8947   <- guard passes, coverage short
n=19   sim 0.8999   (n-1)/(n+1) = 18/20 = 0.9000
n=20   sim 0.9052   (n-1)/(n+1) = 19/21 = 0.9048

matching the end-to-end figures in the original write-up of this PR (0.8750, 0.8939, 0.9005, 0.9043) and the independent reproduction by @LEDazzio01, who suggested the closed form. It is exact only for fixed quantile predictions; with fitted quantile regressors the bounds move and the statement reverts to a simulated one — but that was equally true of the simulation it replaces.

allow_infinite_bounds=True does not help either: the unbounded escape tests alpha_ref >= 1, not the corrected level. At n = 15, alpha = 0.1 the lower side asks for order statistic ceil(0.95 * 16) = 16 out of 15 available. SplitConformalRegressor with AbsoluteConformityScore(sym=False) already refuses the same n, the same alpha, and the same alpha / 2 split.

The fix reads the count from a single row, so an asymmetric quantile score now reports the same effective sample count as every other asymmetric score in MAPIE. That makes the guard reject n <= 19 at alpha = 0.1 — one sample more conservative than the exact requirement n >= 19, which is the same slack AbsoluteConformityScore(sym=False) already carries, rather than a new rule.

The dimension check is load-bearing rather than defensive: AbsoluteQuantileRegressionScore inherits the method and reduces both rows to a single one-dimensional score, so the default path must keep reading the array whole. That path is symmetric and unchanged.

The deeper fix — raising or returning +inf instead of clipping in _compute_regression_quantile — is #974 and out of scope here.

2. sym is not a knob on AbsoluteQuantileRegressionScore (44c4915a)

The score is max(y_pred_lower - y, y - y_pred_upper), a single distribution of absolute distances calibrating both bounds, so it is symmetric by construction. sym was settable anyway, and _check_score accepts any QuantileRegressionScore subclass, so AbsoluteQuantileRegressionScore(sym=False) reached the asymmetric branch of get_bounds. That branch expects one signed distribution per side; given a single non-negative one it takes the beta quantile for the lower bound and adds it, instead of subtracting a 1 - alpha quantile.

Measured on 200 calibration samples at a 0.9 confidence level, with symmetric Gaussian noise:

sym=True   coverage=0.9375  width=3.8130
sym=False  coverage=0.9834  width=5.4773

Conservative rather than under-covering, but 44% wider for no gain, and lopsided — the lower bound widens by 1.3533 against 0.3110 above, an artefact of the branch rather than of the data. QuantileRegressionScore keeps sym for callers who want the asymmetric variant.

Happy to swap this for an explicit ValueError on sym=False if you would rather keep the argument in the signature.

3. HISTORY.md (0482140b)

#958 documented none of its own changes, so this commit originally added an entry for CrossConformalizedQuantileRegressor and the two new quantile conformity scores, plus one for each fix above.

Since then the v1.5.0 pre-release (60031ac) restructured HISTORY.md, opening an empty 1.x.x section and writing its own CrossConformalizedQuantileRegressor entry under 1.5.0. On merging that in I dropped my now-redundant feature entry and kept only the two entries for the fixes in this PR, both under 1.x.x. The same restructure appears to supersede #979, which I will close unless a maintainer sees something left in it.

4. Inverted shapes in the get_signed_conformity_scores docstring (5b9e4675)

Flagged by @LEDazzio01 in the same comment. QuantileRegressionScore.get_signed_conformity_scores documented y as (2, n_samples), y_pred as (n_samples,) and the return as (n_samples, 2). All three are inverted — the method indexes y_pred[0] and y_pred[1] and vstacks the rows:

documented shapes  y (2, 6)  y_pred (6,)    -> returns (4, 6)
actual shapes      y (6,)    y_pred (2, 6)  -> returns (2, 6)

In scope because that row layout is precisely what commit 1 turns on, and commit 1 adds a docstring on get_effective_calibration_samples correctly describing the scores as (2, n_samples) — immediately below one asserting the opposite.

Type of change

  • Bug fix (non-breaking change which fixes an issue)

Removing sym from AbsoluteQuantileRegressionScore.__init__ is technically a signature change, but the class is new in #958 and absent from v1.4.1, so no released configuration can start raising on upgrade and no deprecation cycle is needed.

How Has This Been Tested?

Two new tests in mapie/tests/test_quantile_regression.py, both of which fail on master and pass with commit 1:

  • test_effective_calibration_samples_are_counted_per_side — pins the arithmetic: a (2, n) asymmetric layout must report the same effective count as a one-dimensional asymmetric score for the same n. Before: 15; after: 7, matching AbsoluteConformityScore(sym=False).
  • test_asymmetric_score_rejects_a_calibration_set_too_small_for_alpha_over_two — pins the user-visible behaviour end-to-end at n = 15, confidence_level = 0.9. Before: a finite interval built on a clipped level; after: ValueError.

The range is the point here. Every calibration set in #958's new tests sits well outside the silent-clip windows above — _asymmetric_scores uses 500, test_asymmetric_bounds_split_the_miscoverage_between_both_sides uses 5000, the end-to-end tests use X[:100] — which is why this fix moves no existing test result across the +2118 lines that PR added to the file.

  • pytest mapie/tests/test_quantile_regression.py mapie/tests/test_conformity_scores_bounds.py304 passed.
  • Full suite: make coverage exits clean, 100% coverage, with mapie/conformity_scores/bounds/quantile.py at 100% and all 4 branches covered.
  • Every coverage and width figure above was re-measured against this branch rather than quoted from the earlier review, and the degeneracy table was re-derived from the real _compute_regression_quantile rather than from algebra alone.
  • Merged current master (60031ac, the v1.5.0 pre-release) into the branch, resolving the HISTORY.md and AUTHORS.md conflicts against the restructured files. Merged rather than rebased so the commit hashes referenced above stay valid — happy to rebase and force-push instead if you prefer a linear history.

Checklist

Guidelines

Quality Checks

  • Linting passes successfully: make lint
  • Typing passes successfully: make type-check (green in the Quality checks job; on my machine mypy fails inside the numpy stubs for Python 3.13 before reaching mapie — verified clean with mypy mapie --python-version 3.13, which reports the same 5 pre-existing errors as master and none in the changed files)
  • Unit tests pass successfully: make tests
  • Coverage is 100%: make coverage
  • When updating documentation: doc builds successfully and without warnings: make doc — n/a, no changes under doc/; the only prose changes are class docstrings and HISTORY.md
  • When updating documentation: code examples in doc run successfully: make doctest — n/a, no doctest added or changed (--doctest-modules runs as part of make coverage and is green)

Contribution Documentation

LLM Usage

  • I used a Large Language Model (LLM) for this contribution.

  • I carefully reviewed and verified all LLM-generated content.

  • LLM used: Claude (Claude Code)

  • Purpose: Very limited, mechanical use of Claude Code to assist with tracing the conformity-score and quantile call paths, scaffolding the two tests, and drafting the boilerplate for this description. The core analysis and logic were human-driven. Every number here was produced by executing the code on this branch, and every claim was rigorously checked against the failing-before/passing-after tests rather than taken on the model's word.

`QuantileRegressionScore.get_signed_conformity_scores` returns one row per
side, of shape (2, n_samples). `BaseRegressionScore.get_effective_calibration_
samples` counts non-NaN entries and halves the result for asymmetric scores,
so the two-row layout returned 2n // 2 == n. The halving that accounts for
each side being calibrated at `alpha / 2` was exactly cancelled.

`_check_alpha_and_n_samples` therefore validated `alpha` against n instead of
n // 2, and accepted calibration sets in which the requested order statistic
does not exist. `_compute_regression_quantile` then clipped the corrected
level to 1.0 and both bounds fell back to the extreme order statistic, with no
warning. `allow_infinite_bounds=True` did not help either: the `unbounded`
escape tests `alpha_ref >= 1`, not the corrected level.

At n = 15 and alpha = 0.1 the lower side asks for order statistic
ceil(0.95 * 16) = 16 out of 15 available; isolating the conformal step on
exchangeable draws measured 0.8725 two-sided coverage against a requested
0.90. `SplitConformalRegressor` with `AbsoluteConformityScore(sym=False)`
already refuses the same n, alpha and alpha/2 split.

Reading the count from a single row restores the halving, so an asymmetric
quantile score now reports the same effective sample count as every other
asymmetric score. The default `AbsoluteQuantileRegressionScore` is symmetric
and reduces both rows to a single one-dimensional score, so it is unaffected.

The deeper fix — raising or returning +inf instead of clipping in
`_compute_regression_quantile` — is tracked in scikit-learn-contrib#974 and out of scope here.
The score is `max(y_pred_lower - y, y - y_pred_upper)`, a single distribution
of absolute distances that calibrates both bounds, so it is symmetric by
construction. `sym` was settable anyway, and `_check_score` accepts any
`QuantileRegressionScore` subclass, so `AbsoluteQuantileRegressionScore(sym=
False)` reached the asymmetric branch of `get_bounds`. That branch expects one
signed distribution per side; given a single non-negative one it takes the
`beta` quantile for the lower bound and adds it, instead of subtracting a
`1 - alpha` quantile.

Measured on 200 calibration samples at a 0.9 confidence level, with symmetric
Gaussian noise:

    sym=True   coverage=0.9375  width=3.8130
    sym=False  coverage=0.9834  width=5.4773

Conservative rather than under-covering, but 44% wider for no gain, and
lopsided — the lower bound widens by 1.3533 against 0.3110 above, an artefact
of the branch rather than of the data.

Both scores are new in scikit-learn-contrib#958 and unreleased as of v1.4.1, so no released
configuration can break. `QuantileRegressionScore` keeps `sym` for callers who
want the asymmetric variant.
np.zeros(n_calib)
)

assert asymmetric == n_calib // 2
)

assert asymmetric == n_calib // 2
assert symmetric == n_calib
scikit-learn-contrib#958 documented none of its own changes, so `CrossConformalizedQuantileRegressor`
and the two new quantile conformity scores had no entry. Adds one for them and
one for each fix in this branch.

The 27 duplicated entries scikit-learn-contrib#958 also introduced are removed separately in scikit-learn-contrib#979,
to keep this branch to the quantile scores.
@shivamlalakiya
shivamlalakiya force-pushed the fix/quantile-score-effective-calibration-samples branch from 85d76f8 to 0482140 Compare August 4, 2026 17:48
@LEDazzio01

LEDazzio01 commented Aug 5, 2026

Copy link
Copy Markdown
Contributor

Hi @shivamlalakiya — following on from my comment on #973, same disclaimer: outside contributor, so this is a comment and not a review. @allglc / @GBrelurut own the call.

I checked commit 1 rather than read it, and it holds up. Two things to add: a closed form that I think can replace the simulation, and a docstring defect sitting right next to the code being fixed.

Commit 1 verified

The cancellation is exactly as described. get_signed_conformity_scores ends in

return np.vstack((y - y_pred[0], y - y_pred[1]))

so the array is (2, n), and BaseRegressionScore.get_effective_calibration_samples does np.sum(~np.isnan(scores))2n — before n //= 2. The halving that should account for each side being calibrated at alpha / 2 is spent undoing the row layout instead.

I reproduced your coverage numbers by reimplementing the asymmetric branch of get_bounds against _compute_regression_quantile, 40k trials, alpha = 0.1:

n clipped my sim yours
15 yes 0.8776 0.8750
18 yes 0.8950 0.8939
19 no 0.9025 0.9005
20 no 0.9046 0.9043
200 no 0.8982 0.8991

Same picture, within Monte Carlo error.

The clip-window coverage has a closed form

With fixed quantile predictions — your setup — the clipped case collapses analytically. When ceil((1 - alpha/2)(n+1)) / n > 1 the corrected level is clipped to 1.0 on both sides, and:

  • upper: signed = +1, so the quantile is max(y_cal - q_hi), giving bound_up = max(y_cal)
  • lower: signed = -1, so the quantile is min(y_cal - q_lo), giving bound_low = min(y_cal)

Both q_lo and q_hi cancel. The interval is exactly [min(y_cal), max(y_cal)], so for exchangeable continuous draws

coverage = 1 - 2/(n+1) = (n-1)/(n+1)

distribution-free, no simulation required. At n = 15 that is 14/16 = 0.8750 — your reported figure to four decimals — and at n = 18, 17/19 = 0.8947.

It also gives your threshold directly:

(n-1)/(n+1) < 1 - alpha   <=>   n < (2 - alpha)/alpha

which is the "correctness needs n >=" column of your table, derived rather than tabulated. I checked the smallest n satisfying the closed form against (2-alpha)/alpha for alpha ∈ {0.05, 0.1, 0.2, 0.3} and got 39, 19, 9, 6 against 39, 19, 9, 5.667 — exact.

Worth swapping in, I think. A 200k-trial simulation invites "is the harness right?", whereas (n-1)/(n+1) is checkable by eye, and it makes the failure mode legible: in the window, the interval degenerates to the range of the calibration sample. That is a sharper statement of the bug than a coverage deficit.

The caveat is that this is exact only for fixed predictions. With fitted quantile regressors the bounds move and you would be back to simulation — but that is equally true of the numbers already in the description.

Docstring defect in the code being fixed

QuantileRegressionScore.get_signed_conformity_scores currently documents:

y: NDArray[float] of shape (2, n_samples)
y_pred: NDArray[float] of shape (n_samples,)
Returns: NDArray[float] of shape (n_samples, 2)

All three are wrong. The code indexes y_pred[0] and y_pred[1], so y_pred is the two-row array and y is one-dimensional; the return is (2, n_samples), not (n_samples, 2). Feeding it the documented shapes returns (4, n):

docstring shapes  y (2, 6)  y_pred (6,)    -> returns (4, 6)
actual shapes     y (6,)    y_pred (2, 6)  -> returns (2, 6)

This seems in scope: the bug you are fixing is a misreading of that exact row layout, and commit 1 adds a docstring that correctly says (2, n_samples) while the method directly above still says the opposite. Fixing them together stops the next reader hitting the same trap.

On commit 2

No opinion on removing sym versus raising ValueError — that is a maintainer call. The reasoning that max(y_pred_lower - y, y - y_pred_upper) is symmetric by construction looks right to me, and the lopsided widening you measured (1.3533 below against 0.3110 above) is what you would expect from feeding one non-negative distribution into a branch that wants two signed ones.

Thanks for the detailed write-up — the review-thread permalinks and the split into independent commits made this easy to check piece by piece.


Analysis assisted by Claude (Claude Code); I verified the numbers and reasoning.

`get_signed_conformity_scores` documented `y` as `(2, n_samples)`,
`y_pred` as `(n_samples,)` and the return as `(n_samples, 2)`. All three
are inverted: the method indexes `y_pred[0]` and `y_pred[1]` and vstacks
the two rows, so `y_pred` is the two-row array, `y` is one-dimensional
and the return is `(2, n_samples)`. Passing the documented shapes
`y (2, 6)`, `y_pred (6,)` returns `(4, 6)`; the real shapes return `(2, 6)`.

The row layout is exactly what the effective-calibration-sample fix in
14e0e07 turns on, and that commit added a docstring on
`get_effective_calibration_samples` correctly describing the scores as
`(2, n_samples)` — directly below one asserting the opposite.
…e-effective-calibration-samples

# Conflicts:
#	AUTHORS.md
#	HISTORY.md
@shivamlalakiya

Copy link
Copy Markdown
Contributor Author

Thanks @LEDazzio01 — both points taken, and both are now in the branch.

The closed form replaces the simulation. You're right that (n-1)/(n+1) is the better statement, and the sharper framing is yours: in the window the interval degenerates to the range of the calibration sample. I've rewritten section 1 around it and derived the "correctness needs n >=" column from n < (2 - alpha)/alpha instead of tabulating it.

Two things I checked while adopting it, one of which widens the claim:

Both sides really do clip together, which the degeneracy argument needs. It isn't obvious from the code — alpha_low = beta with reverse=True and alpha_up = 1 - alpha + beta look like different levels — but reverse=True maps alpha_low to alpha_ref = 1 - alpha/2, which is exactly alpha_up. Same reference level, so one threshold governs both bounds.

The degeneracy region is wider than the clip window. The interval is [min(y_cal), max(y_cal)] whenever the corrected level reaches 1, i.e. ceil((1 - alpha/2)(n+1)) >= n, not only when it exceeds 1. At alpha = 0.1 that is every n <= 38, not just the 10..18 window. Running the real _compute_regression_quantile through the asymmetric branch with fixed predictions, bound_low == min(y_cal) and bound_up == max(y_cal) holds for exactly the n that condition predicts and fails at the first n past it, on all four alpha:

alpha interval == [min, max] for n <= coverage (n-1)/(n+1) short of 1 - alpha for n < guard passes from n >=
0.05 78 39 20
0.10 38 19 10
0.20 18 9 5
0.30 12 6 4

So the bug window is the intersection of the last two columns — the table in the description is unchanged — but for 19 <= n <= 38 at alpha = 0.1 the interval is still the calibration range, merely conservative rather than invalid. That was worth knowing and I would not have looked without your comment.

The 1,000,000-trial check lands on the closed form: n = 15 -> 0.8752 against 14/16 = 0.8750, n = 18 -> 0.8948 against 17/19 = 0.8947. And your caveat is in the description: exact only for fixed predictions.

The docstring is fixed in 5b9e467, as a separate commit. Confirmed your reproduction — documented shapes y (2, 6), y_pred (6,) return (4, 6); the real shapes return (2, 6). Agreed it belongs here: the two-row layout is exactly what commit 1 turns on, and commit 1 was adding a correct (2, n_samples) docstring directly below one saying the opposite.

Also merged current master in — the v1.5.0 pre-release restructured HISTORY.md, which resolved the conflict and, incidentally, looks to have superseded #979.


Analysis assisted by Claude (Claude Code); every number above was produced by executing the code and I verified them.

Matches the entry added on the branch for scikit-learn-contrib#973, so the two do not conflict
when both land.
@allglc

allglc commented Aug 11, 2026

Copy link
Copy Markdown
Collaborator

@GBrelurut this change is about your cross-CQR PR, I will let you check it it makes sense. Also maybe check the comments on the cross-CQR PR discussion.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants