#!/usr/bin/env python3 """ protocol_engine.py -- protocol design engine. Python standard library only. Input : one JSON config describing a registered construct test (stages/arms, effect thresholds, alpha, target power, gate rules). Output: a markdown design report (stdout) containing, per arm: 1. required sample size, closed form; 2. Monte Carlo power at that size (cross-check of the closed form); 3. power at the currently registered size; plus a falsifier-fire audit, an input-separation gate budget, a scale table, and a registration skeleton filled from the config. Deterministic: one seed in the config; every simulation stream derives from it. Re-running the same command reproduces every number. Run: python3 protocol_engine.py config_summarization_bias.json """ import json import math import random import sys from statistics import NormalDist ND = NormalDist() # ---------------------------------------------------------------------- # Helpers # ---------------------------------------------------------------------- def z_alpha(alpha, sides): return ND.inv_cdf(1 - alpha / sides) def z_power(power): return ND.inv_cdf(power) def cohens_h(p1, p2): return 2 * math.asin(math.sqrt(p1)) - 2 * math.asin(math.sqrt(p2)) def rate_from_h(p_base, h): """Rate p2 > p_base such that cohens_h(p2, p_base) == h.""" return math.sin(math.asin(math.sqrt(p_base)) + h / 2.0) ** 2 def rate_from_or(p_base, odds_ratio): odds = odds_ratio * p_base / (1.0 - p_base) return odds / (1.0 + odds) def binom_pmf(k, n, p): return math.comb(n, k) * p ** k * (1.0 - p) ** (n - k) def binom_cdf(k, n, p): if k < 0: return 0.0 if k >= n: return 1.0 return sum(binom_pmf(i, n, p) for i in range(k + 1)) def binom_draw(rng, n, p): if hasattr(rng, "binomialvariate"): # Python >= 3.12 return rng.binomialvariate(n, p) return sum(1 for _ in range(n) if rng.random() < p) def fmt(x, digits=3): return f"{x:.{digits}f}" # ---------------------------------------------------------------------- # Closed-form sample sizes # ---------------------------------------------------------------------- def n_two_proportions(h, alpha, sides, power): """Per-group n for a two-proportion test, arcsine (Cohen's h) form.""" z = z_alpha(alpha, sides) + z_power(power) return math.ceil(2.0 * (z / h) ** 2) def n_one_proportion(p0, p1, alpha, sides, power): za = z_alpha(alpha, sides) zb = z_power(power) num = za * math.sqrt(p0 * (1 - p0)) + zb * math.sqrt(p1 * (1 - p1)) return math.ceil((num / abs(p1 - p0)) ** 2) def _beta_cf(a, b, x): """Continued fraction for the incomplete beta function.""" qab, qap, qam = a + b, a + 1.0, a - 1.0 c = 1.0 d = 1.0 - qab * x / qap if abs(d) < 1e-300: d = 1e-300 d = 1.0 / d h = d for m in range(1, 201): m2 = 2 * m aa = m * (b - m) * x / ((qam + m2) * (a + m2)) d = 1.0 + aa * d if abs(d) < 1e-300: d = 1e-300 c = 1.0 + aa / c if abs(c) < 1e-300: c = 1e-300 d = 1.0 / d h *= d * c aa = -(a + m) * (qab + m) * x / ((a + m2) * (qap + m2)) d = 1.0 + aa * d if abs(d) < 1e-300: d = 1e-300 c = 1.0 + aa / c if abs(c) < 1e-300: c = 1e-300 d = 1.0 / d delta = d * c h *= delta if abs(delta - 1.0) < 3e-16: break return h def _reg_incomplete_beta(x, a, b): if x <= 0.0: return 0.0 if x >= 1.0: return 1.0 bt = math.exp(math.lgamma(a + b) - math.lgamma(a) - math.lgamma(b) + a * math.log(x) + b * math.log(1.0 - x)) if x < (a + 1.0) / (a + b + 2.0): return bt * _beta_cf(a, b, x) / a return 1.0 - bt * _beta_cf(b, a, 1.0 - x) / b def t_cdf(t, df): """CDF of the central Student t distribution (incomplete beta).""" x = df / (df + t * t) ib = _reg_incomplete_beta(x, df / 2.0, 0.5) return 1.0 - 0.5 * ib if t >= 0 else 0.5 * ib def t_crit(alpha, sides, df): target = 1.0 - alpha / sides lo, hi = 0.0, 1000.0 while hi - lo > 1e-12: mid = (lo + hi) / 2.0 if t_cdf(mid, df) < target: lo = mid else: hi = mid return (lo + hi) / 2.0 def _chi_pdf(v, df): if v <= 0.0: return 0.0 k = df / 2.0 return math.exp((k - 1.0) * math.log(v) - v / 2.0 - k * math.log(2.0) - math.lgamma(k)) def _adaptive_simpson(f, a, b, eps, depth): c = (a + b) / 2.0 fa, fb, fc = f(a), f(b), f(c) whole = (b - a) * (fa + 4.0 * fc + fb) / 6.0 def step(a, b, eps, whole, fa, fb, fc, depth): c = (a + b) / 2.0 lm, rm = (a + c) / 2.0, (c + b) / 2.0 flm, frm = f(lm), f(rm) left = (c - a) * (fa + 4.0 * flm + fc) / 6.0 right = (b - c) * (fc + 4.0 * frm + fb) / 6.0 delta = left + right - whole if depth <= 0 or abs(delta) <= 15.0 * eps: return left + right + delta / 15.0 return (step(a, c, eps / 2.0, left, fa, fc, flm, depth - 1) + step(c, b, eps / 2.0, right, fc, fb, frm, depth - 1)) return step(a, b, eps, whole, fa, fb, fc, depth) def nct_cdf(t, df, ncp): """CDF of the noncentral t distribution, from its definition T = (Z + ncp) / sqrt(V/df): P(T <= t) = E_V[Phi(t sqrt(V/df) - ncp)] over the chi-square density of V, by adaptive Simpson quadrature.""" if t < 0: return 1.0 - nct_cdf(-t, df, -ncp) f = lambda v: ND.cdf(t * math.sqrt(v / df) - ncp) * _chi_pdf(v, df) hi = df + 20.0 * math.sqrt(2.0 * df) + 60.0 return _adaptive_simpson(f, 0.0, hi, 1e-12, 40) def paired_t_power_exact(n, dz, alpha, sides): """Exact power of the paired (one-sample) t test on differences: noncentral t at df = n-1, noncentrality dz*sqrt(n).""" df = float(n - 1) ncp = dz * math.sqrt(n) crit = t_crit(alpha, sides, df) power = 1.0 - nct_cdf(crit, df, ncp) if sides == 2: power += nct_cdf(-crit, df, ncp) return power def n_paired_dz(dz, alpha, sides, power): """Smallest n whose exact paired-t power reaches the target -- an exact search over the noncentral t power, not the normal closed form.""" n = 2 while paired_t_power_exact(n, dz, alpha, sides) < power: n += 1 return n # ---------------------------------------------------------------------- # Two-proportion z test: p-value, Monte Carlo power, exact power # ---------------------------------------------------------------------- def two_prop_z_pvalue(x1, n1, x2, n2): pooled = (x1 + x2) / (n1 + n2) se = math.sqrt(pooled * (1 - pooled) * (1.0 / n1 + 1.0 / n2)) if se == 0.0: return 1.0 z = (x1 / n1 - x2 / n2) / se return 2.0 * (1.0 - ND.cdf(abs(z))) def mc_power_two_prop(n, p1, p2, alpha, reps, seed): rng = random.Random(seed) hits = 0 for _ in range(reps): x = binom_draw(rng, n, p1) y = binom_draw(rng, n, p2) if two_prop_z_pvalue(x, n, y, n) < alpha: hits += 1 p = hits / reps return p, math.sqrt(p * (1 - p) / reps) def exact_power_two_prop_z(n, p1, p2, alpha): """Exact power of the pooled z test by enumerating all outcome pairs. Practical for small n (used at the registered n = 20).""" probs1 = [binom_pmf(x, n, p1) for x in range(n + 1)] probs2 = [binom_pmf(y, n, p2) for y in range(n + 1)] total = 0.0 for x in range(n + 1): for y in range(n + 1): if two_prop_z_pvalue(x, n, y, n) < alpha: total += probs1[x] * probs2[y] return total def fisher_pvalue(x, n, y, m): """Two-sided Fisher exact p-value for [[x, n-x], [y, m-y]].""" s = x + y lo = max(0, s - m) hi = min(n, s) denom = math.comb(n + m, s) def hyper(k): return math.comb(n, k) * math.comb(m, s - k) / denom p_obs = hyper(x) return sum(hyper(k) for k in range(lo, hi + 1) if hyper(k) <= p_obs * (1 + 1e-12)) def exact_power_fisher(n, p1, p2, alpha): probs1 = [binom_pmf(x, n, p1) for x in range(n + 1)] probs2 = [binom_pmf(y, n, p2) for y in range(n + 1)] total = 0.0 for x in range(n + 1): for y in range(n + 1): if fisher_pvalue(x, n, y, n) < alpha: total += probs1[x] * probs2[y] return total # ---------------------------------------------------------------------- # One-sample binomial test: exact power # ---------------------------------------------------------------------- def one_sample_critical(n, p0, alpha, sides): """Equal-tail exact critical counts. Reject if X <= lo or X >= hi.""" if sides == 2: lo = -1 for k in range(n + 1): if binom_cdf(k, n, p0) <= alpha / 2: lo = k else: break hi = n + 1 for k in range(n, -1, -1): if 1.0 - binom_cdf(k - 1, n, p0) <= alpha / 2: hi = k else: break return lo, hi hi = n + 1 for k in range(n, -1, -1): if 1.0 - binom_cdf(k - 1, n, p0) <= alpha: hi = k else: break return -1, hi def exact_power_one_sample(n, p0, p1, alpha, sides): lo, hi = one_sample_critical(n, p0, alpha, sides) power = 0.0 if lo >= 0: power += binom_cdf(lo, n, p1) if hi <= n: power += 1.0 - binom_cdf(hi - 1, n, p1) return power # ---------------------------------------------------------------------- # Paired-difference test (Stage 2): Monte Carlo power # ---------------------------------------------------------------------- def mc_power_paired_dz(n, dz, alpha, sides, reps, seed): """Independent check of the exact paired-t size: simulate the paired t test itself. Differences D ~ Normal(dz, 1) are drawn, the variance is ESTIMATED from each sample (sample sd, ddof=1), and the statistic is compared with the t critical value -- the t behaviour the closed form must answer to. Because the check shares no formula with n_paired_dz, it can falsify it: run at the old normal-approximation size n = 32 (dz = 0.5, alpha 0.05) it returns ~0.78, below the 0.80 target the old size claimed.""" rng = random.Random(seed) crit = t_crit(alpha, sides, float(n - 1)) hits = 0 for _ in range(reps): diffs = [rng.gauss(dz, 1.0) for _ in range(n)] mean = sum(diffs) / n var = sum((d - mean) ** 2 for d in diffs) / (n - 1) if var <= 0.0: continue t = mean / math.sqrt(var / n) if (sides == 2 and abs(t) > crit) or (sides == 1 and t > crit): hits += 1 p = hits / reps return p, math.sqrt(p * (1 - p) / reps) # ---------------------------------------------------------------------- # Input-separation gate budget # ---------------------------------------------------------------------- def gate_build_count(n_required, fail_rate, assurance): """Smallest build count B with P(Binomial(B, 1-fail) >= n_required) >= assurance. Exact search over the binomial CDF.""" p_pass = 1.0 - fail_rate B = n_required while True: if 1.0 - binom_cdf(n_required - 1, B, p_pass) >= assurance: return B B += 1 def mc_gate_check(B, n_required, fail_rate, reps, seed): rng = random.Random(seed) hits = sum(1 for _ in range(reps) if binom_draw(rng, B, 1.0 - fail_rate) >= n_required) p = hits / reps return p, math.sqrt(p * (1 - p) / reps) def exact_search_one_sample(p0, p1, alpha, sides, target, n_start): """Smallest n >= n_start whose exact binomial power reaches target. The exact test is conservative at small n, so this sits at or above the closed-form value. Search is capped at n_start + 120.""" n = n_start while n <= n_start + 120: if exact_power_one_sample(n, p0, p1, alpha, sides) >= target: return n n += 1 return n_start # ---------------------------------------------------------------------- # Standards # standard : the config's alpha and target power. # decisive : alpha = 0.01, power = 0.95 -- the scale at which a verdict, # either way, is hard to argue with. # ---------------------------------------------------------------------- DECISIVE_ALPHA = 0.01 DECISIVE_POWER = 0.95 def run(cfg): alpha = cfg["alpha"] sides = cfg.get("sides", 2) power_t = cfg["target_power"] reps = cfg.get("mc_reps", 20000) seed0 = cfg.get("seed", 1) out = [] out.append(f"# Design report -- {cfg['construct_name']}") out.append("") out.append(f"Source protocol: {cfg['source_protocol']}") out.append("") out.append(f"- Standard design point: alpha = {alpha} ({sides}-sided), " f"power = {power_t}") out.append(f"- Decisive design point: alpha = {DECISIVE_ALPHA} " f"({sides}-sided), power = {DECISIVE_POWER}") out.append(f"- Monte Carlo reps per cell = {reps}; base seed = {seed0}; " f"one stream per row. Re-running reproduces every number.") out.append("") # ---- Table 1: per-arm sizing at both design points ---------------- out.append("## 1. Required sample size per arm") out.append("") out.append("| Arm | Test | Effect as specified | n at standard point | " "MC power at that n | n at decisive point | Current n | " "Power at current n |") out.append("|---|---|---|---|---|---|---|---|") res = {} # Stage 3 (core): two-proportion, model vs human told-preference rate. s3 = cfg["stage3"] p_base = s3["human_told_preference_rate"] p_h = rate_from_h(p_base, s3["registered_thresholds"]["cohens_h"]) p_or = rate_from_or(p_base, s3["registered_thresholds"]["odds_ratio"]) h_or = abs(cohens_h(p_or, p_base)) h_gov = min(s3["registered_thresholds"]["cohens_h"], h_or) p_alt = p_h if s3["registered_thresholds"]["cohens_h"] <= h_or else p_or n_std = n_two_proportions(h_gov, alpha, sides, power_t) n_dec = n_two_proportions(h_gov, DECISIVE_ALPHA, sides, DECISIVE_POWER) mc_std, se_std = mc_power_two_prop(n_std, p_alt, p_base, alpha, reps, seed0 + 101) n_cur = s3["current_pairs"] mc_cur, se_cur = mc_power_two_prop(n_cur, p_alt, p_base, alpha, reps, seed0 + 102) res["stage3"] = dict(n_std=n_std, n_dec=n_dec, p_alt=p_alt, p_base=p_base, h=h_gov) out.append( f"| Stage 3, evaluative (core) | two-proportion z, pooled | " f"h = {s3['registered_thresholds']['cohens_h']} -> rate {p_base} vs " f"{fmt(p_h)}; OR = {s3['registered_thresholds']['odds_ratio']} -> rate " f"{fmt(p_or)} (h = {fmt(h_or)}); governing h = {fmt(h_gov)} | " f"{n_std} per group | {fmt(mc_std)} +/- {fmt(se_std)} | {n_dec} per group | " f"{n_cur} per group | {fmt(mc_cur)} +/- {fmt(se_cur)} |") # Stage 2: paired SI difference. s2 = cfg["stage2"] dz = s2["planning_effect_dz"] n2_std = n_paired_dz(dz, alpha, sides, power_t) n2_dec = n_paired_dz(dz, DECISIVE_ALPHA, sides, DECISIVE_POWER) mc2, se2 = mc_power_paired_dz(n2_std, dz, alpha, sides, reps, seed0 + 201) mc2c, se2c = mc_power_paired_dz(s2["current_matrices"], dz, alpha, sides, reps, seed0 + 202) res["stage2"] = dict(n_std=n2_std, n_dec=n2_dec) out.append( f"| Stage 2, generative | paired t on SI differences | dz = {dz} " f"(worked value -- the paper fixes no SI effect size) | {n2_std} pairs | " f"{fmt(mc2)} +/- {fmt(se2)} | {n2_dec} pairs | " f"{s2['current_matrices']} pairs | {fmt(mc2c)} +/- {fmt(se2c)} |") # Density calibration arm: one-sample vs chance. da = cfg["density_arm"] p1d = rate_from_h(da["null_preference_rate"], da["planning_effect_h_vs_null"]) nd_cf = n_one_proportion(da["null_preference_rate"], p1d, alpha, sides, power_t) nd_std = exact_search_one_sample(da["null_preference_rate"], p1d, alpha, sides, power_t, nd_cf) nd_cf_dec = n_one_proportion(da["null_preference_rate"], p1d, DECISIVE_ALPHA, sides, DECISIVE_POWER) nd_dec = exact_search_one_sample(da["null_preference_rate"], p1d, DECISIVE_ALPHA, sides, DECISIVE_POWER, nd_cf_dec) exd = exact_power_one_sample(nd_std, da["null_preference_rate"], p1d, alpha, sides) exdd = exact_power_one_sample(nd_dec, da["null_preference_rate"], p1d, DECISIVE_ALPHA, sides) exdc = exact_power_one_sample(20, da["null_preference_rate"], p1d, alpha, sides) res["density"] = dict(n_std=nd_std, n_dec=nd_dec) out.append( f"| Density calibration arm | one-sample exact binomial vs chance | " f"h = {da['planning_effect_h_vs_null']} -> rate " f"{da['null_preference_rate']} vs {fmt(p1d)} (exact search on " f"the unrounded rate {fmt(p1d, 5)}) | {nd_std} pairs " f"(exact power {fmt(exd)}; closed form {nd_cf}) | -- | {nd_dec} pairs " f"(exact power {fmt(exdd)}; closed form {nd_cf_dec}) | 20 pairs " f"(exact power {fmt(exdc)}) |") # Reader validation: one-sample convergence vs chance. rv = cfg["reader_validation"] nr_cf = n_one_proportion(rv["null_convergence_rate"], rv["planning_target_rate"], alpha, sides, power_t) nr_std = exact_search_one_sample(rv["null_convergence_rate"], rv["planning_target_rate"], alpha, sides, power_t, nr_cf) nr_cf_dec = n_one_proportion(rv["null_convergence_rate"], rv["planning_target_rate"], DECISIVE_ALPHA, sides, DECISIVE_POWER) nr_dec = exact_search_one_sample(rv["null_convergence_rate"], rv["planning_target_rate"], DECISIVE_ALPHA, sides, DECISIVE_POWER, nr_cf_dec) exr = exact_power_one_sample(nr_std, rv["null_convergence_rate"], rv["planning_target_rate"], alpha, sides) exrd = exact_power_one_sample(nr_dec, rv["null_convergence_rate"], rv["planning_target_rate"], DECISIVE_ALPHA, sides) exrc = exact_power_one_sample(20, rv["null_convergence_rate"], rv["planning_target_rate"], alpha, sides) res["reader"] = dict(n_std=nr_std, n_dec=nr_dec) out.append( f"| SI reader validation | one-sample exact binomial vs chance | " f"rate {rv['null_convergence_rate']} vs " f"{rv['planning_target_rate']} (worked values -- chance depends on the " f"coding scheme) | {nr_std} items (exact power {fmt(exr)}; closed " f"form {nr_cf}) | -- | {nr_dec} items (exact power {fmt(exrd)}; " f"closed form {nr_cf_dec}) | 20 items " f"(exact power {fmt(exrc)}) |") # Register arm: two-proportion across registers, per stratum. ra = cfg["register_arm"] p_alt_r = rate_from_h(ra["human_baseline_rate"], ra["planning_effect_h_between_registers"]) nreg_std = n_two_proportions(ra["planning_effect_h_between_registers"], alpha, sides, power_t) nreg_dec = n_two_proportions(ra["planning_effect_h_between_registers"], DECISIVE_ALPHA, sides, DECISIVE_POWER) mcr, ser = mc_power_two_prop(nreg_std, p_alt_r, ra["human_baseline_rate"], alpha, reps, seed0 + 301) mcrc, serc = mc_power_two_prop(20, p_alt_r, ra["human_baseline_rate"], alpha, reps, seed0 + 302) res["register"] = dict(n_std=nreg_std, n_dec=nreg_dec) out.append( f"| Register arm | two-proportion z per register stratum | " f"h = {ra['planning_effect_h_between_registers']} -> rate " f"{ra['human_baseline_rate']} vs {fmt(p_alt_r)} | {nreg_std} contents " f"per stratum ({nreg_std * 4} versions) | {fmt(mcr)} +/- {fmt(ser)} | " f"{nreg_dec} contents per stratum | 20 contents per stratum | " f"{fmt(mcrc)} +/- {fmt(serc)} |") out.append("") # ---- Section 2: falsifier-fire audit, core arm at registered size -- approx = ND.cdf(h_gov * math.sqrt(n_cur / 2) - z_alpha(alpha, sides)) mc_f, mc_f_se = mc_power_two_prop(n_cur, p_alt, p_base, alpha, reps, seed0 + 999) ex_z = exact_power_two_prop_z(n_cur, p_alt, p_base, alpha) fisher = exact_power_fisher(n_cur, p_alt, p_base, alpha) one_s = exact_power_one_sample(n_cur, p_base, p_or, alpha, sides) paired_mc, _ = mc_power_paired_dz(n_cur, s3["registered_thresholds"]["cohens_h"], alpha, sides, reps, seed0 + 998) out.append("## 2. Falsifier-fire audit -- Stage 3 at the registered size") out.append("") out.append(f"Setup: {n_cur} pairs per group; true effect set exactly at the " f"registered threshold (h = {fmt(h_gov)}, rates {p_base} vs " f"{fmt(p_alt)}). 'Fires' = the test is not significant, so the " f"registered falsifier triggers against a true effect.") out.append("") out.append("| Reading | Test | Power | P(falsifier fires) |") out.append("|---|---|---|---|") out.append(f"| Registered (model rate vs human rate) | normal approx | " f"{fmt(approx)} | {fmt(1 - approx)} |") out.append(f"| Registered | pooled z, Monte Carlo | {fmt(mc_f)} +/- " f"{fmt(mc_f_se)} | {fmt(1 - mc_f)} |") out.append(f"| Registered | pooled z, exact enumeration | {fmt(ex_z)} | " f"{fmt(1 - ex_z)} |") out.append(f"| Registered | Fisher exact, exact enumeration | " f"{fmt(fisher)} | {fmt(1 - fisher)} |") out.append(f"| Alternate (model rate vs chance, OR threshold) | exact " f"binomial | {fmt(one_s)} | {fmt(1 - one_s)} |") out.append(f"| Alternate (paired per-item contrast) | paired t, Monte " f"Carlo | {fmt(paired_mc)} | {fmt(1 - paired_mc)} |") out.append("") fires = [1 - approx, 1 - mc_f, 1 - ex_z, 1 - fisher] out.append(f"Registered reading, across test variants: P(fire) = " f"{fmt(min(fires))} to {fmt(max(fires))}.") out.append("") # ---- Section 3: gate budget at both design points ------------------ gate = cfg["gate"] for label, target in (("standard", n_std), ("decisive", n_dec)): out.append(f"## 3{label == 'decisive' and 'b' or 'a'}. " f"{gate['name']} -- build budget ({label} point)") out.append("") out.append(f"Target: >= {target} pairs pass the gate. Assurance: " f"P(passing >= {target}) >= {gate['assurance']}. " f"Exact binomial search; last column re-checks by Monte Carlo.") out.append("") out.append("| Pair failure rate | Pairs to build (exact) | " "Expected passing | MC assurance |") out.append("|---|---|---|---|") for j, f in enumerate(gate["failure_rate_scenarios"]): B = gate_build_count(target, f, gate["assurance"]) mc_a, mc_a_se = mc_gate_check(B, target, f, reps, seed0 + 5000 + (0 if label == "standard" else 100) + j) out.append(f"| {f} | {B} | {fmt(B * (1 - f), 1)} | " f"{fmt(mc_a)} +/- {fmt(mc_a_se)} |") out.append("") # ---- Section 4: scale table ---------------------------------------- sc = cfg["scale"] judges = sc["model_families"] + sc["human_raters_per_pair"] out.append("## 4. Scale table") out.append("") out.append(f"Per pair: {sc['versions_per_pair']} stimulus versions " f"({sc['versions_per_pair'] * 2} with the register arm), " f"{judges} Stage 3 judge calls ({sc['model_families']} model " f"families + {sc['human_raters_per_pair']} human rater), " f"SI coding passes = versions x 2 raters. Stage 4 = two further " f"Stage 3 passes (length-mismatched, framing-varied).") out.append("") out.append("| Pairs | Versions | Versions with register arm | " "Stage 3 judge calls | Stage 4 judge calls | SI coding passes |") out.append("|---|---|---|---|---|---|") for P in sc["pair_counts"]: v = P * sc["versions_per_pair"] calls = P * judges out.append(f"| {P} | {v} | {v * 2} | {calls} | {calls * 2} | {v * 2} |") out.append("") # ---- Section 5: registration skeleton ------------------------------- out.append("## 5. Registration skeleton (filled from this config)") out.append("") out.append("1. Construct: summarization bias -- the systematic tendency of a " "model, in generation or in evaluation, to replace or prefer " "reconstructable shown-mode structure with the abstract summary " "label that names its content (paper S4.3).") out.append(f"2. Thresholds, fixed in advance: model told-preference rate " f"exceeds human rate with Cohen's h >= " f"{s3['registered_thresholds']['cohens_h']} or odds ratio >= " f"{s3['registered_thresholds']['odds_ratio']} (paper S5.3). " f"Alpha = {alpha} ({sides}-sided) at the standard point; " f"alpha = {DECISIVE_ALPHA}, power = {DECISIVE_POWER} at the " f"decisive point.") out.append(f"3. Stimuli: build the gate-budget count from Table 3 for the " f"expected failure rate (standard point) -- word count +/-5%, " f"grade matched (paper S5.1).") out.append("4. Coding gate: SI coded by two independent raters per version; " "no category proceeds at Cohen's kappa < 0.60 (paper S5.1).") out.append(f"5. Input-separation gate: a pair enters Stage 3 only if its two " f"members separate on SI by the registered minimum; failed pairs " f"are counted, reported, and replaced from the build budget " f"(thread post #5).") out.append(f"6. Stage 2 (generative): >= {sc['model_families']} model " f"families, shown-mode instruction, SI compared to the human " f"shown target; n = {n2_std} pairs at the standard point, " f"{n2_dec} at the decisive point (Table 1).") out.append(f"7. Stage 3 (evaluative, core): pairs judged by the same models " f"and by a human panel, order counterbalanced; n = {n_std} pairs " f"per group at the standard point, {n_dec} at the decisive " f"point (Table 1).") out.append("8. Stage 4 (discriminant): Stage 3 re-run (a) length-mismatched " "and (b) framing-varied; a preference that vanishes under either " "is reclassified as that known bias (paper S5.4).") out.append(f"9. Added arms, sized in Table 1: density calibration " f"({nd_std} pairs), SI reader validation ({nr_std} items), " f"register arm ({nreg_std} contents per stratum). No arm is " f"registered without its own size from this table.") out.append("10. Decision rules (paper S5.5): Stage 3 supports and Stage 4 " "survives -> confirmed at the scope tested. Stage 3 falsifies " "-> construct withdrawn, published with the same prominence as " "a confirmation. Post-hoc qualifiers are prohibited.") out.append("") return "\n".join(out) + "\n" KNOWN_SCENARIOS = ("null", "registered_alternative") def check_scenario(cfg): """Refuse an unknown truth.scenario by name, on any config that carries one -- an unknown string must never silently select the world a verdict is about. Configs without a truth block (the sizing configs this program normally reads) are unaffected.""" truth = cfg.get("truth") if isinstance(truth, dict) and "scenario" in truth: if truth["scenario"] not in KNOWN_SCENARIOS: sys.exit(f"config refused: truth.scenario " f"{truth['scenario']!r} is not a known scenario " f"(expected 'null' or 'registered_alternative')") def main(): if len(sys.argv) != 2: sys.exit("usage: python3 protocol_engine.py ") with open(sys.argv[1], "r", encoding="utf-8") as fh: cfg = json.load(fh) check_scenario(cfg) if "versions_per_pair" not in cfg.get("scale", {}): cfg["scale"]["versions_per_pair"] = 2 if "mc_reps" not in cfg: cfg["mc_reps"] = 20000 sys.stdout.write(run(cfg)) if __name__ == "__main__": main()