cup26/scripts/buildBacktest.ts

// Walk-forward backtest of the Elo + Dixon-Coles model on historical
// internationals, with an honest train/test split (params fit on pre-2018,
// tested out-of-sample on 2018→now). Outputs proper scoring metrics (Brier,
// log-loss, RPS, accuracy), baseline comparisons, and a calibration/reliability
// analysis → public/data/backtest.json, shown on the Methodology page.
import { readFileSync, writeFileSync, mkdirSync } from 'node:fs';
import { fileURLToPath } from 'node:url';
import { dirname, join } from 'node:path';
import { canonicalTeam } from '../src/lib/teams';
import { HOME_ADV_ELO, importanceWeight, eloDelta } from '../src/lib/model/elo';
import { lambdasFromElo, scoreMatrix, outcomeProbs, poissonPmf, type ModelParams } from '../src/lib/model/poisson';

const ROOT = join(dirname(fileURLToPath(import.meta.url)), '..');
const CSV = join(ROOT, 'data', 'raw', 'international-results.csv');
const OUT_DIR = join(ROOT, 'public', 'data');

const START_ELO = 1500;
const PARAM_FROM = '2006-01-01';
const TRAIN_END = '2018-01-01'; // params fit on [PARAM_FROM, TRAIN_END)
const TEST_TO = '2026-06-01'; // exclude the 2026 World Cup itself

interface Row { date: string; home: string; away: string; hs: number; as: number; tournament: string; neutral: boolean }
type Probs = { h: number; d: number; a: number };
type Outcome = 'h' | 'd' | 'a';

function parse(): Row[] {
  const lines = readFileSync(CSV, 'utf8').split('\n');
  const rows: Row[] = [];
  for (let i = 1; i < lines.length; i++) {
    const c = lines[i]?.split(',');
    if (!c || c.length < 9) continue;
    const hs = Number(c[3]); const as = Number(c[4]);
    if (!Number.isFinite(hs) || !Number.isFinite(as)) continue;
    rows.push({ date: c[0]!, home: c[1]!, away: c[2]!, hs, as, tournament: c[5]!, neutral: c[8]!.trim().toUpperCase() === 'TRUE' });
  }
  rows.sort((a, b) => a.date.localeCompare(b.date));
  return rows;
}

function dcTau(h: number, a: number, lh: number, la: number, rho: number): number {
  if (h === 0 && a === 0) return 1 - lh * la * rho;
  if (h === 0 && a === 1) return 1 + lh * rho;
  if (h === 1 && a === 0) return 1 + la * rho;
  if (h === 1 && a === 1) return 1 - rho;
  return 1;
}

const outcomeOf = (hs: number, as: number): Outcome => (hs > as ? 'h' : hs < as ? 'a' : 'd');

// ---- scoring ----
function brier(p: Probs, o: Outcome): number {
  return (p.h - (o === 'h' ? 1 : 0)) ** 2 + (p.d - (o === 'd' ? 1 : 0)) ** 2 + (p.a - (o === 'a' ? 1 : 0)) ** 2;
}
function logloss(p: Probs, o: Outcome): number {
  return -Math.log(Math.max(1e-12, p[o]));
}
function rps(p: Probs, o: Outcome): number {
  // ordered H, D, A
  const y = { h: o === 'h' ? 1 : 0, d: o === 'd' ? 1 : 0, a: o === 'a' ? 1 : 0 };
  const c1p = p.h, c1y = y.h;
  const c2p = p.h + p.d, c2y = y.h + y.d;
  return ((c1p - c1y) ** 2 + (c2p - c2y) ** 2) / 2;
}
const argmax = (p: Probs): Outcome => (p.h >= p.d && p.h >= p.a ? 'h' : p.d >= p.a ? 'd' : 'a');

function scoreSet(preds: { p: Probs; o: Outcome }[]) {
  let b = 0, l = 0, r = 0, correct = 0;
  for (const { p, o } of preds) { b += brier(p, o); l += logloss(p, o); r += rps(p, o); if (argmax(p) === o) correct++; }
  const n = preds.length;
  return { brier: b / n, logloss: l / n, rps: r / n, accuracy: correct / n };
}

function main(): void {
  const rows = parse();
  const elo = new Map<string, number>();
  const get = (t: string) => elo.get(t) ?? START_ELO;

  // ---- fit params on the training window ----
  let sumTotal = 0, nCal = 0, sxy = 0, sxx = 0;
  const calib: { d: number; h: number; a: number }[] = [];
  // we need a first pass to fit params using pre-match Elo; do it inline while walking
  // (Elo is updated every match; calibration accumulates only within the train window)

  const testPreds: { p: Probs; o: Outcome }[] = [];
  const baseUniform: { p: Probs; o: Outcome }[] = [];
  const baseRate: { p: Probs; o: Outcome }[] = [];
  const baseElo: { p: Probs; o: Outcome }[] = [];
  const reliabilityRaw: { p: number; y: number }[] = [];
  let wcCorrect = 0, wcN = 0;

  // outcome base rates over the test window (computed in a quick pre-scan)
  let th = 0, td = 0, ta = 0, tn = 0;
  for (const r of rows) {
    if (r.date < TRAIN_END || r.date >= TEST_TO) continue;
    const o = outcomeOf(r.hs, r.as); if (o === 'h') th++; else if (o === 'd') td++; else ta++; tn++;
  }
  const rate: Probs = { h: th / tn, d: td / tn, a: ta / tn };

  // params get finalized at TRAIN_END; predictions before that use a rolling fit.
  let params: ModelParams = { goalsPerElo: 0.0057, avgGoals: 2.7, rho: -0.05, homeAdvElo: HOME_ADV_ELO };
  let paramsFixed = false;

  const finalizeParams = (): ModelParams => {
    const avgGoals = sumTotal / Math.max(1, nCal);
    const goalsPerElo = sxy / Math.max(1e-9, sxx);
    let bestRho = -0.05, bestLL = -Infinity;
    for (let rho = -0.2; rho <= 0.1 + 1e-9; rho += 0.02) {
      let ll = 0;
      for (const m of calib) {
        const sup = goalsPerElo * m.d;
        const lh = Math.min(8, Math.max(0.15, avgGoals / 2 + sup / 2));
        const la = Math.min(8, Math.max(0.15, avgGoals / 2 - sup / 2));
        ll += Math.log(Math.max(1e-12, poissonPmf(lh, m.h) * poissonPmf(la, m.a) * dcTau(m.h, m.a, lh, la, rho)));
      }
      if (ll > bestLL) { bestLL = ll; bestRho = rho; }
    }
    return { goalsPerElo, avgGoals, rho: Math.round(bestRho * 100) / 100, homeAdvElo: HOME_ADV_ELO };
  };

  const predict = (eh: number, ea: number, homeAdv: number): Probs => {
    const { lambdaHome, lambdaAway } = lambdasFromElo(eh, ea, params, homeAdv);
    const p = outcomeProbs(scoreMatrix(lambdaHome, lambdaAway, params.rho));
    return { h: p.home, d: p.draw, a: p.away };
  };

  for (const r of rows) {
    if (!paramsFixed && r.date >= TRAIN_END) { params = finalizeParams(); paramsFixed = true; }

    const home = canonicalTeam(r.home), away = canonicalTeam(r.away);
    const eh = get(home), ea = get(away);
    const homeAdv = r.neutral ? 0 : HOME_ADV_ELO;
    const o = outcomeOf(r.hs, r.as);

    // accumulate calibration only in the param-fitting window
    if (r.date >= PARAM_FROM && r.date < TRAIN_END) {
      const effDiff = eh - ea + homeAdv;
      sumTotal += r.hs + r.as; nCal++; sxy += effDiff * (r.hs - r.as); sxx += effDiff * effDiff;
      calib.push({ d: effDiff, h: r.hs, a: r.as });
    }

    // score out-of-sample
    if (r.date >= TRAIN_END && r.date < TEST_TO) {
      const p = predict(eh, ea, homeAdv);
      testPreds.push({ p, o });
      baseUniform.push({ p: { h: 1 / 3, d: 1 / 3, a: 1 / 3 }, o });
      baseRate.push({ p: rate, o });
      // elo-only: expected score → W/A split, fixed draw rate
      const e = 1 / (1 + 10 ** ((ea - eh - homeAdv) / 400));
      baseElo.push({ p: { h: e * (1 - rate.d), d: rate.d, a: (1 - e) * (1 - rate.d) }, o });
      reliabilityRaw.push({ p: p.h, y: o === 'h' ? 1 : 0 }, { p: p.d, y: o === 'd' ? 1 : 0 }, { p: p.a, y: o === 'a' ? 1 : 0 });
      if (r.tournament.includes('FIFA World Cup') && !r.tournament.includes('qualification')) {
        wcN++; if (argmax(p) === o) wcCorrect++;
      }
    }

    // update Elo (always)
    const k = importanceWeight(r.tournament);
    const d = eloDelta(r.hs, r.as, eh, ea, k, homeAdv);
    elo.set(home, eh + d); elo.set(away, ea - d);
  }

  // reliability bins (10) + ECE
  const bins = Array.from({ length: 10 }, () => ({ sp: 0, sy: 0, n: 0 }));
  for (const { p, y } of reliabilityRaw) {
    const idx = Math.min(9, Math.floor(p * 10));
    bins[idx]!.sp += p; bins[idx]!.sy += y; bins[idx]!.n++;
  }
  const reliability = bins.map((b) => ({ predicted: b.n ? b.sp / b.n : 0, observed: b.n ? b.sy / b.n : 0, count: b.n }));
  const totalN = reliabilityRaw.length;
  const ece = bins.reduce((s, b) => s + (b.n ? (b.n / totalN) * Math.abs(b.sp / b.n - b.sy / b.n) : 0), 0);

  const out = {
    generatedAt: new Date().toISOString(),
    paramFrom: PARAM_FROM,
    trainEnd: TRAIN_END,
    testTo: TEST_TO,
    tested: testPreds.length,
    params,
    model: scoreSet(testPreds),
    baselines: {
      uniform: scoreSet(baseUniform),
      baseRate: scoreSet(baseRate),
      eloOnly: scoreSet(baseElo),
    },
    reliability,
    ece: Math.round(ece * 1000) / 1000,
    worldCup: { tested: wcN, accuracy: wcN ? wcCorrect / wcN : 0 },
  };

  mkdirSync(OUT_DIR, { recursive: true });
  writeFileSync(join(OUT_DIR, 'backtest.json'), JSON.stringify(out));
  console.log(`backtest → public/data/backtest.json  (${out.tested} test matches, ${TRAIN_END}→${TEST_TO})`);
  console.log(`  model:   brier ${out.model.brier.toFixed(3)}  logloss ${out.model.logloss.toFixed(3)}  rps ${out.model.rps.toFixed(3)}  acc ${(out.model.accuracy * 100).toFixed(1)}%`);
  console.log(`  uniform: brier ${out.baselines.uniform.brier.toFixed(3)}  logloss ${out.baselines.uniform.logloss.toFixed(3)}  rps ${out.baselines.uniform.rps.toFixed(3)}`);
  console.log(`  baseRate:brier ${out.baselines.baseRate.brier.toFixed(3)}  rps ${out.baselines.baseRate.rps.toFixed(3)}  | eloOnly rps ${out.baselines.eloOnly.rps.toFixed(3)}`);
  console.log(`  ECE ${out.ece}  |  WC accuracy ${(out.worldCup.accuracy * 100).toFixed(1)}% (${out.worldCup.tested})`);
}

main();