Isotonic regression answers one question: of all non-decreasing sequences, which is closest to this data in least squares? Not approximately ā exactly, in one left-to-right pass.
š Live, runs in your browser: https://dev48.in/ml/day82-isotonic-regression.html
Its main job in practice is calibration: turning scores that rank well into numbers that mean "probability".
It cannot touch your ranking metrics
Calibration is a monotone map. A monotone map cannot reorder anything. So AUC is untouched by construction:
| scores | AUC | Brier, raw | Brier, calibrated | improvement |
|---|---|---|---|---|
| squash 1.00 | 0.8186 | 0.1729 | 0.1739 | ā0.6% |
| squash 0.60 | 0.8186 | 0.1844 | 0.1739 | 5.7% |
| squash 0.35 | 0.8186 | 0.2047 | 0.1739 | 15.1% |
| squash 0.20 | 0.8186 | 0.2217 | 0.1739 | 21.6% |
| squash 0.10 | 0.8186 | 0.2351 | 0.1739 | 26.0% |
0.8186 in every row, to four decimals, while Brier moves 26%. The AUC column is not a coincidence to be checked, it is a theorem to be relied on.
This cuts both ways. If you are reporting AUC to show that calibration helped, you are reporting the one number that cannot respond. And the squash-1.00 row shows calibration slightly hurting when the scores were already fine.
It overfits, and nothing in it stops that
This is the part the tutorials leave out. Isotonic regression is an exact projection onto a convex cone. It is not a smoother. Nothing in it regularises.
| training size | blocks fitted | Brier in-sample | Brier out-of-sample | ratio |
|---|---|---|---|---|
| 40 | 3 | 0.1350 | 0.1650 | 1.22Ć |
| 100 | 7 | 0.1588 | 0.1961 | 1.23Ć |
| 400 | 10 | 0.1671 | 0.1764 | 1.06Ć |
| 2000 | 17 | 0.1745 | 0.1595 | 0.91Ć |
Given a small validation set it will happily fit a step per point. That is not a bug ā it is the method doing exactly what it was asked to do.
Optimality is characterised, not just computed
The useful thing about a projection is that you can check the answer without re-running the algorithm. Two conditions:
# 1. inside every level set, the fit IS the weighted mean of the data there
# (otherwise moving that block toward the mean lowers the error
# without breaking monotonicity - so it was not optimal)
# 2. the residual is orthogonal to the fit: sum_i w_i (y_i - f_i) f_i = 0
Both hold exactly, and neither calls PAVA. If your implementation passes these on your data, it is right on your data ā independent of whether you trust the merge logic.
Verified by a completely different algorithm
PAVA is a merge procedure: a stack of blocks, pooled left to right. The verifier never pools anything. It uses the other classical characterisation ā the isotonic fit is the slope of the greatest convex minorant of the cumulative-sum diagram:
// points P_0 = (0,0), P_i = (W_i, S_i); take the lower convex hull;
// the fit is the slope of the hull segment covering each index
for (let i = 0; i <= n; i++) {
while (hull.length >= 2 && cross(hull[at(-2)], hull[at(-1)], i) <= 0) hull.pop();
hull.push(i);
}
A convex hull, not a merge. Both were then checked against a brute-force search over all 2^(nā1) block partitions on 55 sequences ā including constants, strict reversals, singletons, negatives and weighted cases.
Largest SSE gap vs brute force: 0.00e+0.
One assertion worth keeping: pooling on >= rather than > is why four equal points become one block and not four. The fit is identical either way ā but the block count is what the overfitting table reads, so the wrong comparison gives you a correct curve with a meaningless x-axis.
357 verifier asserts, 0 failures. 16 in-page checks.
Top comments (0)