calibration_report.py

#!/usr/bin/env python3
"""Calibration report — aggregate analysis over results/games/*.json (the historical tracker).

This is the M5 calibration loop: after each settled game is added, re-run to see where the model
beats / trails the crowd by category, the running RBP, and the biggest hits/misses — the signal
that drives model tuning. Stdlib-only.

Usage:
    python scripts/calibration_report.py            # full report over all settled games
    python scripts/calibration_report.py --json      # machine-readable summary
"""

from __future__ import annotations

import glob
import json
import sys
from pathlib import Path

REPO_ROOT = Path(__file__).resolve().parent.parent
RESULTS_DIR = REPO_ROOT / "results" / "games"
LEADERBOARD = REPO_ROOT / "results" / "leaderboard.json"


def load_games() -> list[dict]:
    games = []
    for f in sorted(glob.glob(str(RESULTS_DIR / "*.json"))):
        try:
            d = json.load(open(f, encoding="utf-8"))
            if d.get("status") == "settled":
                games.append(d)
        except (OSError, json.JSONDecodeError):
            continue
    return games


def _brier(p: float, o: float) -> float:
    return ((p - o) / 100.0) ** 2


def murphy_decomposition(rows: list[dict], n_bins: int = 3) -> dict:
    """Murphy 3-component Brier decomposition: BS = reliability - resolution + uncertainty.

    Uses few bins (default 3) because the settled sample is small — fine-grained bins would be
    noise. Reliability low = well-calibrated; resolution high = forecasts separate outcomes.
    """
    pts = [(r["you"] / 100.0, r["outcome"] / 100.0) for r in rows]
    if not pts:
        return {}
    n = len(pts)
    obar = sum(o for _, o in pts) / n
    bins: list[list] = [[] for _ in range(n_bins)]
    for p, o in pts:
        idx = min(n_bins - 1, int(p * n_bins))
        bins[idx].append((p, o))
    reliability = resolution = 0.0
    for b in bins:
        if not b:
            continue
        nb = len(b)
        fbar = sum(p for p, _ in b) / nb
        obar_k = sum(o for _, o in b) / nb
        reliability += nb * (fbar - obar_k) ** 2
        resolution += nb * (obar_k - obar) ** 2
    return {
        "reliability": round(reliability / n, 4),   # lower better (calibration error)
        "resolution": round(resolution / n, 4),     # higher better (discrimination)
        "uncertainty": round(obar * (1 - obar), 4),  # intrinsic
        "n_bins": n_bins,
    }


def temperature_estimate(rows: list[dict]) -> dict:
    """Best single temperature T (logit(p)/T) minimizing Brier on settled data.

    T<1 => model underconfident (should extremize); T>1 => overconfident (should shrink).
    Reported as a DIAGNOSTIC only — do not apply until the sample is large (>=50) to avoid overfit.
    """
    import math
    pts = [(max(0.01, min(0.99, r["you"] / 100.0)), r["outcome"] / 100.0) for r in rows]
    if len(pts) < 5:
        return {"T": None, "note": "need >=5 points"}

    def brier_at(T: float) -> float:
        s = 0.0
        for p, o in pts:
            lp = math.log(p / (1 - p)) / T
            pc = 1 / (1 + math.exp(-lp))
            s += (pc - o) ** 2
        return s / len(pts)

    best_T, best_b = 1.0, brier_at(1.0)
    T = 0.5
    while T <= 2.5:
        b = brier_at(T)
        if b < best_b:
            best_T, best_b = T, b
        T += 0.05
    signal = "extremize (underconfident)" if best_T < 0.95 else (
        "shrink (overconfident)" if best_T > 1.05 else "well-calibrated")
    return {"T": round(best_T, 2), "brier_at_T": round(best_b, 4), "signal": signal,
            "note": "diagnostic only; apply T globally once >=50 settled predictions"}


def summarize(games: list[dict]) -> dict:
    rows = []
    for g in games:
        for q in g.get("questions", []):
            # Skip voided/unscored markets (e.g. a player prop removed when the player isn't in the
            # final squad): these carry null you/outcome and were never scored by the platform.
            if q.get("you") is None or q.get("outcome") is None:
                continue
            rows.append({**q, "game": g.get("game"), "mult": g.get("stage_multiplier", 1)})

    n = len(rows)
    by_cat: dict[str, list] = {}
    beat = below = neutral = 0
    no_rbp = yes_rbp = 0.0
    no_n = yes_n = 0
    you_brier = []
    crowd_brier = []
    for r in rows:
        by_cat.setdefault(r["category"], []).append(r)
        if r.get("beat_crowd") == "beat":
            beat += 1
        elif r.get("beat_crowd") == "below":
            below += 1
        elif r.get("beat_crowd") == "neutral":
            neutral += 1
        rbp = r.get("rbp")
        if isinstance(rbp, (int, float)):
            if r["outcome"] == 0:
                no_rbp += rbp
                no_n += 1
            else:
                yes_rbp += rbp
                yes_n += 1
        if isinstance(r.get("crowd"), (int, float)):
            you_brier.append(_brier(r["you"], r["outcome"]))
            crowd_brier.append(_brier(r["crowd"], r["outcome"]))

    cat_summary = {}
    for c, rs in by_cat.items():
        rbps = [x["rbp"] for x in rs if isinstance(x.get("rbp"), (int, float))]
        gaps = [x["crowd"] - x["you"] for x in rs if isinstance(x.get("crowd"), (int, float))]
        cat_summary[c] = {
            "n": len(rs),
            "total_rbp": round(sum(rbps), 2) if rbps else None,
            "avg_rbp": round(sum(rbps) / len(rbps), 2) if rbps else None,
            "avg_crowd_minus_you": round(sum(gaps) / len(gaps), 1) if gaps else None,
        }

    total_rbp = sum(x["rbp"] for x in rows if isinstance(x.get("rbp"), (int, float)))
    return {
        "games": len(games),
        "questions": n,
        "total_rbp": round(total_rbp, 2),
        "avg_rbp_per_q": round(total_rbp / n, 2) if n else None,
        "beat_crowd": beat,
        "below_crowd": below,
        "neutral": neutral,
        "you_mean_brier": round(sum(you_brier) / len(you_brier), 4) if you_brier else None,
        "crowd_mean_brier": round(sum(crowd_brier) / len(crowd_brier), 4) if crowd_brier else None,
        "no_avg_rbp": round(no_rbp / no_n, 2) if no_n else None,
        "yes_avg_rbp": round(yes_rbp / yes_n, 2) if yes_n else None,
        "by_category": dict(sorted(cat_summary.items(),
                                   key=lambda kv: (kv[1]["total_rbp"] is None, -(kv[1]["total_rbp"] or 0)))),
        "brier_decomposition": murphy_decomposition(rows),
        "temperature": temperature_estimate(rows),
        "rows": rows,
    }


def print_report(s: dict) -> None:
    print(f"== Calibration report: {s['games']} games, {s['questions']} questions ==")
    print(f"Total RBP: {s['total_rbp']:+.2f}  |  avg/question: {s['avg_rbp_per_q']:+.2f}")
    print(f"Beat crowd: {s['beat_crowd']}  |  below: {s['below_crowd']}  |  neutral: {s['neutral']}")
    if s["you_mean_brier"] is not None:
        edge = s["crowd_mean_brier"] - s["you_mean_brier"]
        print(f"Mean Brier — you {s['you_mean_brier']} vs crowd {s['crowd_mean_brier']}  (you {'better' if edge>0 else 'worse'} by {abs(edge):.4f})")
    print(f"NO-outcome avg RBP: {s['no_avg_rbp']}  |  YES-outcome avg RBP: {s['yes_avg_rbp']}")
    print("\nBy category (sorted by total RBP):")
    print(f"  {'category':14} {'n':>2}  {'total':>7}  {'avg':>6}  crowd-you")
    for c, v in s["by_category"].items():
        print(f"  {c:14} {v['n']:>2}  {str(v['total_rbp']):>7}  {str(v['avg_rbp']):>6}  {v['avg_crowd_minus_you']}")

    rows = [r for r in s["rows"] if isinstance(r.get("rbp"), (int, float))]
    rows.sort(key=lambda r: r["rbp"])
    print("\nBiggest misses:")
    for r in rows[:4]:
        print(f"  {r['rbp']:+6.2f}  you={r['you']} crowd={r.get('crowd')} out={r['outcome']}  {r['question'][:50]}")
    print("Biggest wins:")
    for r in rows[-4:][::-1]:
        print(f"  {r['rbp']:+6.2f}  you={r['you']} crowd={r.get('crowd')} out={r['outcome']}  {r['question'][:50]}")

    bd = s.get("brier_decomposition") or {}
    if bd:
        print(f"\nBrier decomposition ({bd.get('n_bins')} bins): reliability {bd.get('reliability')} "
              f"(lower=better calib) | resolution {bd.get('resolution')} (higher=better) | "
              f"uncertainty {bd.get('uncertainty')}")
    t = s.get("temperature") or {}
    if t.get("T") is not None:
        print(f"Temperature estimate: T={t['T']} -> {t.get('signal')}  [{t.get('note')}]")

    if LEADERBOARD.exists():
        lb = json.load(open(LEADERBOARD, encoding="utf-8"))
        print("\nLeaderboard history:")
        for h in lb.get("history", []):
            print(f"  after game {h['after_game']}: {h['cumulative_points']} pts, "
                  f"position {h['position']}/{h.get('field_size','?')} (game RBP {h.get('game_rbp')})")


def main() -> None:
    games = load_games()
    if not games:
        print("No settled games in results/games/. Add one to start tracking.")
        return
    s = summarize(games)
    if "--json" in sys.argv:
        s.pop("rows", None)
        print(json.dumps(s, indent=2))
    else:
        print_report(s)


if __name__ == "__main__":
    main()