diff --git a/decisiongrounding/schema/paired_stats.schema.json b/decisiongrounding/schema/paired_stats.schema.json index 7ed3ed6..6056002 100644 --- a/decisiongrounding/schema/paired_stats.schema.json +++ b/decisiongrounding/schema/paired_stats.schema.json @@ -54,6 +54,8 @@ "n_discordant": { "type": "integer", "minimum": 0 }, "statistic": { "type": "integer", "minimum": 0 }, "p_value": { "type": "number", "minimum": 0, "maximum": 1 }, + "p_value_holm": { "type": "number", "minimum": 0, "maximum": 1 }, + "family": { "enum": ["confirmatory", "secondary", "exploratory"] }, "method": { "const": "mcnemar-exact-binomial" }, "degenerate": { "const": true } } diff --git a/decisiongrounding/scoring/crossover.py b/decisiongrounding/scoring/crossover.py index 3b7215e..3c40dd6 100644 --- a/decisiongrounding/scoring/crossover.py +++ b/decisiongrounding/scoring/crossover.py @@ -29,7 +29,7 @@ from scenarios.loader import Scenario from scoring.metrics import summarize from scoring.scorer import score -from scoring.stats import stats_by_n +from scoring.stats import annotate_holm_family, stats_by_n from util.io import atomic_write_text _TOKEN = re.compile(r"[a-z0-9]+") @@ -176,6 +176,17 @@ def _point(n, adhered, total, retrieved_flags, token_estimates, usages, return point +def finalize_crossover_stats(cells) -> dict | None: + """The crossover `stats` block: per-N paired significance plus the + pre-registered Holm annotation on the H1 family (see + `scoring.stats.annotate_holm_family`).""" + if not cells: + return None + stats = stats_by_n(cells) + annotate_holm_family(stats) + return stats + + def _envelope(discriminating, use_real, pool, seed, ns, model_version, embedder_spec, pool_dir, scenarios_dir, arms, points, per_scenario, errors, cells): """The dataset dict both builders return.""" @@ -208,7 +219,7 @@ def _envelope(discriminating, use_real, pool, seed, ns, model_version, embedder_ # paired analysis (spec/analysis-plan-amendment-1.md) runs on. The # aggregated fractions above are for display; these are for inference. "cells": cells, - "stats": stats_by_n(cells) if cells else None, + "stats": finalize_crossover_stats(cells), } @@ -809,7 +820,7 @@ def _aggregate_seeds(per_seed, arms, ns, pairs) -> dict: # cross-seed statistics is scenario x seed (common random numbers). cells = [c for _, ds in per_seed for c in (ds.get("cells") or [])] base["cells"] = cells - base["stats"] = stats_by_n(cells) if cells else None + base["stats"] = finalize_crossover_stats(cells) paired = _paired(per_seed, ns, pairs) if paired: base["paired"] = paired @@ -885,7 +896,7 @@ def merge_seed_datasets(existing: dict, new_per_seed, arms, ns, pair) -> dict: if old_cells is not None: cells = list(old_cells) + [c for _, ds in add for c in (ds.get("cells") or [])] base["cells"] = cells - base["stats"] = stats_by_n(cells) if cells else None + base["stats"] = finalize_crossover_stats(cells) else: base["cells"] = None base["stats"] = None diff --git a/decisiongrounding/scoring/stats.py b/decisiongrounding/scoring/stats.py index 7de1437..c9846f9 100644 --- a/decisiongrounding/scoring/stats.py +++ b/decisiongrounding/scoring/stats.py @@ -22,6 +22,15 @@ # pins which pair carries the confirmatory weight. CONFIRMATORY_PAIR = ("rac", "naive_rag") +# The pre-registered H1 family (spec/analysis-plan-amendment-1.md): the +# confirmatory test is rac-vs-naive_rag adherence at N=300, reported +# UNcorrected; the secondary cells at N in {50, 150} of the same pair/outcome +# are Holm-corrected; N=10 is descriptive and everything else is exploratory, +# reported raw and labelled as such. +CONFIRMATORY_OUTCOME = "adherent" +CONFIRMATORY_N = 300 +SECONDARY_NS = (50, 150) + # Join fields for pairing records across arms. `example_id` is the alias the # gitchameleon resolution records use for `scenario_id`; `seed` and `N` scope # crossover cells so a pair is always the same scenario under the same @@ -98,12 +107,17 @@ def risk_difference(b: int, c: int, n: int) -> dict: d = (b - c) / n over all n pairs; the paired SE uses only the discordant structure: sqrt(b + c - (b - c)^2 / n) / n. + + A risk difference is bounded [-1, 1], but the Wald interval is a normal + approximation that can spill past that (e.g. b large, c = 0, small n). The + CI is clipped to the parameter space so it never reports an impossible + rate difference. """ if n <= 0: return {"diff": 0.0, "ci": [0.0, 0.0], "degenerate": True} diff = (b - c) / n se = ((b + c - (b - c) ** 2 / n) ** 0.5) / n - return {"diff": diff, "ci": [diff - Z95 * se, diff + Z95 * se]} + return {"diff": diff, "ci": [max(-1.0, diff - Z95 * se), min(1.0, diff + Z95 * se)]} def paired_odds_ratio(b: int, c: int) -> dict: @@ -223,6 +237,55 @@ def stats_by_n( } +def holm_adjust(pvalues: list[float]) -> list[float]: + """Holm-Bonferroni step-down adjusted p-values, returned in input order. + + Controls the family-wise error rate across `m = len(pvalues)` tests with no + independence assumption. The k-th smallest raw p (k = 1..m) is scaled by + (m - k + 1); a running maximum keeps the adjusted sequence monotonic and + each value is capped at 1. `m <= 1` returns the inputs unchanged (nothing + to correct). Pure and deterministic. + """ + m = len(pvalues) + if m <= 1: + return [min(1.0, p) for p in pvalues] + order = sorted(range(m), key=lambda i: pvalues[i]) + adj = [0.0] * m + running = 0.0 + for rank, idx in enumerate(order): + running = max(running, min(1.0, (m - rank) * pvalues[idx])) + adj[idx] = running + return adj + + +def annotate_holm_family(stats_by_n_result: dict) -> None: + """Tag the pre-registered H1 family in a crossover `stats` block, in place. + + For the confirmatory pair/outcome (rac vs naive_rag, adherent), the cells + at `SECONDARY_NS` get a Holm-adjusted `p_value_holm` (Holm across whichever + of {50, 150} are present) and `family = "secondary"`; the `CONFIRMATORY_N` + cell is tagged `family = "confirmatory"` and left UNcorrected — the single + pre-registered primary test. Every other cell is untouched (exploratory, + reported raw and labelled). See spec/analysis-plan-amendment-1.md. + """ + key = f"{CONFIRMATORY_PAIR[0]}_vs_{CONFIRMATORY_PAIR[1]}" + + def _mcnemar_at(n: int) -> dict | None: + block = stats_by_n_result.get(n) + if not block or block.get("outcome") != CONFIRMATORY_OUTCOME: + return None + return (block.get("pairs") or {}).get(key, {}).get("mcnemar") + + secondary = [mc for mc in (_mcnemar_at(n) for n in SECONDARY_NS) if mc is not None] + if secondary: + for mc, adj in zip(secondary, holm_adjust([mc["p_value"] for mc in secondary])): + mc["p_value_holm"] = adj + mc["family"] = "secondary" + conf = _mcnemar_at(CONFIRMATORY_N) + if conf is not None: + conf["family"] = "confirmatory" + + def records_from_runs(runs: list[dict], *, outcome: str = "adherent") -> list[dict]: """Flatten run-report entries (runner/cli.py `runs`) into stats records.""" return [ diff --git a/decisiongrounding/scripts/paper_figs.py b/decisiongrounding/scripts/paper_figs.py index 8df0c4d..467b97a 100644 --- a/decisiongrounding/scripts/paper_figs.py +++ b/decisiongrounding/scripts/paper_figs.py @@ -59,6 +59,15 @@ def _fmt_p_tex(p: float) -> str: return "$<$0.0001" if p < 1e-4 else f"{p:.4f}" +def _holm_tex(mc: dict) -> str: + """The Holm-p / family cell for the paper table.""" + if "p_value_holm" in mc: + return _fmt_p_tex(mc["p_value_holm"]) + if mc.get("family") == "confirmatory": + return "conf.\\ (uncorr.)" + return "--" + + def write_stats_table(ds: dict, outdir: Path) -> Path | None: """Emit `stats_table.tex` — a booktabs tabular of the pre-registered paired analysis per N, \\input by the paper so the numbers cannot drift from the @@ -69,9 +78,9 @@ def write_stats_table(ds: dict, outdir: Path) -> Path | None: rows = [ "% Auto-generated by scripts/paper_figs.py from the crossover dataset.", "% Do not edit: regenerate with `make paper-figs`.", - r"\begin{tabular}{llrrrrr}", + r"\begin{tabular}{llrrrrrr}", r"\toprule", - r"$N$ & pair & $b$ & $c$ & McNemar $p$ & risk diff [95\% CI] & OR [95\% CI] \\", + r"$N$ & pair & $b$ & $c$ & McNemar $p$ & Holm $p$ / family & risk diff [95\% CI] & OR [95\% CI] \\", r"\midrule", ] for n in sorted(stats, key=lambda k: int(k)): @@ -82,10 +91,16 @@ def write_stats_table(ds: dict, outdir: Path) -> Path | None: else f"{orr['or']:.2f} [{orr['ci'][0]:.2f}, {orr['ci'][1]:.2f}]") pair_tex = name.replace("_vs_", " vs ").replace("_", r"\_") rows.append( - f"{n} & {pair_tex} & {mc['b']} & {mc['c']} & {p_txt} & " + f"{n} & {pair_tex} & {mc['b']} & {mc['c']} & {p_txt} & {_holm_tex(mc)} & " f"{rd['diff']:+.2f} [{rd['ci'][0]:+.2f}, {rd['ci'][1]:+.2f}] & {or_txt} \\\\" ) - rows += [r"\bottomrule", r"\end{tabular}"] + rows += [ + r"\bottomrule", + r"\end{tabular}", + r"% Holm column: rac-vs-naive\_rag adherence secondary cells N$\in\{50,150\}$ are", + r"% Holm-corrected (pre-registered); N=300 is the uncorrected confirmatory test;", + r"% all other cells are exploratory and uncorrected (--).", + ] out = outdir / "stats_table.tex" out.write_text("\n".join(rows) + "\n", encoding="utf-8") return out diff --git a/decisiongrounding/scripts/report.py b/decisiongrounding/scripts/report.py index 3796a71..66de176 100644 --- a/decisiongrounding/scripts/report.py +++ b/decisiongrounding/scripts/report.py @@ -31,7 +31,7 @@ from scoring.cost import cost_by_arm, dollars # noqa: E402 from scoring.health import scenario_health # noqa: E402 from scoring.snr import signal_to_noise # noqa: E402 -from scoring.stats import paired_significance # noqa: E402 +from scoring.stats import SECONDARY_NS, paired_significance # noqa: E402 _ARM_DESC = { "context_dump": "pastes the entire corpus into the prompt (the no-retrieval ceiling)", @@ -235,11 +235,20 @@ def _fmt_p(p: float) -> str: return "<0.0001" if p < 1e-4 else f"{p:.4f}" +def _holm_cell(mc: dict) -> str: + """The family / Holm-adjusted-p cell for one McNemar result.""" + if "p_value_holm" in mc: + return f"{_fmt_p(mc['p_value_holm'])} (Holm)" + if mc.get("family") == "confirmatory": + return "confirmatory (uncorr.)" + return "—" + + def _stats_pair_rows(pairs_block: dict) -> list[str]: """Markdown rows for one paired-significance block (per arm pair).""" rows = [ - "| pair | n | b | c | McNemar p | risk diff [95% CI] | odds ratio [95% CI] |", - "|---|--:|--:|--:|--:|--:|--:|", + "| pair | n | b | c | McNemar p | Holm p / family | risk diff [95% CI] | odds ratio [95% CI] |", + "|---|--:|--:|--:|--:|--:|--:|--:|", ] for name, st in pairs_block.items(): mc, rd, orr = st["mcnemar"], st["risk_difference"], st["odds_ratio"] @@ -248,6 +257,7 @@ def _stats_pair_rows(pairs_block: dict) -> list[str]: else f"{orr['or']:.2f} [{orr['ci'][0]:.2f}, {orr['ci'][1]:.2f}]") rows.append( f"| `{name}` | {st['n_pairs']} | {mc['b']} | {mc['c']} | {p_txt} | " + f"{_holm_cell(mc)} | " f"{rd['diff']:+.2f} [{rd['ci'][0]:+.2f}, {rd['ci'][1]:+.2f}] | {or_txt} |" ) return rows @@ -277,6 +287,12 @@ def _stats_section(run: dict, dataset: dict | None) -> list[str]: "`b` = first arm adherent where the second is not; `c` = the reverse. " "Analysis only — these statistics never gate anything (ADR-066/ADR-097).", "", + "**Multiple comparisons.** The confirmatory test is `rac`-vs-`naive_rag` " + "adherence at **N=300**, reported uncorrected. Its secondary cells at " + "**N∈{50,150}** are **Holm-corrected** (the `Holm p / family` column). " + "Every other pair and N is exploratory context, reported raw and " + "uncorrected (`—`).", + "", ] return intro + parts @@ -547,9 +563,24 @@ def delta(s): return "\n".join(lines) +def _secondary_holm_line(dataset: dict) -> str | None: + """One-line summary of the Holm-corrected H1 secondary cells (N∈{50,150}).""" + ds_stats = dataset.get("stats") or {} + bits = [] + for n in SECONDARY_NS: + block = ds_stats.get(n) or ds_stats.get(str(n)) or {} + mc = (block.get("pairs") or {}).get("rac_vs_naive_rag", {}).get("mcnemar") + if mc and "p_value_holm" in mc: + bits.append(f"N={n}: Holm p={_fmt_p(mc['p_value_holm'])}") + if not bits: + return None + return "Secondary (Holm-corrected across N∈{50,150}): " + "; ".join(bits) + "." + + def _mcnemar_line(dataset: dict, top) -> str | None: """The pre-registered confirmatory sentence for rac vs naive_rag at the top - N, when the dataset carries per-cell stats.""" + N, when the dataset carries per-cell stats — followed by the Holm-corrected + secondary summary when present.""" ds_stats = dataset.get("stats") or {} block = ds_stats.get(top) or ds_stats.get(str(top)) or {} st = (block.get("pairs") or {}).get("rac_vs_naive_rag") @@ -557,13 +588,16 @@ def _mcnemar_line(dataset: dict, top) -> str | None: return None mc = st["mcnemar"] if mc.get("degenerate"): - return (f"Pre-registered confirmatory test at N={top}: exact McNemar is " - "degenerate (the arms never disagreed on any paired cell) — " - "no discordant evidence either way.") - return (f"Pre-registered confirmatory test at N={top}: exact McNemar " - f"b={mc['b']}, c={mc['c']}, p={_fmt_p(mc['p_value'])}; paired risk " - f"difference {st['risk_difference']['diff']:+.2f} " - f"[{st['risk_difference']['ci'][0]:+.2f}, {st['risk_difference']['ci'][1]:+.2f}].") + line = (f"Pre-registered confirmatory test at N={top} (uncorrected): exact " + "McNemar is degenerate (the arms never disagreed on any paired " + "cell) — no discordant evidence either way.") + else: + line = (f"Pre-registered confirmatory test at N={top} (uncorrected): exact " + f"McNemar b={mc['b']}, c={mc['c']}, p={_fmt_p(mc['p_value'])}; paired " + f"risk difference {st['risk_difference']['diff']:+.2f} " + f"[{st['risk_difference']['ci'][0]:+.2f}, {st['risk_difference']['ci'][1]:+.2f}].") + secondary = _secondary_holm_line(dataset) + return f"{line}\n\n{secondary}" if secondary else line _PROSE_PIPELINE = """\ diff --git a/decisiongrounding/tests/test_report_holm.py b/decisiongrounding/tests/test_report_holm.py new file mode 100644 index 0000000..85e496d --- /dev/null +++ b/decisiongrounding/tests/test_report_holm.py @@ -0,0 +1,48 @@ +"""Holm correction + family labels render in the report and the paper table.""" + +from pathlib import Path + +from scoring.stats import annotate_holm_family, stats_by_n +from scripts.paper_figs import write_stats_table +from scripts.report import _mcnemar_line, _stats_pair_rows, _stats_section + + +def _annotated_stats(): + cells = [] + for n in (50, 150, 300): + for i in range(6): + cells.append({"seed": 0, "N": n, "arm": "rac", "scenario_id": f"s{i}", "adherent": True}) + cells.append({"seed": 0, "N": n, "arm": "naive_rag", "scenario_id": f"s{i}", "adherent": i < 2}) + stats = stats_by_n(cells) + annotate_holm_family(stats) + return stats + + +def test_pair_rows_show_holm_and_family(): + stats = _annotated_stats() + sec = "\n".join(_stats_pair_rows(stats[50]["pairs"])) + assert "Holm p / family" in sec + assert "(Holm)" in sec # secondary cell corrected + conf = "\n".join(_stats_pair_rows(stats[300]["pairs"])) + assert "confirmatory (uncorr.)" in conf # N=300 uncorrected, tagged + + +def test_stats_section_has_multiple_comparisons_note(): + ds = {"stats": _annotated_stats()} + out = "\n".join(_stats_section({}, ds)) + assert "Multiple comparisons" in out and "Holm-corrected" in out + assert "N=300" in out + + +def test_mcnemar_line_is_uncorrected_with_secondary_holm_summary(): + ds = {"stats": _annotated_stats()} + line = _mcnemar_line(ds, 300) + assert "uncorrected" in line + assert "Holm-corrected across N" in line # the secondary companion line + + +def test_paper_table_has_holm_column(tmp_path): + out = write_stats_table({"stats": _annotated_stats()}, tmp_path) + tex = Path(out).read_text() + assert "Holm $p$ / family" in tex + assert "conf.\\ (uncorr.)" in tex # the N=300 confirmatory cell diff --git a/decisiongrounding/tests/test_stats.py b/decisiongrounding/tests/test_stats.py index eefcf18..3a1f346 100644 --- a/decisiongrounding/tests/test_stats.py +++ b/decisiongrounding/tests/test_stats.py @@ -8,8 +8,10 @@ import pytest from scoring.stats import ( + annotate_holm_family, binom_two_sided_p, collect_pairs, + holm_adjust, mcnemar_exact, pair_stats, paired_odds_ratio, @@ -80,6 +82,63 @@ def test_risk_difference_known_value(): assert risk_difference(0, 0, 0)["degenerate"] is True +def test_risk_difference_ci_is_clamped_to_parameter_space(): + # b=9,c=0,n=10: diff=0.9, Wald upper = 0.9 + 1.96*sqrt(0.9)/10 ~= 1.086 -> clip to 1. + hi = risk_difference(9, 0, 10) + assert hi["diff"] == pytest.approx(0.9) + assert hi["ci"][1] == 1.0 and hi["ci"][0] >= -1.0 + # symmetric: lower bound clips to -1. + lo = risk_difference(0, 9, 10) + assert lo["ci"][0] == -1.0 and lo["ci"][1] <= 1.0 + + +# --- Holm correction ---------------------------------------------------------- + +def test_holm_adjust_known_vectors(): + assert holm_adjust([0.01, 0.04]) == pytest.approx([0.02, 0.04]) + assert holm_adjust([]) == [] + assert holm_adjust([0.03]) == pytest.approx([0.03]) # m=1 identity + assert holm_adjust([0.6, 0.7]) == pytest.approx([1.0, 1.0]) # capped at 1 + # order preserved; monotonic when sorted (the larger raw p can't adjust below + # the smaller's adjusted value). + adj = holm_adjust([0.04, 0.01]) + assert adj == pytest.approx([0.04, 0.02]) + + +def _stats_block(pvals_by_n, pair="rac_vs_naive_rag", outcome="adherent"): + return {n: {"outcome": outcome, + "pairs": {pair: {"mcnemar": {"p_value": p}}}} + for n, p in pvals_by_n.items()} + + +def test_annotate_holm_family_corrects_secondary_leaves_confirmatory(): + stats = _stats_block({10: 0.9, 50: 0.01, 150: 0.04, 300: 0.001}) + annotate_holm_family(stats) + mc = lambda n: stats[n]["pairs"]["rac_vs_naive_rag"]["mcnemar"] + # secondary {50,150} Holm-corrected across the family of 2 + assert mc(50)["p_value_holm"] == pytest.approx(0.02) and mc(50)["family"] == "secondary" + assert mc(150)["p_value_holm"] == pytest.approx(0.04) and mc(150)["family"] == "secondary" + # confirmatory N=300 tagged but UNcorrected (no holm key) + assert mc(300)["family"] == "confirmatory" and "p_value_holm" not in mc(300) + # N=10 (descriptive) untouched + assert "family" not in mc(10) and "p_value_holm" not in mc(10) + + +def test_annotate_holm_family_ignores_non_confirmatory_pairs(): + stats = _stats_block({50: 0.01, 150: 0.04}, pair="context_dump_vs_no_grounding") + annotate_holm_family(stats) + mc = stats[50]["pairs"]["context_dump_vs_no_grounding"]["mcnemar"] + assert "family" not in mc and "p_value_holm" not in mc + + +def test_annotate_holm_family_partial_grid_is_identity(): + # only N=50 present -> family of 1 -> Holm is the identity + stats = _stats_block({50: 0.03}) + annotate_holm_family(stats) + mc = stats[50]["pairs"]["rac_vs_naive_rag"]["mcnemar"] + assert mc["p_value_holm"] == pytest.approx(0.03) and mc["family"] == "secondary" + + def test_paired_odds_ratio_zero_cell_is_degenerate_not_fudged(): orr = paired_odds_ratio(5, 0) assert orr["or"] is None and orr["ci"] is None and orr["degenerate"] is True @@ -182,11 +241,25 @@ def test_paired_significance_validates_against_schema(): import jsonschema # type: ignore except ImportError: pytest.skip("jsonschema not installed") + schema = json.loads(_STATS_SCHEMA.read_text()) out = paired_significance(_designed_records()) - jsonschema.validate(out, json.loads(_STATS_SCHEMA.read_text())) + jsonschema.validate(out, schema) report = json.loads(_HEADLINE.read_text()) real = paired_significance(records_from_runs(report["runs"])) - jsonschema.validate(real, json.loads(_STATS_SCHEMA.read_text())) + jsonschema.validate(real, schema) + # a Holm-annotated crossover block (p_value_holm + family) still validates + cells = [] + for n in (50, 150, 300): + for i in range(6): + cells.append({"seed": 0, "N": n, "arm": "rac", + "scenario_id": f"s{i}", "adherent": True}) + cells.append({"seed": 0, "N": n, "arm": "naive_rag", + "scenario_id": f"s{i}", "adherent": i < 2}) + annotated = stats_by_n(cells) + annotate_holm_family(annotated) + assert "p_value_holm" in annotated[50]["pairs"]["rac_vs_naive_rag"]["mcnemar"] + for block in annotated.values(): + jsonschema.validate(block, schema) # --- crossover integration ----------------------------------------------------