Inverse maps¶

Data these snippets assume

Blocks on this page use a fixed set of names for held-out data instead of re-deriving it each time. To run them, define the names first:

# docs: no-run — the definitions the snippets on this page assume
import numpy as np
from probcal import BetaCalibrator, Masterscale, make_pd_portfolio
from probcal.monitor import CalibrationMonitor

cal_set, new_set = make_pd_portfolio(n=3000, random_state=0), make_pd_portfolio(n=1000, random_state=1)
s_cal, y_cal = cal_set.scores, cal_set.y          # held-out calibration scores and outcomes
w_cal = np.ones_like(y_cal)                       # uniform sample weights
s_new = new_set.scores                            # scores of new obligors
ms = Masterscale.from_edges([0.01, 0.05], names=["G1", "G2", "G3"])   # the masterscale
p_cal = BetaCalibrator().fit(s_cal, y_cal).predict_proba(s_cal)         # calibrated PDs
grades = ms.assign(p_cal)                                              # rating labels
segments = np.array(["seg-a", "seg-b", "seg-c"])[np.arange(len(s_cal)) % 3]  # segment labels
model = ...   # any object with predict_proba(X); the docs use a stub that reads X[:, 0]
mon = CalibrationMonitor(alpha=0.05)              # with a few batches of a calibrated forecast applied

Every block also names, in its first comment line, which of these it uses.

Once decisions are made on calibrated probabilities, every consumer that lives on the raw score needs a translation service. A credit policy says "approve below 2% PD"; the deployed model emits raw margins; the cutoff engine, the masterscale, and the counterfactual generator all operate upstream of the calibrator. The capability to translate calibrated intervals back to raw-score intervals belongs to the calibrator itself, since only it knows its functional form, its block structure, and its output range. probcal builds it into the base contract: interval_inverse(lo, hi, *, space, buffer_logit), with thresholds.py providing thin functional wrappers.

The preimage identity¶

For a monotone calibration map \( g \), inverting a decision interval commutes with composition:

\[ \{\,x : g(f(x)) \in [\mathrm{lo}, \mathrm{hi}]\,\} \;=\; \{\,x : f(x) \in [\,g^{-1}(\mathrm{lo}),\; g^{-1}(\mathrm{hi})\,]\,\}. \]

The consequence reaches further than notation: calibration does not change counterfactual geometry, only the target interval. A counterfactual engine searching for the smallest feature change that lands \( f(x) \) in a raw interval works unchanged after calibration is deployed; it just needs the translated interval. This is the designed hand-off to treecf: probcal computes lo_z, hi_z = cal.interval_inverse(0.0, 0.02, space="logit") and the engine consumes Target.raw(range=(lo_z, hi_z)). The division of responsibilities is deliberate: probcal owns the capability and publishes the duck-typed protocol (interval_inverse plus is_monotone_); target-construction ergonomics belong to the consumer.

One trap deserves its single sentence here and a longer one in the FAQ: after calibration is deployed, a consumer-side Target.probability(...) silently inverts the model's own sigmoid link rather than the calibrator, and therefore targets the uncalibrated probability.

Generalized inverses and plateau semantics¶

Step calibrators have no literal inverse, so the contract is defined by generalized inverses: for non-decreasing \( g \),

\[ \mathrm{raw\_lo} = \inf\{\, s : g(s) \ge \mathrm{lo} \,\}, \qquad \mathrm{raw\_hi} = \sup\{\, s : g(s) \le \mathrm{hi} \,\}, \]

with \( \mathrm{lo} = 0 \) and \( \mathrm{hi} = 1 \) mapping to \( -\infty \) and \( +\infty \). For an isotonic map this means the preimage of a target begins exactly at the left edge of the first qualifying block, a semantics the tests pin down explicitly. Implementations follow the map's structure: closed form for Platt, temperature and the offset; monotone bisection for beta and monotone splines; searchsorted on the block structure for isotonic, CIR, binning, scaling-binning and the scalarized Venn–Abers (whose monotonicity the distribution-free chapter establishes). The wrapper composes right-to-left through the pipeline: the offset inverts first, subtracting \( \delta \) on the logit, then the calibrator's inverse applies.

Two refusals are part of the contract. If \( [\mathrm{lo}, \mathrm{hi}] \) does not intersect the calibrator's output range (routine for isotonic maps on low-PD data, whose range is the span of block means), the method raises UnattainableTargetError naming both intervals; silent clamping would convert a policy error into a wrong cutoff. And a non-monotone calibrator (is_monotone_ = False, ENIR being the resident example) has preimages that may be unions of intervals; the method raises NotImplementedError with an explanation, and the documented recommendation is simply to use a monotone calibrator when recourse or thresholding is downstream.

Robustness: plateaus and drift¶

A counterfactual that lands just past a block edge of a step calibrator is fragile twice over: the calibrated value jumps discretely at the edge, and any refit moves the edge. The first fragility argues for continuous calibrators, beta or CIR, wherever recourse is in scope. The second is addressed by buffer_logit: the calibrated interval is shrunk by a margin in logit space before inversion, so the produced raw interval is conservative by that margin. The guarantee it buys is concrete because the offset is a pure logit translation: a central-tendency update of magnitude at most \( m \) cannot invalidate a counterfactual built with buffer_logit = m. Larger buffers give tighter raw intervals, a monotonicity the tests check, and the trade is explicit: resistance to recalibration drift, paid in recourse difficulty. The concern is the one the algorithmic-recourse literature formalizes: recourse recommendations invalidated by model updates and distribution shift (Rawal, Kamar and Lakkaraju, 2020), and the case for building recourse robust to such shifts (Upadhyay, Joshi and Lakkaraju, 2021).

The masterscale workflow¶

Rating systems define grades on calibrated PD, a masterscale \( \{\text{grade } j : [\mathrm{lo}_j, \mathrm{hi}_j)\} \), while the model emits raw margins, with calibration sitting between. calibrated_bands_to_raw(calibrator, bands, ...) maps the entire masterscale to raw intervals in one call, which is the canonical workflow: grade edges translated once per recalibration, consumed by cutoff engines and counterfactual targeting alike (the output plugs directly into band-style targets on the raw scale). Since grade boundaries are policy artifacts that outlive model versions, the translation, not the masterscale, is what changes when the calibrator is refitted, and the audit story stays clean: policy fixed, mapping versioned.

In probcal¶

# s_cal, y_cal: held-out calibration scores and outcomes
from probcal import BetaCalibrator, UnattainableTargetError, calibrated_bands_to_raw

cal = BetaCalibrator().fit(s_cal, y_cal)

lo_s, hi_s = cal.interval_inverse(0.0, 0.02)               # "PD <= 2%" in score space
lo_z, hi_z = cal.interval_inverse(0.0, 0.02, space="logit")  # ... in raw margins

# Robust to the next re-anchoring of magnitude <= 0.1 log-odds:
lo_z, hi_z = cal.interval_inverse(0.0, 0.02, space="logit", buffer_logit=0.1)

masterscale = {"A": (0.0, 0.01), "B": (0.01, 0.03), "C": (0.03, 1.0)}
raw_bands = calibrated_bands_to_raw(cal, masterscale, space="logit")

try:
    cal.interval_inverse(0.95, 1.0)     # beyond an isotonic map's output range
except UnattainableTargetError as err:
    print(err)                           # named intervals — never a silent clamp

Point inverses: an exact preimage for the parametric family¶

interval_inverse answers "what raw range maps into this calibrated interval"; a second question is narrower and, for the affine and beta families, answerable exactly: "what raw score maps to this one calibrated probability". point_inverse(p, *, space) is that exact preimage: a single closed-form (or fixed, certified-precision) computation, not bisection.

It is defined on:

BaseCalibrator (inherited by every calibrator): the affine-logit path. Any calibrator whose affine_logit_coeffs_ is not None (PlattCalibrator, TemperatureCalibrator, and BetaCalibrator's tied "a"/"ab" variants) inverts logit g(s) = a * logit(s) + b in one line: z = (logit(p) - b) / a.
BetaCalibrator (all three variants, overriding the base method): the full "abm" map is not affine on the logit scale, so it gets its own layered exact construction; see below.
LogitOffset: z = logit(p) - delta, the same closed form as its interval_inverse.

A non-monotone calibrator, or a monotone one with no affine-logit or beta closed form (isotonic/CIR/binning/scaling-binning, splines, Venn–Abers, ENIR), raises NotImplementedError naming interval_inverse, which is the right tool wherever the map is a step function (plateaus have no single preimage) or an otherwise non-affine monotone curve (only a generalized inverse is well-defined).

The step-calibrator contract (what a recourse engine may rely on)¶

For the step calibrators (IsotonicCalibrator, HistogramBinningCalibrator), interval_inverse implements the exact generalized-inverse rule below, pinned by regression tests (tests/test_calibrator_protocol.py) because counterfactual engines (treecf) treat it as a contract:

Left inverse inf{s : g(s) ≥ lo}: the left edge of the first block/bin whose level reaches lo.
Right inverse sup{s : g(s) ≤ hi}: the boundary after the last block/bin whose level stays within hi (or the raw score 1.0 when that is the last block).
Closed intervals. A target equal to a plateau level qualifies the whole plateau: interval_inverse(t, t) at a block level t returns the block's full raw extent, closed on the left.
Constant extension. Beyond the outer block levels the map is constant, so a one-sided target (lo = 0 or hi ≥ the top level) maps to the full raw range on that side: 0.0/1.0 in probability space, ∓inf in logit space. An interval that misses the output range entirely raises UnattainableTargetError, never a clamp.

Worked example: an isotonic fit with block levels (0.01, 0.04, 0.04, 0.11) over blocks [0, 0.2) [0.2, 0.5) [0.5, 0.8) [0.8, 1]. The target t = 0.04 gives interval_inverse(0.04, 0.04) = (0.2, 0.8): the left edge of the first 0.04 block and the boundary after the second. A counterfactual that must land at 0.04 calibrated therefore lands anywhere in [0.2, 0.8), and the cheapest feasible point is the block boundary nearest the applicant. interval_inverse(0.12, 1.0) raises: 0.12 exceeds every level.

CenteredIsotonicCalibrator is not a step function; it interpolates linearly between block centers, so a block level is attained at a single point unless adjacent block means tie: interval_inverse(t, t) returns a point preimage. That distinction is also pinned by test.

The probability boundary¶

Two refusals keep point_inverse inside the no-silent-clamp doctrine. Targets are validated strictly inside (0, 1) before any clipping: p = 0 and p = 1 are not attained by any finite raw score, so the whole call raises UnattainableTargetError (naming the offending values) rather than clipping to [1e-12, 1 − 1e-12] and returning a finite "inverse". That is the same all-or-nothing convention the degenerate beta branches already used. And when space="probability", any root with |z| > logit(1 − 1e-12) ≈ 27.631 is refused with a pointer at space="logit": σ(z) beyond that bound rounds into the clipping zone, predict_proba clips it back, and the round trip breaks silently (target 0.998 → raw 1.0 → forward 0.17 was the failure this closes). On the logit scale the answer is exact, so nothing is lost; only the lossy representation is refused. For affine-logit maps the bound aligns point_inverse with interval_inverse exactly: p > g(1 − 1e-12) ⟺ z > logit(1 − 1e-12).

The beta inverse: a layered exact construction¶

BetaCalibrator.point_inverse solves, with z = logit(s) and K = logit(p) - c,

\[ h(z) = a z + (b - a)\,\mathrm{softplus}(z) = K, \qquad \mathrm{softplus}(z) = \ln(1 + e^{z}), \]

which has no elementary closed form when a != b (it is transcendental; see the Lambert-W connection below). The construction layers three ideas, each exact in a different regime:

Layer 1: minimax hyperbola seed. Replacing softplus(z) by the minimax hyperbola (z + sqrt(z^2 + kappa)) / 2 with kappa = 1.524 (max pointwise deviation 0.076) turns the equation quadratic, with admissible root

\[ z_0 = \frac{(a+b) K - (b-a)\sqrt{K^2 + \kappa a b}}{2 a b}, \]

exact at a = b (collapses to K / a) and in both tails, with error bounded by 0.076 * |b - a| / min(a, b) elsewhere.

Layer 2: certified Halley correction. Up to 4 fixed Newton–Halley steps refine z_0 to machine precision, exiting as soon as the residual certificate |h(z) - K| <= 1e-13 * max(1, |K|) is met: a bounded, finite expression, not open-ended iteration. The residual is itself a certificate: |z - z*| <= |h(z) - K| / min(a, b).

Layer 3: the Lambert-W tail form (theory, not shipped code). As |z| -> infinity the equation admits a closed form in the Lambert-W function. For the left tail (z -> -infinity, softplus(z) -> 0):

\[ z = \frac{K}{a} - W_0\!\left(\frac{b-a}{a}\, e^{K/a}\right), \]

with the symmetric right-tail form (z -> +infinity) obtained by the a <-> b, z -> -z, K -> -K swap. Both are exact only in the limit; away from the tails the correct branch and argument regime are case-dependent, so probcal ships the seed-plus-Halley construction above rather than a Lambert-W evaluator. It is recorded here as the closed form the numerics are approximating, not as an implementation.

Degenerate exponents fall back to their own closed forms: a == b uses the affine formula directly; a == 0 (h ranges over (0, infinity), so p above sigma(c) only) and b == 0 (range (-infinity, 0), p below sigma(c) only) invert via expm1/log, and raise UnattainableTargetError outside the attainable range; a == b == 0 (a constant map) raises NotImplementedError, since no point has a well-defined preimage.

References¶

Rawal, K., Kamar, E., Lakkaraju, H. (2020). "Algorithmic Recourse in the Wild: Understanding the Impact of Data and Model Shifts." arXiv:2012.11788.
Upadhyay, S., Joshi, S., Lakkaraju, H. (2021). "Towards Robust and Reliable Algorithmic Recourse." Advances in Neural Information Processing Systems 34 (NeurIPS 2021).