#!/usr/bin/env python3 """主分析：按 claims.json 预注册判据做配对检验。用 scipy 实现。 计数式判据一律绑连续量：每个「N/12 个种子上 X 更好」旁边必须同时给出 配对 Wilcoxon 双侧 p、配对差值中位数、rank-biserial 效应量。 """ from __future__ import annotations import json from pathlib import Path import numpy as np from scipy.stats import wilcoxon PKG = Path(__file__).resolve().parents[2] RAW = json.loads((PKG / "work/results/eval_raw.json").read_text(encoding="utf-8")) SEEDS = [str(s) for s in RAW["eval_seeds"]] def series(block, arm, key): return np.array([block[arm][s][key] for s in SEEDS], dtype=float) def paired(better, worse, label, metric): """diff = worse - better；diff>0 表示 better 那一臂的极化方差更低（更好）。""" d = worse - better n_better = int((d > 0).sum()) res = wilcoxon(worse, better, alternative="two-sided", method="exact") r = np.abs(d) order = np.argsort(r, kind="stable") ranks = np.empty(len(d)) ranks[order] = np.arange(1, len(d) + 1) wp = float(ranks[d > 0].sum()) wn = float(ranks[d < 0].sum()) return { "comparison": label, "metric": metric, "n_seeds": len(d), "n_seeds_better": n_better, "count_statement": f"{n_better}/{len(d)} 个评估种子上更好", "median_paired_diff": float(np.median(d)), "mean_paired_diff": float(d.mean()), "wilcoxon_statistic": float(res.statistic), "wilcoxon_p_two_sided_exact": float(res.pvalue), "rank_biserial": float((wp - wn) / (wp + wn)) if (wp + wn) else 0.0, "paired_diffs": [float(v) for v in d], } def main() -> int: real = RAW["real"] arms = ["null_no_intervention", "random_budget_matched", "degree_top_b", "extremity_top_b", "extremity_double_budget", "offline_rl_fqi", "random_variant_a", "random_variant_b"] out = {"eval_windows": RAW["eval_windows"], "n_eval_windows": RAW["n_eval_windows"], "n_nodes": RAW["n_nodes"], "eval_seeds": RAW["eval_seeds"], "budget_total_standard_arms": RAW["budget_total_standard_arms"]} null_end = series(real, "null_no_intervention", "P_end") out["per_arm"] = {} for a in arms: pe, pa = series(real, a, "P_end"), series(real, a, "P_auc") spend = [real[a][s]["spend"] for s in SEEDS] budget = real[a][SEEDS[0]]["budget"] out["per_arm"][a] = { "P_end_median": float(np.median(pe)), "P_end_mean": float(pe.mean()), "P_end_std": float(pe.std(ddof=1)), "P_auc_median": float(np.median(pa)), "reduction_vs_null_median": float(np.median(null_end - pe)), "relative_reduction_vs_null_median_pct": float(np.median((null_end - pe) / null_end) * 100), "budget_declared": budget, "spend_min": int(min(spend)), "spend_max": int(max(spend)), "budget_respected": bool(max(spend) <= budget), } tests = {} tests["H1_primary_rl_vs_extremity_P_end"] = paired( series(real, "offline_rl_fqi", "P_end"), series(real, "extremity_top_b", "P_end"), "offline_rl_fqi 对 extremity_top_b（主对照）", "P_end") tests["H1_secondary_rl_vs_extremity_P_auc"] = paired( series(real, "offline_rl_fqi", "P_auc"), series(real, "extremity_top_b", "P_auc"), "offline_rl_fqi 对 extremity_top_b（次要终点 AUC）", "P_auc") for a in ["random_budget_matched", "degree_top_b", "extremity_top_b", "offline_rl_fqi"]: tests[f"H2_sanity_{a}_vs_null_P_end"] = paired( series(real, a, "P_end"), null_end, f"{a} 对零模型 null_no_intervention", "P_end") tests["H3_ruler_positive_double_vs_single_budget"] = paired( series(real, "extremity_double_budget", "P_end"), series(real, "extremity_top_b", "P_end"), "extremity_double_budget 对 extremity_top_b（正对照：预算翻倍，方向已知）", "P_end") tests["H4_ruler_negative_random_a_vs_b"] = paired( series(real, "random_variant_a", "P_end"), series(real, "random_variant_b", "P_end"), "random_variant_a 对 random_variant_b（负对照：应无系统差异）", "P_end") out["tests"] = tests verdict = {} h1 = tests["H1_primary_rl_vs_extremity_P_end"] verdict["H1_primary"] = ("supported" if h1["wilcoxon_p_two_sided_exact"] < 0.05 and h1["median_paired_diff"] > 0 else "refuted") verdict["H2_sanity"] = "supported" if all( tests[f"H2_sanity_{a}_vs_null_P_end"]["wilcoxon_p_two_sided_exact"] < 0.05 and tests[f"H2_sanity_{a}_vs_null_P_end"]["median_paired_diff"] > 0 for a in ["random_budget_matched", "degree_top_b", "extremity_top_b", "offline_rl_fqi"] ) else "refuted" h3 = tests["H3_ruler_positive_double_vs_single_budget"] verdict["H3_ruler_positive"] = ("supported" if h3["wilcoxon_p_two_sided_exact"] < 0.05 and h3["median_paired_diff"] > 0 else "refuted") h4 = tests["H4_ruler_negative_random_a_vs_b"] verdict["H4_ruler_negative"] = ("supported" if h4["wilcoxon_p_two_sided_exact"] >= 0.05 else "refuted") out["verdict"] = verdict sh = RAW["sign_shuffled_graph_null_model"] sh_null, sh_ext = series(sh, "null_no_intervention", "P_end"), series(sh, "extremity_top_b", "P_end") out["graph_null_model_sign_shuffled"] = { "what": "保持拓扑、时间戳次序、正负边总数，只把全部边的 ±1 符号全局随机置换", "real_null_P_end_median": float(np.median(null_end)), "shuffled_null_P_end_median": float(np.median(sh_null)), "real_extremity_reduction_median": float(np.median(null_end - series(real, "extremity_top_b", "P_end"))), "shuffled_extremity_reduction_median": float(np.median(sh_null - sh_ext)), "null_level_shift_test": paired(sh_null, null_end, "符号打乱图 对 真实图（零模型臂的极化水平）", "P_end"), } out["budget_accounting"] = { a: {"declared": out["per_arm"][a]["budget_declared"], "max_spend_over_seeds": out["per_arm"][a]["spend_max"], "respected": out["per_arm"][a]["budget_respected"]} for a in arms} out["all_budgets_respected"] = all(v["respected"] for v in out["budget_accounting"].values()) (PKG / "work/results/analysis.json").write_text( json.dumps(out, ensure_ascii=False, indent=1), encoding="utf-8") print(json.dumps({"verdict": verdict, "H1": {k: h1[k] for k in ("count_statement", "median_paired_diff", "wilcoxon_p_two_sided_exact", "rank_biserial")}, "budgets_ok": out["all_budgets_respected"]}, ensure_ascii=False, indent=1)) return 0 if __name__ == "__main__": raise SystemExit(main())