Prediction reports: training-set membership, novelty and conformal intervals (v4.4.0)

CHANGELOG.md

@@ -6,6 +6,49 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to
 [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [4.4.0] - 2026-08-31
+
+### Added
+
+- `prediction_report`: predictions with provenance and uncertainty per PSM. Returns a
+  DataFrame with, next to `predicted_rt`: a conformal prediction interval (`ci_lower`,
+  `ci_upper`) at a chosen coverage, exact-match membership against the calibration reference
+  (`in_reference`) and the Levenshtein distance to the closest reference sequence
+  (`dist_to_reference`); with a training index also membership in the corpus the bundled
+  multitask model was trained on (`in_training`), membership within the training sets of the
+  setups the calibration selected (`in_selected_heads_training`) and the distance to the
+  closest training sequence (`dist_to_training`, exact up to 10 and capped beyond).
+
+  The interval is cross-fitted split-conformal on the reference: the reference is split into
+  folds, each fold is predicted by a calibration fitted on the other folds, and the half-width
+  is a finite-sample quantile of those honest residuals per predicted-RT bin. On eight PRIDE
+  setups no DeepLC model was trained on, the empirical coverage of the 90 % interval was 0.88
+  to 0.97 per setup (median 0.91), with widths from 4 % of the gradient on well-behaved setups
+  to 79 % on a run that pools several fractions. Coverage is marginal, not per-peptide.
+
+  With a multi-head calibration the width is also **per peptide** (`per_peptide_width`, on by
+  default): the residuals are divided by how far the combined setup heads lie apart for that
+  peptide before the quantile is taken, and multiplied by it again at prediction time. Two
+  peptides predicted at the same retention time therefore no longer share one interval. On the
+  six held-out setups this raised the worst conditional slice from 0.851 to 0.882 and the
+  Spearman correlation between width and error from 0.15 to 0.25, for 10 % wider intervals;
+  the largest gains are on the setups where the RT-only width was weakest. Set
+  `per_peptide_width=False` for widths that depend on the predicted retention time alone.
+
+- `Calibration.disagreement`, the per-input uncertainty a calibration can report, implemented
+  by `MultiHeadRidgeCalibration` as the ridge-weighted spread of its calibrated head
+  estimates and returning None elsewhere.
+
+- `TrainingIndex`: an index of the multitask training corpus (10,105,640 canonical
+  peptidoform keys, their 65,139,832 setup observations, 6,157,558 unique stripped sequences).
+  Distributed separately from the package as a single 105 MB `.dlcidx` file: an LZMA zip
+  holding 40-bit key hashes in a bucketed layout (false positive about once per 100,000
+  membership queries, irrelevant for a provenance flag), per-key setup lists and the unique
+  sequences. A raw memory-mapped directory form with exact 64-bit hashes is read as well.
+  `prediction_report` takes either as an optional argument and works without one.
+
+- Dependency: `rapidfuzz` (Levenshtein distances).
+
 ## [4.3.0] - 2026-09-02
 
 ### Changed

deeplc/__init__.py

@@ -11,9 +11,12 @@
     save_model,
     train,
 )
+from deeplc.report import TrainingIndex, prediction_report
 
 __version__: str = version("deeplc")
 __all__: list[str] = [
+    "TrainingIndex",
+    "prediction_report",
     "calibrate",
     "predict",
     "predict_and_calibrate",

deeplc/calibration.py

@@ -46,6 +46,16 @@ def transform(self, source: np.ndarray) -> np.ndarray:
         """Transform source values into the calibrated target space."""
         ...
 
+    def disagreement(self, source: np.ndarray) -> np.ndarray | None:  # noqa: ARG002
+        """
+        Per-input uncertainty score, or None when the calibration has none.
+
+        A calibration that combines several estimates of the same retention time can report
+        how far they lie apart for each input, which :func:`deeplc.report.prediction_report`
+        uses to scale its prediction intervals per peptide.
+        """
+        return None
+
 
 class IdentityCalibration(Calibration):
     """No calibration; returns inputs unchanged."""
@@ -444,13 +454,40 @@ def transform(self, source: np.ndarray) -> np.ndarray:
             )
         if source.shape[0] == 0:
             return np.array([])
-        calibrated = np.column_stack(
+        return np.asarray(self._ridge.predict(self._calibrated_columns(source)), dtype=np.float64)
+
+    def _calibrated_columns(self, source: np.ndarray) -> np.ndarray:
+        """Each selected head's own estimate of the retention time, in the reference's unit."""
+        head_idx = cast(np.ndarray, self._head_idx)
+        return np.column_stack(
             [
                 np.asarray(cal.transform(source[:, head].astype(np.float32)), dtype=np.float64)
                 for cal, head in zip(self._head_calibrations, head_idx, strict=True)
             ]
         )
-        return np.asarray(self._ridge.predict(calibrated), dtype=np.float64)
+
+    def disagreement(self, source: np.ndarray) -> np.ndarray | None:
+        """
+        How far the combined setup heads lie apart for each input, in the reference's unit.
+
+        Every selected head estimates the retention time of the same peptide, so the spread of
+        those estimates, weighted by the ridge weight each head received, is an uncertainty
+        that varies from peptide to peptide rather than only along the gradient. Returns None
+        while the calibration is unfitted or combines a single head, which carries no spread.
+        """
+        if not self.is_fitted:
+            return None
+        columns = np.asarray(source, dtype=np.float64)
+        if columns.ndim == 1:
+            columns = columns[:, None]
+        if columns.shape[0] == 0 or len(cast(np.ndarray, self._head_idx)) < 2:
+            return None
+        weights = np.abs(np.asarray(self._ridge.coef_, dtype=np.float64).ravel())
+        total = weights.sum()
+        weights = weights / total if total > 0 else np.full(len(weights), 1 / len(weights))
+        estimates = self._calibrated_columns(columns)
+        mean = estimates @ weights
+        return np.sqrt(((estimates - mean[:, None]) ** 2) @ weights)
 
 
 def _rank_heads_by_correlation(source: np.ndarray, target: np.ndarray) -> np.ndarray: