feat(#505): q95 and E[y|y>0] as selector frames — and the separation that keeps them out of the contract

docs/ADRs/023_collapsing_posterior_draws.md

@@ -75,10 +75,13 @@ rented GPU — the numpy comes down first, and the conversion is re-runnable fro
 
 ### 2. The method is the one the model declares. It is never hard-coded.
 
-`collapse_origin(..., aggregate_method=...)` implements `arithmetic_mean` and `median` — the same
+`collapse_origin(..., estimator=...)` implements `arithmetic_mean` and `median` — the same
 two names `vhy_021` defines — and **refuses any other string** rather than falling
 back to a default. The converter's default is `arithmetic_mean` because all eight declare it.
 
+*(Amended 2026-09-29: the parameter is now `estimator`, and the tool computes two more of them.
+The set a model may **declare** is still exactly these two. See Amendment 1.)*
+
 **That default and the models' declaration are two places stating one fact, so they are pinned
 together.** `tests/test_roster_conformance.py::test_collapse_declaration_matches_the_converter`
 imports `DEFAULT_AGGREGATE_METHOD` and asserts every roster member matches it, and that every
@@ -207,8 +210,9 @@ python -m tools.collapse.collapse_predictions models/<model> --run-type calibrat
 python -m tools.collapse.plot_collapse_audit <parquet> --draws-dir <origin_i> --out audit.png
 ```
 
-`--aggregate-method` exists and must be passed if a model ever stops declaring `arithmetic_mean`;
-the roster test names it when that happens.
+`--estimator` exists and must be passed if a model ever stops declaring `arithmetic_mean`;
+the roster test names it when that happens. `--aggregate-method` remains accepted as the name it
+had before Amendment 1, so runbooks and saved command lines keep working.
 
 **The plot script is part of the procedure, not a debugging aid.** The tests prove the arithmetic;
 they cannot tell you the field stopped looking like conflict. A person looks at the panels before
@@ -252,6 +256,65 @@ a real run.
 
 ---
 
+## Amendment 1 (2026-09-29) — an estimator this tool can compute is not a method a model may declare
+
+**Status:** Accepted. Extends §2; overrides nothing.
+
+### What changed
+
+`tools/collapse` gained two point estimators beyond the two `vhy_021` defines:
+
+| name | what it is |
+|---|---|
+| `q95` | the 95th percentile across draws, numpy's linear interpolation |
+| `conditional_mean` | `E[y|y>0]` — the mean over the positive draws, equal to `E[y] / P(y>0)` |
+
+The parameter is renamed `aggregate_method` → `estimator`, because for these two the old name
+was false: they are not aggregate methods and `vhy_021` does not define them.
+
+### Why the two sets stay separate
+
+`AGGREGATE_METHODS` is **contract vocabulary** — the only names a model may put in its
+`aggregate_method`, checked against every roster member by
+`tests/test_roster_conformance.py::test_collapse_declaration_matches_the_converter`.
+`ESTIMATORS` is **what this tool can compute**, a strict superset.
+
+Merging them would be the whole defect: a config could declare `q95`, and the roster test — whose
+job is to catch exactly that — would wave it through as a legitimate ADR-021 method. The
+separation is pinned by `test_the_experiment_estimators_are_not_declarable_aggregate_methods`,
+which asserts `AGGREGATE_METHODS` is a *strict* subset and names the difference.
+
+**The standing rule: adding an estimator means adding a key to `ESTIMATORS` and nothing else.
+`AGGREGATE_METHODS` changes only when views-hydranet's ADR-021 changes.**
+
+### What the two extra estimators are for
+
+views-models#505's selector experiment. One posterior cube, collapsed three ways, submitted as
+three apparent models (`_I` mean, `_II` q95, `_III` conditional mean), so the ensemble selector's
+own criteria choose between a mean that under-predicts total fatalities ~5x and two estimators
+that do not. Measured on the 2026-09-28 calibration run, summed over 13 origins and all 7 landed
+models, `q95` and `conditional_mean` run **~4x the mean** — the right order to close that gap.
+
+They are **not** candidates for delivery. No model declares them, `arithmetic_mean` remains the
+default and the delivered estimator, and the mean frames this amendment's code produces are
+**byte-exact against the parquets delivered before it** — checked for all 7 models x 13 origins as
+the acceptance condition for the change.
+
+### What is deliberately not recorded in the parquet
+
+Which estimator produced a frame. The output carries `month_id`, `priogrid_id` and three `pred_*`
+columns and nothing else, by §5, and adding provenance to the frame would change a shape
+`ensemble-updater` already reads. **Provenance lives in the output directory the caller picks, and
+nowhere else** — so a frame moved out of its directory is a frame whose estimator is unknowable.
+That is a real sharp edge and the reason the directories are named `I_mean`, `II_q95` and
+`III_conditional_mean` rather than anything shorter.
+
+### Caveat carried forward
+
+At the roster's `D x K = 16`, `q95` sits between the 15th and 16th order statistics — a coarse
+tail estimate that one draw moves. It is reported, not corrected for. Reading either new estimator
+as a calibrated forecast rather than as a probe of the selector's criteria would be a misuse.
+
 ## References
 
 - `tools/collapse/` — the converter and the audit plots; `tools/collapse/__init__.py` has the usage

tests/test_collapse_predictions.py

@@ -21,6 +21,7 @@
 from tools.collapse.collapse_predictions import (
     AGGREGATE_METHODS,
     DEFAULT_AGGREGATE_METHOD,
+    ESTIMATORS,
     TARGETS,
     CollapseError,
     collapse_origin,
@@ -255,20 +256,30 @@ def test_the_default_method_is_the_one_the_roster_declares():
     assert set(AGGREGATE_METHODS) == {"arithmetic_mean", "median"}
 
 
+def test_the_experiment_estimators_are_not_declarable_aggregate_methods():
+    """The separation is the point. `q95` and `conditional_mean` are estimators this tool can
+    compute; they are NOT names a model may declare. Collapsing the two sets would let a config
+    declare `q95` and have the roster test wave it through as an ADR-021 aggregate method."""
+    assert set(AGGREGATE_METHODS) < set(ESTIMATORS), "aggregate methods must be a strict subset"
+    assert set(ESTIMATORS) - set(AGGREGATE_METHODS) == {"q95", "conditional_mean"}
+    for experiment_only in ("q95", "conditional_mean"):
+        assert experiment_only not in AGGREGATE_METHODS
+
+
 def test_median_is_available_and_is_not_the_mean(origin):
     """views-hydranet ADR-021 allows median. If a model ever declares it we must honour it."""
     raw = np.load(origin / "lr_sb_best" / "y_pred.npy")
-    df = collapse_origin(origin, aggregate_method="median")
+    df = collapse_origin(origin, estimator="median")
     np.testing.assert_allclose(
         df["pred_lr_sb_best"].to_numpy(), np.median(raw.astype("float64"), axis=1), rtol=1e-12
     )
     assert not np.allclose(df["pred_lr_sb_best"], raw.astype("float64").mean(axis=1))
 
 
-def test_mutation_unknown_aggregate_method_is_refused(origin):
+def test_mutation_unknown_estimator_is_refused(origin):
     """A typo must not fall back to the mean and ship an estimator nobody chose."""
-    with pytest.raises(CollapseError, match="unknown aggregate_method"):
-        collapse_origin(origin, aggregate_method="geometric_mean")
+    with pytest.raises(CollapseError, match="unknown estimator"):
+        collapse_origin(origin, estimator="geometric_mean")
 
 
 def test_mutation_log_space_in_ONE_target_only_is_refused(tmp_path):
@@ -319,3 +330,113 @@ def test_mutation_duplicate_identifier_rows_are_refused(origin):
         np.savez(d / "identifiers.npz", time=month, unit=unit)
     with pytest.raises(CollapseError, match="duplicate"):
         collapse_origin(origin)
+
+
+# ── the views-models#505 selector experiment: q95 and E[y|y>0] ────────────────────
+
+
+def test_q95_is_available_and_is_not_the_mean(origin):
+    """Frame II. numpy's linear interpolation, pinned — a different interpolation would move
+    every value in the frame without changing a single row count."""
+    raw = np.load(origin / "lr_sb_best" / "y_pred.npy").astype("float64")
+    df = collapse_origin(origin, estimator="q95")
+    np.testing.assert_allclose(
+        df["pred_lr_sb_best"].to_numpy(), np.quantile(raw, 0.95, axis=1), rtol=1e-12
+    )
+    assert not np.allclose(df["pred_lr_sb_best"], raw.mean(axis=1))
+
+
+def test_conditional_mean_is_the_mean_of_the_positive_draws(origin):
+    """Frame III, against the definition written the obvious way rather than the vectorised
+    way the implementation uses. If the two agree, the `sum / n_positive` shortcut is sound."""
+    raw = np.load(origin / "lr_sb_best" / "y_pred.npy").astype("float64")
+    df = collapse_origin(origin, estimator="conditional_mean")
+    expected = np.array([row[row > 0].mean() if (row > 0).any() else 0.0 for row in raw])
+    np.testing.assert_allclose(df["pred_lr_sb_best"].to_numpy(), expected, rtol=1e-12)
+
+
+def test_conditional_mean_equals_mean_over_probability_of_positive(origin):
+    """The other half of the identity: `E[y|y>0] == E[y] / P(y>0)`. Stated as a test because
+    the docstring claims it, and a claim in a docstring is not evidence."""
+    raw = np.load(origin / "lr_sb_best" / "y_pred.npy").astype("float64")
+    df = collapse_origin(origin, estimator="conditional_mean")
+    p_positive = (raw > 0).mean(axis=1)
+    live = p_positive > 0
+    np.testing.assert_allclose(
+        df["pred_lr_sb_best"].to_numpy()[live],
+        raw.mean(axis=1)[live] / p_positive[live],
+        rtol=1e-12,
+    )
+
+
+def test_conditional_mean_is_never_below_the_arithmetic_mean(origin):
+    """A property, not an example: dividing by `n_positive <= D` cannot shrink the value. It is
+    what lets the scale guard's threshold stay unchanged for the experiment frames — neither
+    new estimator can drop a target's maximum below one the mean already cleared."""
+    for target in TARGETS:
+        mean = collapse_origin(origin, estimator="arithmetic_mean")[f"pred_{target}"].to_numpy()
+        cond = collapse_origin(origin, estimator="conditional_mean")[f"pred_{target}"].to_numpy()
+        assert (cond >= mean - 1e-12).all()
+
+
+def test_conditional_mean_of_an_all_silent_cell_is_zero_not_nan(tmp_path):
+    """The 99.6% case on real output. NaN here would be joined and scored by
+    `ensemble-updater` rather than refused, so the estimator must commit to 0.0."""
+    o = tmp_path / "predictions_calibration_20260101_000000" / "origin_0"
+    draws = np.zeros((ROWS, DRAWS))
+    draws[0, :] = 50.0  # one live cell, so the scale guard has something to clear
+    for target in TARGETS:
+        _write_target(o, target, draws)
+    df = collapse_origin(o, estimator="conditional_mean")
+    silent = df["pred_lr_sb_best"].to_numpy()[1:]
+    assert np.isfinite(silent).all(), "an all-zero cell must not produce NaN or inf"
+    assert (silent == 0.0).all()
+
+
+def test_conditional_mean_equals_the_mean_when_every_draw_is_positive(tmp_path):
+    """The boundary of the inequality above: at `n_positive == D` the two estimators are the
+    same number. A frame III that differed from frame I everywhere would mean the positive mask
+    was wrong, and this is the case that would catch it."""
+    o = tmp_path / "predictions_calibration_20260101_000000" / "origin_0"
+    rng = np.random.default_rng(7)
+    draws = rng.uniform(1.0, 100.0, size=(ROWS, DRAWS))  # no zeros at all
+    for target in TARGETS:
+        _write_target(o, target, draws)
+    mean = collapse_origin(o, estimator="arithmetic_mean")["pred_lr_sb_best"].to_numpy()
+    cond = collapse_origin(o, estimator="conditional_mean")["pred_lr_sb_best"].to_numpy()
+    np.testing.assert_allclose(cond, mean, rtol=1e-6)
+
+
+def test_every_estimator_produces_the_same_shape_and_identifiers(origin):
+    """Whichever estimator ran, the frame the selector reads is the same shape with the same
+    keys. Only `pred_*` may differ — an estimator that dropped or reordered rows would make the
+    three frames non-comparable, which is the whole point of submitting them together."""
+    frames = {name: collapse_origin(origin, estimator=name) for name in ESTIMATORS}
+    reference = frames[DEFAULT_AGGREGATE_METHOD]
+    for name, df in frames.items():
+        assert list(df.columns) == list(reference.columns), name
+        pd.testing.assert_series_equal(df["month_id"], reference["month_id"], check_names=False)
+        pd.testing.assert_series_equal(
+            df["priogrid_id"], reference["priogrid_id"], check_names=False
+        )
+
+
+def test_convert_model_honours_the_estimator_end_to_end(tmp_path, origin):
+    """The CLI path, not just the function: two estimators over the same source must write
+    parquets that differ in value and agree in every identifier."""
+    model = tmp_path / "m"
+    source = origin.parent
+    generated = model / "data" / "generated"
+    generated.mkdir(parents=True)
+    source.rename(generated / source.name)
+
+    out_i = tmp_path / "frame_I"
+    out_ii = tmp_path / "frame_II"
+    written_i = convert_model(model, out_dir=out_i, estimator="arithmetic_mean")
+    written_ii = convert_model(model, out_dir=out_ii, estimator="q95")
+    assert len(written_i) == len(written_ii) == 1
+
+    a, b = pd.read_parquet(written_i[0]), pd.read_parquet(written_ii[0])
+    pd.testing.assert_series_equal(a["month_id"], b["month_id"])
+    pd.testing.assert_series_equal(a["priogrid_id"], b["priogrid_id"])
+    assert not np.allclose(a["pred_lr_sb_best"], b["pred_lr_sb_best"])