"""生成诊断（方法.md 6.12）的预期值 JSON（diagnosis.json）。 用法（在仓库根目录）： python3 tests/fixtures/rules/gen_diagnosis.py 只用标准库；依赖同目录的 rules_ref.py 里的基础件（折线、有理数分布、Stieltjes 精确积分、比特、惊讶、分歧量）。 重跑结果逐位相同。 独立性：诊断的成分模型、方向收益、信息价值、价值减成本、选实验、结算都在本文件里只照 docs/方法.md 6.12 和 6.7、6.9 的文字写成，没有读过 src/rules/diagnosis.ts 和 src/rules/value.ts。 算法（有理数精确，浮点输入先精确转成有理数）： - 每个成分每个实验的分布：各条解释的折线加这个实验自己的开放预留（6.3"诊断"），换成收益 g = 改进量（x − b，越小越好时 b − x）； - E(g_d)：各成分期望按起始权重平均； - 全做（k = 实验总数）的 VOI：Σ_j w_j · E_j[max(0, max_d g_d)]（rules_ref.expect_max 的 floor0 模式）− max(0, max_d E g_d)； - 只做一个实验 d0 的 VOI（exactK1）：d0 的各成分在相邻折点之间密度是常数，所以后验权重在每一段里是常数， 没做的实验的期望收益 c_e 在每一段里是常数，max(0, v, c_e…) 对 v 的积分有闭式； - 价值 = u·m′·VOI，R̂ = ⌈(余额 − 1) ÷ K⌉，成本 = 2·max(0, ĝ₀)·[(K+1)/(K+2) − K/(K+1)]； - 推出的方向：成分新权重 ∝ 起始权重 × 各有效实验在实测处的密度（所在段斜率，段左闭右开）。 生成的文件：diagnosis.json，格式： {"generator", "tolerance", "scenarios": [{"name", "section", "input", "expected", "note"}]} - 场景 input：{"b", "direction", "lambda", "experiments": ["d1", ...], "explanations": [{"id", "author", "ranges": {"d1": {"q1","q5","q9"}, ...}}]}， 选实验场景另带 "balance"、"batchCap"、"maxExperiments"、"u"、"decay"，结算场景另带 "epsilon" 和 "outcomes"。 - 分数一律是原值；收益、VOI 是收益刻度（改进量）。 """ import json import math import os import sys from fractions import Fraction as Fr from itertools import combinations sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) import rules_ref as R # noqa: E402 HERE = os.path.dirname(os.path.abspath(__file__)) F = float def tri(t): return {"q1": t[0], "q5": t[1], "q9": t[2]} def ids_key(s): return R.factor_key(s) class Model: """诊断的成分模型（有理数）。explanations 按编号排好；成分 = 各解释 + 开放预留。""" def __init__(self, b, direction, lam, experiments, explanations): self.b, self.direction, self.lam = fr(b), direction, fr(lam) self.exps = sorted(experiments, key=ids_key) self.expl = sorted(explanations, key=lambda e: ids_key(e["id"])) authors = [e["author"] for e in self.expl] self.w = R.explanation_weights(authors, lam) + [fr(lam)] # raw[d][j]：分数原值上的点列；gain[d][j]：收益（改进量）上的点列 self.raw, self.gain = {}, {} for d in self.exps: pts = [R.polyline(*e["ranges"][d]) for e in self.expl] lo, hi = R.reserve_bounds(pts) pts.append(R.uniform_pts(lo, hi)) self.raw[d] = pts self.gain[d] = [R.improvement_pts(p, self.b, direction) for p in pts] self.mu = {d: [R.mean(R.single(p)) for p in self.gain[d]] for d in self.exps} self.prior = {d: sum((wj * m for wj, m in zip(self.w, self.mu[d])), Fr(0)) for d in self.exps} def g0(self): return max(Fr(0), max(self.prior.values())) def voi_all(self): total = Fr(0) for j, wj in enumerate(self.w): if wj == 0: continue total += wj * R.expect_max([R.single(self.gain[d][j]) for d in self.exps], floor0=True) return total - self.g0() def voi_k1(self, d0): """只做 d0：每一段里后验权重和没做的实验的期望收益都是常数，对 max(0, v, c…) 精确积分。""" bps = sorted(set().union(*[set(x for x, _ in p) for p in self.gain[d0]])) total = Fr(0) for a, b in zip(bps, bps[1:]): dens = [R.pts_left(p, b) - R.pts_left(p, a) for p in self.gain[d0]] dens = [x / (b - a) for x in dens] marg = sum((wj * x for wj, x in zip(self.w, dens)), Fr(0)) if marg == 0: continue c0 = Fr(0) for e in self.exps: if e == d0: continue num = sum((wj * x * m for wj, x, m in zip(self.w, dens, self.mu[e])), Fr(0)) c0 = max(c0, num / marg) t = min(max(c0, a), b) total += marg * (c0 * (t - a) + (b * b - t * t) / 2) return total - self.g0() def mixture(self, d): """实验 d 的成分混合（分数原值）。""" return R.Dist(list(zip(self.w, self.raw[d]))) def fr(x): return R.fr(x) def model_of(inp): expl = [{"id": e["id"], "author": e["author"], "ranges": {d: (r["q1"], r["q5"], r["q9"]) for d, r in e["ranges"].items()}} for e in inp["explanations"]] return Model(inp["b"], inp["direction"], inp["lambda"], inp["experiments"], expl) def explanation(eid, author, ranges): return {"id": eid, "author": author, "ranges": {d: tri(t) for d, t in ranges.items()}} def m_prime(balance, K, r): R_ = -(-(balance - 1) // K) return R_, sum((fr(r) ** j for j in range(R_)), Fr(0)) def cost_of(g0, K): return 2 * g0 * (Fr(K + 1, K + 2) - Fr(K, K + 1)) S = [] def select_case(name, inp, balance, K, D, u, r, note): inp = dict(inp, balance=balance, batchCap=K, maxExperiments=D, u=u, decay=r) m = model_of(inp) rows = [] if balance >= 2: size = min(len(m.exps), D) R_, mp = m_prime(balance, K, r) c = cost_of(m.g0(), K) for combo in combinations(m.exps, size): if size == len(m.exps): voi = m.voi_all() elif size == 1: voi = m.voi_k1(combo[0]) else: voi = None # 网格近似，没有精确值，只检查结构 if voi is None: rows.append({"experiments": list(combo), "voi": None}) continue value = fr(u) * mp * voi rows.append({"experiments": list(combo), "voi": F(voi), "rHat": R_, "mPrime": F(mp), "value": F(value), "cost": F(c), "net": F(value - c)}) expected = { "priorGains": {d: F(m.prior[d]) for d in m.exps}, "g0": F(m.g0()), "weights": [F(x) for x in m.w], "combos": rows, } if all(r_["voi"] is not None for r_ in rows) and rows: eligible = [r_ for r_ in rows if r_["net"] > 1e-12] if eligible: top = max(r_["net"] for r_ in eligible) near = [r_ for r_ in eligible if r_["net"] >= top - 1e-12] best = min(near, key=lambda r_: [ids_key(x) for x in r_["experiments"]]) expected["diagnose"] = True expected["chosen"] = best["experiments"] expected["chosenNet"] = best["net"] else: expected["diagnose"] = False expected["chosen"] = [] expected["chosenNet"] = 0.0 elif not rows: expected["diagnose"] = False expected["chosen"] = [] expected["chosenNet"] = 0.0 S.append({"name": name, "section": "6.12", "input": inp, "expected": expected, "note": note}) def settle_case(name, inp, outcomes, note): inp = dict(inp, outcomes=outcomes) m = model_of(inp) eps = fr(inp["epsilon"]) delta = F(eps / 100) b = F(m.b) scores = {d: (o["score"] if o["kind"] == "valid" else None) for d, o in outcomes.items()} rows = {} for d in m.exps: dist = m.mixture(d) row = {"outcome": outcomes[d]["kind"]} row["merged"] = {"q1": F(R.quantile(dist, Fr(1, 10))), "q5": F(R.quantile(dist, Fr(1, 2))), "q9": F(R.quantile(dist, Fr(9, 10)))} if scores[d] is not None: x = scores[d] row["gain"] = F(R.improve(x, m.b, m.direction)) row["worldBits"] = R.world_bits(dist, b, delta) sp = R.surprise(dist, x, b, delta) row["surprise"] = sp persons = {} for e in m.expl: persons.setdefault(e["author"], []).append(R.polyline(*e["ranges"][d])) person_dists = [R.Dist([(Fr(1, len(ps)), p) for p in ps]) for _, ps in sorted(persons.items(), key=lambda kv: ids_key(kv[0]))] row["disagreement"] = R.disagreement(person_dists, b, delta)[0] rows[d] = row # 确认与否定：只有有效实验判，落在 [q1, q9]（含两端）确认。 judgments = {} for e in m.expl: judgments[e["id"]] = {d: ("confirmed" if e["ranges"][d][0] <= scores[d] <= e["ranges"][d][2] else "denied") for d in m.exps if scores[d] is not None} # 诊断推出的方向 new = [] for j, wj in enumerate(m.w): v = wj for d in m.exps: if scores[d] is not None: v *= R.pts_density(m.raw[d][j], fr(scores[d])) new.append(v) tot = sum(new, Fr(0)) kept = tot == 0 w2 = m.w if kept else [v / tot for v in new] gains = {} for d in m.exps: if scores[d] is not None: gains[d] = {"gain": F(R.improve(scores[d], m.b, m.direction)), "measured": True} else: gains[d] = {"gain": F(sum((wj * mu for wj, mu in zip(w2, m.mu[d])), Fr(0))), "measured": False} best = m.exps[0] for d in m.exps[1:]: if gains[d]["gain"] > gains[best]["gain"]: best = d expected = {"experiments": rows, "judgments": judgments, "direction": {"weights": [F(x) for x in w2], "weightsKept": kept, "gains": gains, "experiment": best, "gain": gains[best]["gain"], "found": gains[best]["gain"] > 0}} S.append({"name": name, "section": "6.12 结算；6.9", "input": inp, "expected": expected, "note": note}) # ------------------------------------------------------------------ 场景 def base(explanations, experiments, b=0.0, direction="maximize", lam=0.1): return {"b": b, "direction": direction, "lambda": lam, "experiments": experiments, "explanations": explanations} # 1. 两条解释对同一个实验预测相反：A 说涨 0.2 左右，B 说跌 0.2 左右。 OPP = [explanation("e1", "researcher-1", {"d1": (0.15, 0.2, 0.25)}), explanation("e2", "researcher-2", {"d1": (-0.25, -0.2, -0.15)})] select_case("opposite_lambda0", base(OPP, ["d1"], lam=0.0), 7, 2, 4, 0.5, 0.5, "λ = 0，两个成分的范围不相交：做了就知道是哪一条。手算：做了以后 E[max(0, g)] = 0.5 × 0.2 = 0.1，" "不做时 max(0, E g) = max(0, 0.5×0.2 + 0.5×(−0.2)) = 0，VOI = 0.1。") select_case("opposite_lambda01", base(OPP, ["d1"]), 7, 2, 4, 0.5, 0.5, "同上，λ = 0.1，开放预留铺在两端之外。") # 2. 两条解释预测一样：λ = 0 时信息价值为 0；λ = 0.1 时只剩开放预留带来的一点。 SAME = [explanation("e1", "researcher-1", {"d1": (0.1, 0.2, 0.3)}), explanation("e2", "researcher-2", {"d1": (0.1, 0.2, 0.3)})] select_case("same_lambda0", base(SAME, ["d1"], lam=0.0), 7, 2, 4, 0.5, 0.5, "λ = 0，两条解释一样，做不做都不改变看法：VOI = 0（手算：g 的取值范围 0 到 0.4 全在 0 以上，E[max(0, g)] = E g = 0.2）。净值 = −成本 < 0，不诊断。") select_case("same_lambda01", base(SAME, ["d1"]), 7, 2, 4, 0.5, 0.5, "λ = 0.1：解释们的范围都在 0 以上，只有开放预留（铺在 −0.2 到 0.6）能给出负的结果，所以 VOI 只剩预留带来的一小部分。") # 3. 三个实验、三条解释，全做（3 ≤ D）。 MULTI = [explanation("e1", "researcher-1", {"d1": (0.02, 0.05, 0.08), "d2": (-0.05, 0.0, 0.05), "d3": (0.0, 0.1, 0.2)}), explanation("e2", "researcher-2", {"d1": (0.0, 0.03, 0.06), "d2": (0.1, 0.2, 0.3), "d3": (-0.1, -0.05, 0.0)}), explanation("e3", "researcher-3", {"d1": (-0.02, 0.04, 0.1), "d2": (-0.2, -0.1, 0.0), "d3": (0.05, 0.1, 0.15)})] select_case("multi_all", base(MULTI, ["d1", "d2", "d3"]), 9, 2, 4, 0.5, 0.5, "三个实验不超过 D = 4，只有全做一个组合。R̂ = ⌈8 ÷ 2⌉ = 4，m′ = 1.875，成本 = 2ĝ₀ × (3/4 − 2/3)。") # 4. 实验数超过 D：只比较装满的组合。D = 1 时是三个单实验组合（有精确值）；D = 2 时只有三个两实验组合（只查结构）。 select_case("over_cap_d1", base(MULTI, ["d1", "d2", "d3"]), 9, 2, 1, 0.5, 0.5, "三个实验超过 D = 1：只有三个单实验组合，没有空组合；精确值用逐段闭式积分。") select_case("over_cap_d2", base(MULTI, ["d1", "d2", "d3"]), 9, 2, 2, 0.5, 0.5, "三个实验超过 D = 2：只有三个两实验组合，没有单实验组合；只检查组合列表（网格近似没有精确值）。") # 5. 余额不到 2 次不诊断。 select_case("balance_one", base(MULTI, ["d1", "d2", "d3"]), 1, 2, 4, 0.5, 0.5, "余额 1：没有可选组合，不诊断。") # 6. 平票按实验编号：d2 和 d10 的预测完全一样，D = 1 时两个单实验组合净值相同，取 d2（自然序，d2 排在 d10 前面）。 TIE = [explanation("e1", "researcher-1", {"d2": (-0.2, 0.0, 0.2), "d10": (-0.2, 0.0, 0.2)}), explanation("e2", "researcher-2", {"d2": (0.1, 0.3, 0.5), "d10": (0.1, 0.3, 0.5)})] select_case("tie_by_id", base(TIE, ["d10", "d2"]), 7, 2, 1, 0.5, 0.5, "两个实验的预测完全一样，净值相同（差不到 1e-12），取实验编号小的 d2。") # 7. 结算。 def shifted(expls, b): """把以基线为 0 的改进量范围搬成分数原值（基线 b）。""" return [dict(e, ranges={d: {k: v + b for k, v in r.items()} for d, r in e["ranges"].items()}) for e in expls] SETTLE = shifted(MULTI[:2] + [explanation("e3", "researcher-2", {"d1": (-0.02, 0.04, 0.1), "d2": (-0.2, -0.1, 0.0), "d3": (0.05, 0.1, 0.15)})], 1.0) SINP = dict(base(SETTLE, ["d1", "d2", "d3"], b=1.0), epsilon=0.01) settle_case("settle_mixed", SINP, {"d1": {"kind": "valid", "score": 1.045}, "d2": {"kind": "insufficient"}, "d3": {"kind": "untested"}}, "d1 有效，d2 证据不足（不判、不记账），d3 没做。研究者 2 交了两条解释，各占他那份的一半。") settle_case("settle_two_valid", SINP, {"d1": {"kind": "valid", "score": 1.07}, "d2": {"kind": "valid", "score": 1.25}, "d3": {"kind": "untested"}}, "两个有效实验：权重按两处密度连乘更新；d2 的实测 +0.25 落在 e2 的范围里、其余两条之外。") settle_case("settle_no_direction", dict(SINP), {"d1": {"kind": "valid", "score": 0.97}, "d2": {"kind": "valid", "score": 0.8}, "d3": {"kind": "valid", "score": 0.9}}, "全部实验的收益都不为正（−0.03、−0.2、−0.1）：最大收益 −0.03 ≤ 0，记没有指出能涨分的方向。") settle_case("settle_all_insufficient", SINP, {"d1": {"kind": "insufficient"}, "d2": {"kind": "insufficient"}, "d3": {"kind": "insufficient"}}, "没有一个有效实验：不判、不记账，权重保持起始权重，方向按起始权重估的期望收益。") settle_case("settle_outside_all", SINP, {"d1": {"kind": "valid", "score": 3.0}, "d2": {"kind": "untested"}, "d3": {"kind": "untested"}}, "d1 实测 +2.0 落在所有成分（含开放预留）的范围之外，权重无从更新：保持起始权重（weightsKept）。") MIN_EXPL = [explanation("e1", "researcher-1", {"d1": (0.95, 0.9, 0.85)[::-1], "d2": (0.9, 1.0, 1.1)}), explanation("e2", "researcher-2", {"d1": (0.8, 0.85, 0.9), "d2": (0.95, 1.05, 1.15)})] MIN_INP = dict(base(MIN_EXPL, ["d1", "d2"], b=1.0, direction="minimize"), epsilon=0.01) select_case("minimize_all", dict(base(MIN_EXPL, ["d1", "d2"], b=1.0, direction="minimize")), 7, 2, 4, 0.5, 0.5, "越小越好：收益 = b − x。d1 两条解释都预测降到 0.9 左右（收益约 +0.1），d2 预测上升。") settle_case("settle_minimize", MIN_INP, {"d1": {"kind": "valid", "score": 0.88}, "d2": {"kind": "valid", "score": 1.05}}, "越小越好：d1 实测 0.88，收益 = 1.0 − 0.88 = 0.12；d2 实测 1.05，收益 −0.05。") doc = {"generator": "gen_diagnosis.py", "tolerance": { "exact": "有理数精确算出的量（成分权重、事前期望收益、ĝ₀、全做的 VOI、R̂、m′、成本、价值、净值、收益、更新后的权重）：绝对误差 1e-12", "exactK1": "只做一个实验的 VOI 是逐段闭式积分，实现用部分格网格近似：|差| ≤ 1e-9 + 2e-4·|精确值|", "bits": "世界比特、惊讶、分歧量是浮点：绝对误差 1e-9；合并范围是分位点：绝对误差 1e-9", "discrete": "是否诊断、选中的实验、确认与否定、有没有指出方向、实验的结局：必须完全相同"}, "scenarios": S} with open(os.path.join(HERE, "diagnosis.json"), "w", encoding="utf-8") as f: json.dump(doc, f, ensure_ascii=False, indent=1) f.write("\n") print(len(S), "个场景")