Coverage for orchestrant/benchmark/stats.py: 99%
227 statements
« prev ^ index » next coverage.py v7.16.2, created at 2026-10-18 00:48 +0000
« prev ^ index » next coverage.py v7.16.2, created at 2026-10-18 00:48 +0000
1"""Small statistics helpers, so scores are not published as bare fractions."""
3import math
4from fractions import Fraction
7# Two-sided p-value below which a paired difference counts as separable.
8ALPHA = 0.05
11def wilson_interval(successes, trials, z=1.96):
12 """95 % confidence interval for a proportion. Returns (low, high) in 0..1."""
13 if trials <= 0:
14 return (0.0, 1.0)
15 if successes < 0 or successes > trials:
16 raise ValueError(f"successes {successes} out of range for {trials} trials")
17 p = successes / trials
18 denom = 1 + z * z / trials
19 centre = (p + z * z / (2 * trials)) / denom
20 half = z * math.sqrt(p * (1 - p) / trials + z * z / (4 * trials * trials)) / denom
21 return (max(0.0, centre - half), min(1.0, centre + half))
24def intervals_overlap(a_succ, a_tot, b_succ, b_tot, z=1.96):
25 """Do two scores' intervals overlap? If so, they are not separable."""
26 a_lo, a_hi = wilson_interval(a_succ, a_tot, z)
27 b_lo, b_hi = wilson_interval(b_succ, b_tot, z)
28 return a_lo <= b_hi and b_lo <= a_hi
31def format_score(successes, trials, width=None):
32 """'8/12 = 67% [39-86%]' — the score with what it can actually support."""
33 if trials <= 0:
34 return "n/a"
35 lo, hi = wilson_interval(successes, trials)
36 pct = 100 * successes / trials
37 s = f"{successes}/{trials} = {pct:.0f}% [{100 * lo:.0f}-{100 * hi:.0f}%]"
38 return f"{s:{width}}" if width else s
41# Floor for the back-flip rate: seeing none in a few dozen cases does not make it zero.
42DEFAULT_BACK_FLIP_RATE = 0.05
45def _counts(value):
46 """(passes, attempts) from a (passes, attempts) pair or a one-draw bool."""
47 if isinstance(value, bool):
48 return (int(value), 1)
49 passes, attempts = value
50 if not 0 <= passes <= max(attempts, 0):
51 raise ValueError(f"passes {passes} out of range for {attempts} attempts")
52 return passes, attempts
55def clustered_rate(cases, z=1.96):
56 """Pooled pass rate with the CASE, not the draw, as the unit of sampling.
58 Returns None when nothing was attempted, else a dict: rate, se, design_effect,
59 n_eff, low, high, n_cases, passes, attempts.
60 """
61 counts = [c for c in map(_counts, cases.values()) if c[1] > 0]
62 if not counts:
63 return None
64 passes = sum(p for p, _ in counts)
65 attempts = sum(m for _, m in counts)
66 n_cases = len(counts)
67 rate = passes / attempts
68 naive = rate * (1 - rate) / attempts
69 var = None
70 if n_cases > 1:
71 resid = sum((p - rate * m) ** 2 for p, m in counts)
72 # CR1 cluster-robust variance of a ratio estimator (Miller 2024, "Adding Error Bars to Evals").
73 var = n_cases / (n_cases - 1) * resid / attempts**2
74 if var is None or naive == 0:
75 deff = None
76 n_eff = attempts**2 / sum(m * m for _, m in counts)
77 else:
78 deff = var / naive
79 # Floored at 1: repeats are not anti-correlated, so a lower ratio is noise.
80 n_eff = attempts / max(1.0, deff)
81 low, high = wilson_interval(rate * n_eff, n_eff, z)
82 return {
83 "rate": rate,
84 "se": None if var is None else math.sqrt(var),
85 "design_effect": deff,
86 "n_eff": n_eff,
87 "low": low,
88 "high": high,
89 "n_cases": n_cases,
90 "passes": passes,
91 "attempts": attempts,
92 }
95def paired_difference(a_cases, b_cases, z=1.96):
96 """Mean per-case change in pass rate, B minus A, with a case-clustered interval.
98 Returns None when no case is shared, else (mean, low, high, n_cases).
99 """
100 diffs = []
101 for key in a_cases:
102 if key not in b_cases:
103 continue
104 (pa, ma), (pb, mb) = _counts(a_cases[key]), _counts(b_cases[key])
105 if ma > 0 and mb > 0:
106 # Exact: float thirds that cancel printed "paired diff -0pt" for no change.
107 diffs.append(Fraction(pb, mb) - Fraction(pa, ma))
108 if not diffs:
109 return None
110 n = len(diffs)
111 mean = sum(diffs, Fraction(0)) / n
112 if n < 2:
113 return (float(mean), -1.0, 1.0, n)
114 spread = sum(((d - mean) ** 2 for d in diffs), Fraction(0)) / (n * (n - 1))
115 m, se = float(mean), math.sqrt(float(spread))
116 return (m, max(-1.0, m - z * se), min(1.0, m + z * se), n)
119def pass_hat_k(cases, k):
120 """Unbiased chance a case passes k draws out of k (tau-bench's pass^k), or None."""
121 if k < 1:
122 raise ValueError(f"k must be at least 1, got {k}")
123 terms = [
124 math.comb(p, k) / math.comb(m, k)
125 for p, m in map(_counts, cases.values())
126 if m >= k
127 ]
128 return sum(terms) / len(terms) if terms else None
131def _sign_test_rejects(n_cases, alpha):
132 """Per discordant count d, the most `better` cases still read as a regression, or -1."""
133 limits = []
134 for d in range(n_cases + 1):
135 limit, tail, coef, k = -1, 0, 1, 0
136 while 2 * k < d:
137 tail += coef
138 if min(1.0, 2 * (tail / 2**d)) >= alpha:
139 break
140 limit = k
141 coef = coef * (d - k) // (k + 1)
142 k += 1
143 limits.append(limit)
144 return limits
147def _binomial_pmf(n, p):
148 """P(X = k) for k = 0..n, X ~ Binomial(n, p), in log space against underflow."""
149 if p <= 0:
150 return [1.0] + [0.0] * n
151 if p >= 1:
152 return [0.0] * n + [1.0]
153 lp, lq, lg = math.log(p), math.log1p(-p), math.lgamma
154 return [
155 math.exp(lg(n + 1) - lg(k + 1) - lg(n - k + 1) + k * lp + (n - k) * lq)
156 for k in range(n + 1)
157 ]
160def paired_power(
161 n_cases, drop, back_flip_rate=DEFAULT_BACK_FLIP_RATE, alpha=ALPHA, _limits=None
162):
163 """Exact chance paired_sign_test flags `drop` (cases worsen at back_flip_rate + drop)."""
164 p_worse, p_better = back_flip_rate + drop, back_flip_rate
165 if drop < 0 or p_better < 0 or p_worse + p_better > 1:
166 raise ValueError(f"no such case mix: drop {drop}, back-flip {back_flip_rate}")
167 discordant = p_worse + p_better
168 if n_cases <= 0 or discordant == 0:
169 return 0.0
170 limits = _limits or _sign_test_rejects(n_cases, alpha)
171 power = 0.0
172 for d, p_d in enumerate(_binomial_pmf(n_cases, discordant)):
173 if limits[d] >= 0 and p_d >= 1e-15:
174 better = _binomial_pmf(d, p_better / discordant)
175 power += p_d * sum(better[: limits[d] + 1])
176 return power
179def paired_mde(n_cases, back_flip_rate=DEFAULT_BACK_FLIP_RATE, power=0.8, alpha=ALPHA):
180 """Smallest net drop, as a fraction of the cases, caught with `power`, or None."""
181 if n_cases <= 0:
182 return None
183 limits = _sign_test_rejects(n_cases, alpha)
184 low, high = 0.0, 1.0 - 2 * back_flip_rate
185 if high <= 0 or paired_power(n_cases, high, back_flip_rate, alpha, limits) < power:
186 return None
187 while high - low > 1e-4:
188 mid = (low + high) / 2
189 if paired_power(n_cases, mid, back_flip_rate, alpha, limits) >= power:
190 high = mid
191 else:
192 low = mid
193 return high
196def back_flip_estimate(n_cases, back_flips=0):
197 """(rate, source) for paired_mde: back_flips / n_cases, floored at DEFAULT_BACK_FLIP_RATE."""
198 observed = back_flips / n_cases if n_cases > 0 else 0.0
199 if observed >= DEFAULT_BACK_FLIP_RATE:
200 return observed, "observed"
201 return DEFAULT_BACK_FLIP_RATE, f"assumed, {back_flips or 'none'} observed"
204def paired_mde_note(n_cases, back_flips=0, power=0.8, alpha=ALPHA):
205 """One line: the smallest drop a paired comparison of n_cases would catch."""
206 rate, source = back_flip_estimate(n_cases, back_flips)
207 mde = paired_mde(n_cases, rate, power, alpha)
208 if mde is None:
209 return (
210 f"at {n_cases} paired cases no drop is caught with {power:.0%} power "
211 f"(back-flip rate {rate:.0%}, {source}) — 'no regression' here "
212 f"means 'cannot tell'"
213 )
214 return (
215 f"minimum detectable drop at {power:.0%} power: {100 * mde:.0f}pt "
216 f"(~{mde * n_cases:.0f} of {n_cases} paired cases; back-flip rate "
217 f"{rate:.0%}, {source}) — a smaller real drop is missed more than "
218 f"{1 - power:.0%} of the time"
219 )
222def clustered_note(cases):
223 """' clustered [66-88%, deff 2.6]' when some case's repeats disagree, else ''."""
224 if not any(m > 1 and 0 < p < m for p, m in map(_counts, cases.values())):
225 return ""
226 c = clustered_rate(cases)
227 deff = "" if c["design_effect"] is None else f", deff {c['design_effect']:.1f}"
228 return f" clustered [{100 * c['low']:.0f}-{100 * c['high']:.0f}%{deff}]"
231def paired_diff_note(a_cases, b_cases):
232 """'paired diff -22pt [-38, -6]' over the shared cases, or None if none are shared."""
233 paired = paired_difference(a_cases, b_cases)
234 if paired is None:
235 return None
236 mean, low, high, _ = paired
237 return f"paired diff {100 * mean:+.0f}pt [{100 * low:+.0f}, {100 * high:+.0f}]"
240def pass_k_note(a_cases, b_cases, k):
241 """'pass^3 73% -> 73% (...)': an agent that retries nothing needs every draw to pass."""
242 rates = (pass_hat_k(cases, k) for cases in (a_cases, b_cases))
243 shown = " -> ".join("n/a" if v is None else f"{v:.0%}" for v in rates)
244 return f"pass^{k} {shown} (a case counts only when all {k} of its draws pass)"
247def significance_note(a_label, a_succ, a_tot, b_label, b_succ, b_tot):
248 """Plain-language verdict on whether two scores are separable at all."""
249 if a_tot <= 0 or b_tot <= 0:
250 return "one side has no observations — nothing to compare"
251 a_rate, b_rate = a_succ / a_tot, b_succ / b_tot
252 if intervals_overlap(a_succ, a_tot, b_succ, b_tot):
253 better = a_label if a_rate > b_rate else b_label
254 if a_rate == b_rate:
255 return "identical rates"
256 return (
257 f"{better} scores higher, but the 95 % intervals OVERLAP at this "
258 f"sample size — not separable on this evidence alone"
259 )
260 better = a_label if a_rate > b_rate else b_label
261 return f"{better} is higher and the intervals do not overlap"
264def smallest_separable_rate(trials, from_rate=1.0, z=1.96):
265 """The lowest rate still distinguishable from `from_rate` here, or None if none is."""
266 if trials <= 0:
267 return None
268 successes = round(from_rate * trials)
269 for lower in range(successes - 1, -1, -1):
270 if not intervals_overlap(successes, trials, lower, trials, z):
271 return lower / trials
272 return None
275def power_note(trials, from_rate=1.0):
276 """One line stating what a 'no regression' verdict is actually worth."""
277 mde = smallest_separable_rate(trials, from_rate)
278 if mde is None:
279 return (
280 f"at n={trials} this suite cannot prove ANY drop — "
281 f"'no regression' here means 'cannot tell'"
282 )
283 return (
284 f"at n={trials} the smallest provable drop is {100 * from_rate:.0f}% -> "
285 f"{100 * mde:.0f}%; anything subtler passes unnoticed"
286 )
289def diff_interval(a_succ, a_tot, b_succ, b_tot, z=1.96):
290 """Newcombe 95 % interval for b_rate - a_rate; (-1.0, 1.0) when a side has no trials."""
291 if a_tot <= 0 or b_tot <= 0:
292 return (-1.0, 1.0)
293 pa, pb = a_succ / a_tot, b_succ / b_tot
294 la, ua = wilson_interval(a_succ, a_tot, z)
295 lb, ub = wilson_interval(b_succ, b_tot, z)
296 d = pb - pa
297 low = d - math.sqrt((pb - lb) ** 2 + (ua - pa) ** 2)
298 high = d + math.sqrt((ub - pb) ** 2 + (pa - la) ** 2)
299 return (max(-1.0, low), min(1.0, high))
302def paired_sign_test(discordant_a, discordant_b):
303 """Exact two-sided sign test p-value over the cases only A, or only B, got right."""
304 if discordant_a < 0 or discordant_b < 0:
305 raise ValueError("discordant counts cannot be negative")
306 n = discordant_a + discordant_b
307 if n == 0:
308 return 1.0
309 k = min(discordant_a, discordant_b)
310 tail = sum(math.comb(n, i) for i in range(k + 1)) / 2**n
311 return min(1.0, 2 * tail)
314def _rate(value):
315 if isinstance(value, bool):
316 return 1.0 if value else 0.0
317 passes, attempts = value
318 return passes / attempts if attempts else None
321def paired_outcomes(a_cases, b_cases):
322 """Count shared cases where A did better, B did better, or neither; unmeasured are skipped."""
323 a_better = b_better = ties = 0
324 for key in set(a_cases) & set(b_cases):
325 ra, rb = _rate(a_cases[key]), _rate(b_cases[key])
326 if ra is None or rb is None:
327 continue
328 if ra > rb:
329 a_better += 1
330 elif rb > ra:
331 b_better += 1
332 else:
333 ties += 1
334 return a_better, b_better, ties
337def smallest_detectable_flips(alpha=ALPHA):
338 """How many cases must flip ONE way, with none flipping back, to be seen."""
339 k = 1
340 while paired_sign_test(k, 0) >= alpha:
341 k += 1
342 return k
345def paired_power_note(alpha=ALPHA):
346 k = smallest_detectable_flips(alpha)
347 return (
348 f"paired: {k} cases flipping the same way (none flipping back) would "
349 f"be detected; fewer cannot be, whatever the suite size"
350 )
353def tiers(rows, key, alpha=ALPHA):
354 """Group ranked rows whose neighbours the paired sign test cannot separate.
356 `key(row)` returns the row's per-case outcomes; a new tier starts where p < alpha.
357 """
358 groups = []
359 for row in rows:
360 if groups:
361 a_better, b_better, _ = paired_outcomes(key(groups[-1][-1]), key(row))
362 if paired_sign_test(a_better, b_better) < alpha:
363 groups.append([row])
364 continue
365 groups[-1].append(row)
366 else:
367 groups.append([row])
368 return groups