diff --git a/cmd_chat/sor/analysis/confirm.py b/cmd_chat/sor/analysis/confirm.py new file mode 100644 index 0000000..1e9c4e8 --- /dev/null +++ b/cmd_chat/sor/analysis/confirm.py @@ -0,0 +1,204 @@ +"""RQ1/RQ2 confirmatory decision layer — the frozen prereg §6 gates, in code. + +This turns the §6 analysis plan into mechanical, pre-registered decisions for the +four confirmatory tests the **lead paper (G4 + RQ1 + RQ2)** reports: + + * **RQ1-P1 (leak):** bridge-on correlation AUC + BCa 95% CI. Gate = the CI + **excludes 0.5** (D2). Materiality is a *separate label*, not the gate: + CI lower bound ≥ 0.60 ⇒ "material", 0.5 < lo < 0.60 ⇒ "weak-but-real". + * **RQ1-P2 (padding efficacy):** paired ΔAUC = AUC(no-pad) − AUC(pad) over + circuits; effective iff the CI **> 0**. + * **RQ2-P1 (anonymity set, two-sided):** ΔH = H(federated) − H(single, matched + N) with Miller–Madow per-circuit entropy; **grow** if CI > 0, **honest-shrink** + if CI < 0, **inconclusive** if it spans 0. The sign is *not* presumed. + * **RQ2-P3 (mechanism):** Spearman ρ between top-k=3 bridge concentration and + per-circuit H; negative ρ quantifies funnelling. + +Multiplicity: :func:`apply_holm` corrects the reported tests against the **full +frozen family of 7** (§6), not the reported subset — see the module's +``FROZEN_FAMILY`` and the stage-05 Holm clarification. The bootstrap p-values are +used only to *order* the Holm step-down; every reported decision is a CI gate, +never a bare p (§6, rigor-standards §Statistics). + +Pure analysis: no I/O beyond what a caller serialises, no engine, no traffic. All +detectors are the §5-calibrated instruments (`detectors.py`), never re-fit here. +""" + +from __future__ import annotations + +from dataclasses import dataclass +from typing import Dict, List, Sequence, Tuple + +from cmd_chat.sor.analysis.detectors import auc +from cmd_chat.sor.analysis.stats import ( + CIResult, + HolmResult, + bootstrap_ci, + holm_bonferroni, + mean, + miller_madow_entropy_bits, + spearman, + two_sample_diff_ci, + two_sided_bootstrap_p, +) + +# The pre-registered confirmatory family (frozen prereg §6). The lead paper +# reports the first four; RQ3's three are the severable follow-on (D6/§8). The +# family SIZE stays 7 for Holm regardless of how many are reported here. +FROZEN_FAMILY: Tuple[str, ...] = ( + "RQ1-P1", "RQ1-P2", "RQ2-P1", "RQ2-P3", + "RQ3-P1-perf", "RQ3-P1-latency", "RQ3-P2", +) +FROZEN_FAMILY_SIZE = len(FROZEN_FAMILY) # 7 +LEAD_PAPER_TESTS: Tuple[str, ...] = ("RQ1-P1", "RQ1-P2", "RQ2-P1", "RQ2-P3") + +MATERIALITY_FLOOR = 0.60 # §6 [APPROVAL] — a *label*, not the RQ1 gate. +TOP_K = 3 # RQ2-P3 [APPROVAL]. + +# A circuit-pair scored by the correlator: (score, linked?) — linked = same +# circuit (diagonal), unlinked = different (off-diagonal). +Pair = Tuple[float, bool] + + +@dataclass(frozen=True) +class ConfirmTest: + """One confirmatory test's frozen decision: the effect size, its CI, the + pre-registered gate outcome, any secondary label, and the bootstrap p used + only for Holm ordering.""" + + name: str + effect: str # human name of the effect size (e.g. "AUC", "ΔH") + ci: CIResult + decision: str # e.g. "leak", "null", "grow", "shrink", "inconclusive" + label: str # secondary label (materiality / direction); "" if none + p_for_holm: float + + def as_dict(self) -> Dict: + return { + "name": self.name, + "effect": self.effect, + "decision": self.decision, + "label": self.label, + "p_for_holm": self.p_for_holm, + **self.ci.as_dict(), + } + + +def _auc_of_pairs(pairs: Sequence[Pair]) -> float: + pos = [s for s, linked in pairs if linked] + neg = [s for s, linked in pairs if not linked] + return auc(pos, neg) + + +def rq1_p1_leak( + pairs: Sequence[Pair], *, seed: int = 0, n_resamples: int = 10_000, alpha: float = 0.05 +) -> ConfirmTest: + """RQ1-P1: bridge-on correlation AUC with BCa CI over circuit pairs. Gate = CI + excludes 0.5 (a leak is present). Materiality label from the CI lower bound.""" + ci, dist = bootstrap_ci(list(pairs), _auc_of_pairs, n_resamples=n_resamples, + alpha=alpha, seed=seed, return_dist=True) + if ci.excludes(0.5) and ci.strictly_greater(0.5): + decision = "leak" + label = "material" if ci.lo >= MATERIALITY_FLOOR else "weak-but-real" + elif ci.strictly_less(0.5): + decision = "anomaly-below-chance" # AUC < 0.5 (should not happen for a real leak) + label = "" + else: + decision = "null" # CI spans 0.5 — no measurable leak + label = "" + return ConfirmTest("RQ1-P1", "AUC", ci, decision, label, + two_sided_bootstrap_p(dist, 0.5)) + + +@dataclass(frozen=True) +class PairedCircuit: + """One circuit contributing correlator pairs under both conditions — the unit + of the RQ1-P2 *paired* bootstrap.""" + + nopad_pairs: Tuple[Pair, ...] + pad_pairs: Tuple[Pair, ...] + + +def _delta_auc(units: Sequence[PairedCircuit]) -> float: + nopad = [p for u in units for p in u.nopad_pairs] + pad = [p for u in units for p in u.pad_pairs] + return _auc_of_pairs(nopad) - _auc_of_pairs(pad) + + +def rq1_p2_padding( + circuits: Sequence[PairedCircuit], *, seed: int = 0, + n_resamples: int = 10_000, alpha: float = 0.05, +) -> ConfirmTest: + """RQ1-P2: paired ΔAUC = AUC(no-pad) − AUC(pad), resampling circuits (the + paired unit). Padding effective iff the CI > 0.""" + ci, dist = bootstrap_ci(list(circuits), _delta_auc, n_resamples=n_resamples, + alpha=alpha, seed=seed, return_dist=True) + decision = "padding-effective" if ci.strictly_greater(0.0) else "padding-ineffective" + return ConfirmTest("RQ1-P2", "ΔAUC", ci, decision, "", + two_sided_bootstrap_p(dist, 0.0)) + + +def _mean_mm_entropy(circuits: Sequence[Sequence[float]]) -> float: + """Mean per-circuit Miller–Madow entropy (bits) over an arm's circuits.""" + return mean([miller_madow_entropy_bits(c) for c in circuits]) + + +def rq2_p1_delta_h( + federated: Sequence[Sequence[float]], + single_house: Sequence[Sequence[float]], + *, seed: int = 0, n_resamples: int = 10_000, alpha: float = 0.05, +) -> ConfirmTest: + """RQ2-P1 (two-sided): ΔH = mean H(federated) − mean H(single-house, matched + N), Miller–Madow per circuit, BCa CI over circuits. grow / honest-shrink / + inconclusive by the sign of the CI — the design does not presume the sign. + + ``federated`` / ``single_house`` are lists of per-circuit sender-posterior + count vectors. The caller is responsible for the §6 matched-N sizing (the + single-house arm's node count equals the federated arm's total consenting + nodes); this function asserts nothing about N — it reports ΔH honestly.""" + ci, dist = two_sample_diff_ci(list(federated), list(single_house), _mean_mm_entropy, + n_resamples=n_resamples, alpha=alpha, seed=seed, + return_dist=True) + if ci.strictly_greater(0.0): + decision, label = "grow", "federation grows the anonymity set" + elif ci.strictly_less(0.0): + decision, label = "shrink", "honest null — federation shrinks (reported with equal prominence)" + else: + decision, label = "inconclusive", "CI spans 0" + return ConfirmTest("RQ2-P1", "ΔH", ci, decision, label, + two_sided_bootstrap_p(dist, 0.0)) + + +def rq2_p3_funnel( + concentration: Sequence[float], per_circuit_h: Sequence[float], + *, seed: int = 0, n_resamples: int = 10_000, alpha: float = 0.05, +) -> ConfirmTest: + """RQ2-P3 (mechanism): Spearman ρ between top-k=3 bridge concentration and + per-circuit entropy H, with BCa CI over circuits. Negative ρ quantifies + funnelling. ``concentration[i]``/``per_circuit_h[i]`` are the two measurements + for circuit i.""" + units = list(zip(concentration, per_circuit_h)) + stat = lambda us: spearman([x for x, _ in us], [y for _, y in us]) + ci, dist = bootstrap_ci(units, stat, n_resamples=n_resamples, alpha=alpha, + seed=seed, return_dist=True) + if ci.strictly_less(0.0): + decision, label = "funnel", "negative ρ — concentration funnels the anonymity set" + elif ci.strictly_greater(0.0): + decision, label = "anti-funnel", "positive ρ" + else: + decision, label = "inconclusive", "CI spans 0" + return ConfirmTest("RQ2-P3", "spearman_rho", ci, decision, label, + two_sided_bootstrap_p(dist, 0.0)) + + +def apply_holm( + tests: Sequence[ConfirmTest], *, alpha: float = 0.05, + family_size: int = FROZEN_FAMILY_SIZE, +) -> List[HolmResult]: + """Holm–Bonferroni over the reported ``tests``, corrected against the full + frozen family (``family_size`` defaults to 7). Reporting a subset of the + pre-registered family never shrinks the correction to the subset — see the + stage-05 Holm clarification. Ordering uses the bootstrap p-values; the + reported decisions remain the CI gates above.""" + return holm_bonferroni({t.name: t.p_for_holm for t in tests}, + alpha=alpha, family_size=family_size) diff --git a/cmd_chat/sor/analysis/stats.py b/cmd_chat/sor/analysis/stats.py index fbb561b..57d5c16 100644 --- a/cmd_chat/sor/analysis/stats.py +++ b/cmd_chat/sor/analysis/stats.py @@ -128,6 +128,19 @@ def _bca_endpoints( return q_lo, q_hi, "bca" +def two_sided_bootstrap_p(thetas: Sequence[float], null: float) -> float: + """A bootstrap two-sided p-value for H0: θ = ``null`` from a bootstrap + distribution ``thetas`` — 2·min(P(θ* ≤ null), P(θ* ≥ null)), capped at 1.0. + Used only to *order* the Holm family; the pre-registered decision gates are + the CIs, never a p-value alone (§6).""" + b = len(thetas) + if b == 0: + return 1.0 + le = sum(1 for t in thetas if t <= null) / b + ge = sum(1 for t in thetas if t >= null) / b + return min(1.0, 2.0 * min(le, ge)) + + def bootstrap_ci( units: Sequence[T], statistic: Callable[[Sequence[T]], float], @@ -136,16 +149,19 @@ def bootstrap_ci( alpha: float = DEFAULT_ALPHA, seed: int = 0, method: str = "bca", -) -> CIResult: + return_dist: bool = False, +): """One-sample bootstrap CI of ``statistic`` over ``units`` (the unit of analysis — a circuit-pair for RQ1, a circuit for RQ2). Resamples ``units`` with replacement ``n_resamples`` times. ``method="bca"`` applies bias-correction + acceleration (falling back to percentile if degenerate); ``"percentile"`` - forces the plain interval.""" + forces the plain interval. With ``return_dist=True`` also returns the sorted + bootstrap distribution (so a p-value can be derived from the same resamples).""" n = len(units) if n == 0: - return CIResult(float("nan"), float("nan"), float("nan"), - alpha, n_resamples, "empty", seed) + res = CIResult(float("nan"), float("nan"), float("nan"), + alpha, n_resamples, "empty", seed) + return (res, []) if return_dist else res theta_hat = float(statistic(units)) rng = random.Random(seed) thetas: List[float] = [] @@ -160,8 +176,9 @@ def bootstrap_ci( else: q_lo, q_hi, used = alpha / 2.0, 1.0 - alpha / 2.0, "percentile" - return CIResult(theta_hat, _percentile(thetas, q_lo), _percentile(thetas, q_hi), - alpha, n_resamples, used, seed) + res = CIResult(theta_hat, _percentile(thetas, q_lo), _percentile(thetas, q_hi), + alpha, n_resamples, used, seed) + return (res, thetas) if return_dist else res def two_sample_diff_ci( @@ -173,15 +190,18 @@ def two_sample_diff_ci( alpha: float = DEFAULT_ALPHA, seed: int = 0, method: str = "bca", -) -> CIResult: + return_dist: bool = False, +): """Bootstrap CI for the difference ``statistic(A) - statistic(B)`` of two independent arms (RQ2-P1: ΔH = H(federated) − H(single-house, matched N)). Each arm is resampled independently. BCa uses a combined leave-one-out - jackknife across both arms (each point dropped from its own arm).""" + jackknife across both arms (each point dropped from its own arm). With + ``return_dist=True`` also returns the sorted bootstrap distribution.""" na, nb = len(units_a), len(units_b) if na == 0 or nb == 0: - return CIResult(float("nan"), float("nan"), float("nan"), - alpha, n_resamples, "empty", seed) + res = CIResult(float("nan"), float("nan"), float("nan"), + alpha, n_resamples, "empty", seed) + return (res, []) if return_dist else res theta_hat = float(statistic(units_a)) - float(statistic(units_b)) rng = random.Random(seed) thetas: List[float] = [] @@ -203,8 +223,9 @@ def two_sample_diff_ci( else: q_lo, q_hi, used = alpha / 2.0, 1.0 - alpha / 2.0, "percentile" - return CIResult(theta_hat, _percentile(thetas, q_lo), _percentile(thetas, q_hi), - alpha, n_resamples, used, seed) + res = CIResult(theta_hat, _percentile(thetas, q_lo), _percentile(thetas, q_hi), + alpha, n_resamples, used, seed) + return (res, thetas) if return_dist else res def mean(xs: Sequence[float]) -> float: diff --git a/tests/test_sor_confirm.py b/tests/test_sor_confirm.py new file mode 100644 index 0000000..088b61c --- /dev/null +++ b/tests/test_sor_confirm.py @@ -0,0 +1,149 @@ +"""RQ1/RQ2 confirmatory gates — validated on synthetic ground truth only. + +Each test constructs a distribution whose §6 verdict is known by construction +(a clear leak / a null; padding that works / doesn't; federation that grows / +shrinks / is flat; a funnelling mechanism) and asserts the frozen gate fires the +right way. No confirmatory-cell data is involved. The Holm test pins that the +reported RQ1/RQ2 subset is corrected against the full frozen family of 7. +""" + +from cmd_chat.sor.analysis.confirm import ( + FROZEN_FAMILY_SIZE, + PairedCircuit, + apply_holm, + rq1_p1_leak, + rq1_p2_padding, + rq2_p1_delta_h, + rq2_p3_funnel, +) + + +class _LCG: + def __init__(self, seed=0x2468ACE0): + self.s = seed + def u(self): + self.s = (1103515245 * self.s + 12345) & 0x7FFFFFFF + return self.s / 0x7FFFFFFF + + +def _pairs(hi_mean, lo_mean, n=80, rng=None): + """n linked pairs near hi_mean + n unlinked near lo_mean, with jitter.""" + rng = rng or _LCG() + out = [] + for _ in range(n): + out.append((hi_mean + 0.1 * rng.u(), True)) + out.append((lo_mean + 0.1 * rng.u(), False)) + return out + + +# --------------------------------------------------------------------------- # +# RQ1-P1 — leak gate + materiality label. +# --------------------------------------------------------------------------- # +def test_rq1p1_material_leak(): + t = rq1_p1_leak(_pairs(0.9, 0.1), seed=1, n_resamples=2000) + assert t.decision == "leak" + assert t.label == "material" # CI lower bound >= 0.60 + assert t.ci.excludes(0.5) + + +def test_rq1p1_null_when_scores_overlap(): + # Linked and unlinked drawn from the same band -> AUC ~ 0.5, CI spans it. + t = rq1_p1_leak(_pairs(0.5, 0.5), seed=1, n_resamples=2000) + assert t.decision == "null" + assert not t.ci.excludes(0.5) + + +# --------------------------------------------------------------------------- # +# RQ1-P2 — padding efficacy (paired ΔAUC). +# --------------------------------------------------------------------------- # +def test_rq1p2_padding_effective_when_pad_lowers_auc(): + rng = _LCG(0x1111) + circuits = [] + for _ in range(30): + nopad = tuple(_pairs(0.9, 0.1, n=4, rng=rng)) # strong linkage + pad = tuple(_pairs(0.5, 0.5, n=4, rng=rng)) # padding blurs it (AUC~0.5) + circuits.append(PairedCircuit(nopad, pad)) + t = rq1_p2_padding(circuits, seed=2, n_resamples=2000) + assert t.decision == "padding-effective" + assert t.ci.strictly_greater(0.0) + + +def test_rq1p2_padding_ineffective_when_no_change(): + rng = _LCG(0x2222) + circuits = [] + for _ in range(30): + nopad = tuple(_pairs(0.7, 0.3, n=4, rng=rng)) + pad = tuple(_pairs(0.7, 0.3, n=4, rng=rng)) # identical regime + circuits.append(PairedCircuit(nopad, pad)) + t = rq1_p2_padding(circuits, seed=2, n_resamples=2000) + assert t.decision == "padding-ineffective" + assert not t.ci.strictly_greater(0.0) + + +# --------------------------------------------------------------------------- # +# RQ2-P1 — ΔH two-sided (grow / shrink / inconclusive). +# --------------------------------------------------------------------------- # +def _uniform_circuits(n_circuits, n_senders, per=50): + return [[per] * n_senders for _ in range(n_circuits)] + + +def _skewed_circuits(n_circuits, n_senders): + # One dominant sender -> low entropy. + return [[1000] + [1] * (n_senders - 1) for _ in range(n_circuits)] + + +def test_rq2p1_grow_when_federation_is_more_uniform(): + fed = _uniform_circuits(30, 8) # high H (~3 bits) + single = _skewed_circuits(30, 8) # low H + t = rq2_p1_delta_h(fed, single, seed=3, n_resamples=2000) + assert t.decision == "grow" + assert t.ci.strictly_greater(0.0) + + +def test_rq2p1_honest_shrink_reported(): + fed = _skewed_circuits(30, 8) # federation funnels -> low H + single = _uniform_circuits(30, 8) # matched-N single house, high H + t = rq2_p1_delta_h(fed, single, seed=3, n_resamples=2000) + assert t.decision == "shrink" + assert t.ci.strictly_less(0.0) + + +def test_rq2p1_inconclusive_when_arms_match(): + fed = _uniform_circuits(30, 8) + single = _uniform_circuits(30, 8) + t = rq2_p1_delta_h(fed, single, seed=3, n_resamples=2000) + assert t.decision == "inconclusive" + assert not t.ci.excludes(0.0) + + +# --------------------------------------------------------------------------- # +# RQ2-P3 — funnelling mechanism (Spearman). +# --------------------------------------------------------------------------- # +def test_rq2p3_funnel_negative_rho(): + # Higher top-3 concentration -> lower per-circuit entropy. + conc = [i / 20.0 for i in range(20)] + h = [3.0 - c for c in conc] + t = rq2_p3_funnel(conc, h, seed=4, n_resamples=2000) + assert t.decision == "funnel" + assert t.ci.strictly_less(0.0) + + +# --------------------------------------------------------------------------- # +# Holm over the full frozen family (size 7) while reporting 4. +# --------------------------------------------------------------------------- # +def test_apply_holm_corrects_against_family_of_seven(): + t1 = rq1_p1_leak(_pairs(0.95, 0.05), seed=1, n_resamples=1500) # tiny p + t2 = rq1_p2_padding( + [PairedCircuit(tuple(_pairs(0.9, 0.1, n=4)), tuple(_pairs(0.5, 0.5, n=4))) + for _ in range(30)], seed=2, n_resamples=1500) + fed, single = _uniform_circuits(30, 8), _skewed_circuits(30, 8) + t3 = rq2_p1_delta_h(fed, single, seed=3, n_resamples=1500) + t4 = rq2_p3_funnel([i / 20.0 for i in range(20)], + [3.0 - i / 20.0 for i in range(20)], seed=4, n_resamples=1500) + holm = apply_holm([t1, t2, t3, t4]) + assert len(holm) == 4 + # Smallest-p test gets the full-family multiplier of 7, not 4. + top = min(holm, key=lambda h: h.rank) + assert top.multiplier == FROZEN_FAMILY_SIZE == 7 + mults = sorted(h.multiplier for h in holm) + assert mults == [4, 5, 6, 7]