# Independent re-derivation (reviewer-written): forward DP over success counts with exact
# Fraction probabilities and p-values computed via scipy-free cumulative sums; compares every
# row of the bundle's results/R1.json.
import json, sys
from fractions import Fraction
from math import comb
R = json.load(open(sys.argv[1]))
def pval(n, k):  # exact one-sided upper tail under p=1/2
    return Fraction(sum(comb(n, j) for j in range(k, n + 1)), 2 ** n)
def rej(rule, n, k, H):
    if rule == 'likelihood_ratio':
        return Fraction(3 ** k, 2 ** n) >= 20
    if rule == 'fixed_binomial' and n != H:
        return False
    return pval(n, k) <= Fraction(1, 20)
def calc(rule, H, p):
    alive = {0: Fraction(1)}; P = Fraction(0); E = Fraction(0)
    for n in range(1, H + 1):
        nxt = {}
        for k, w in alive.items():
            nxt[k] = nxt.get(k, 0) + w * (1 - p)
            nxt[k + 1] = nxt.get(k + 1, 0) + w * p
        alive = {}
        for k, w in nxt.items():
            if rej(rule, n, k, H):
                P += w; E += n * w
            else:
                alive[k] = w
    E += H * sum(alive.values())
    return P, E
bad = 0
for r in R['rows']:
    a, d = map(int, r['p'].split('/'))
    P, E = calc(r['rule'], r['horizon'], Fraction(a, d))
    ok = (f'{P.numerator}/{P.denominator}' == r['probability_exact'] and float(P) == r['rejection_probability']
          and float(E) == r['expected_observations'])
    bad += not ok
    if not ok: print('MISMATCH', r)
print(json.dumps({'rows_checked': len(R['rows']), 'mismatches': bad,
                  'null_h100': R['null_horizon_100'], 'oracle_checks_passed': R['oracle_checks_passed']}))