#!/usr/bin/env python3 """Batched (GPU) backtest for the Stage-3 flash-crash detector. The full cascade scores Stage-3 once per tick (a TCN forward per tick), which is correct in production but makes offline backtests impractically slow. This script precomputes every stride-1 200-tick window score in one batched GPU pass, then replays the (fixed) cascade decision cheaply: alert(t) <= s3(t) >= stage3_threshold AND Stage-5 posterior >= alert_threshold Stage-5 fusion uses the trained Stage-2 score (cheap) + neutral Stage-4 (0.5). Stage-1 is advisory (never a veto). Crash windows use the same 2%-drop-in-60s label rule as the cascade backtest, and alerts are coalesced by --cooldown-ms. Usage: PYTHONPATH=ml python scripts/run_backtest_batched.py \ --data data/parquet/LUNAUSDT_2022-05-10.parquet --model models/stage3_tcn_prod.pt """ import argparse import logging import sys from pathlib import Path import numpy as np import pandas as pd import torch ML_DIR = Path(__file__).resolve().parent.parent / "ml" sys.path.insert(0, str(ML_DIR)) from flash_crash_watchdog.data.historical_loader import df_to_ticks # noqa: E402 from flash_crash_watchdog.data.labels import label_crashes # noqa: E402 from flash_crash_watchdog.features import FEATURE_NAMES, FeatureExtractor # noqa: E402 from flash_crash_watchdog.models.stage3_tcn import normalize_z, TCNConfig, TCNDetector # noqa: E402 torch.serialization.add_safe_globals([TCNConfig]) from flash_crash_watchdog.models.stage2_isolation_forest import Stage2IsolationForest # noqa: E402 from flash_crash_watchdog.models.stage5_bayesian import Stage5Config, Stage5Bayesian # noqa: E402 logger = logging.getLogger(__name__) TCN_FEATURES = FEATURE_NAMES[:17] W = 200 # trained window length NORM_WINDOW = 500 # rolling-z window (matches training / Stage3TCN.feed) CH = 4096 # GPU batch per forward chunk def device() -> str: return "cuda" if torch.cuda.is_available() else "cpu" def stage3_scores(model: TCNDetector, F: np.ndarray, dev: str) -> np.ndarray: """Every stride-1 200-window Stage-3 score (batched GPU).""" if len(F) <= W: return np.zeros(0, dtype=np.float32) norm = normalize_z(F, NORM_WINDOW) n = len(norm) - W + 1 out = np.empty(n, dtype=np.float32) for s in range(0, n, CH): e = min(s + CH, n) win = np.stack([norm[i:i + W] for i in range(s, e)]) x = torch.from_numpy(np.ascontiguousarray(win)).permute(0, 2, 1).float().to(dev) with torch.no_grad(): out[s:e] = model(x)[:, -1].cpu().numpy() return out def score_day(model: TCNDetector, dev: str, parquet: str | Path, max_ticks: int, crash_pct: float): """Feature-extract + score one day. Returns (scores, times, ticks, mids, crash_ranges).""" df = pd.read_parquet(parquet) if max_ticks > 0 and len(df) > max_ticks: df = df.iloc[: max_ticks] ticks = list(df_to_ticks(df, symbol="DAY")) extractor = FeatureExtractor() F = np.zeros((len(ticks), 17), dtype=np.float32) times = np.empty(len(ticks), dtype=np.int64) mids = np.full(len(ticks), np.nan, dtype=np.float64) for i, t in enumerate(ticks): fd = extractor.extract(t) F[i] = [float(fd.get(k, 0.0)) or 0.0 for k in TCN_FEATURES] times[i] = t.book.timestamp_ms mp = t.book.mid_price if mp is not None: mids[i] = mp F = np.nan_to_num(F, nan=0.0, posinf=0.0, neginf=0.0) filled = pd.Series(mids).ffill().bfill().to_numpy() mids = np.where(np.isnan(filled), 1.0, filled) scores = stage3_scores(model, F, dev) # scores[j] => window ending at tick j+W-1 crashes = label_crashes(ticks, drop_threshold_pct=crash_pct) crash_ranges = [(c.start_ts, c.end_ts) for c in crashes] return scores, times, ticks, mids, crash_ranges def evaluate(scores, times, ticks, s2, s5, threshold, crash_ranges, cooldown_ms, mids=None, confirmation_drop_pct=0.0, drop_lookback=0, min_tv_bps=0.0): """Replay + metrics for one threshold. Returns (alerts,tp,fp,fn,prec,rec,f1,med_ttd,rate). When confirmation_drop_pct > 0, an alert fires only if the mid is in a REAL sustained decline: mids[alert_tick] <= mids[alert_tick - drop_lookback] * (1 - confirmation_drop_pct/100). This rejects transient normal-day spikes while keeping crash-day sustained drops. drop_lookback==0 => W//2 (100 ticks). When min_tv_bps > 0, an alert fires only if the window's trailing realized volatility (std/mean of the mid over the 200-tick window, in bps) is above the threshold. Normal-day false positives are on ~0-vol flat windows; crash onsets are on high-vol windows — so this suppresses calm-day chatter while keeping crash detections (and it is observable in real time). """ n_scores = len(scores) if drop_lookback <= 0: drop_lookback = W // 2 alerts = [] last_ts = None for j, tick_idx in enumerate(np.arange(n_scores) + W - 1): s3v = float(scores[j]) if s3v < threshold: continue tick = ticks[tick_idx] s2v, _ = s2.score(tick) alert = s5.aggregate(tick, s2v, s3v, 0.5) if alert is None: continue ts = times[tick_idx] if min_tv_bps > 0: if mids is None or tick_idx < W: continue seg = mids[tick_idx - W + 1:tick_idx + 1] mn = float(seg.mean()) tv = float(seg.std() / mn) * 10000.0 if mn > 0 else 0.0 if tv < min_tv_bps: continue if confirmation_drop_pct > 0: if mids is None or tick_idx < drop_lookback: continue if not (mids[tick_idx] <= mids[tick_idx - drop_lookback] * (1.0 - confirmation_drop_pct / 100.0)): continue if last_ts is not None and ts - last_ts < cooldown_ms: continue last_ts = ts alerts.append(ts) matched = set() ttd = [] tp = 0 for ts in alerts: hit = [i for i, (a, b) in enumerate(crash_ranges) if a <= ts <= b] if hit: tp += 1 matched.add(hit[0]) ttd.append(crash_ranges[hit[0]][1] - ts) fn = len(crash_ranges) - len(matched) fp = len(alerts) - tp prec = tp / (tp + fp) if tp + fp else 0.0 rec = tp / (tp + fn) if tp + fn else 0.0 f1 = 2 * prec * rec / (prec + rec) if prec + rec else 0.0 return (len(alerts), tp, fp, fn, prec, rec, f1, float(np.median(ttd)) if ttd else 0.0, 100.0 * len(alerts) / max(1, len(ticks))) def main() -> int: parser = argparse.ArgumentParser(description="Batched GPU Stage-3 backtest") parser.add_argument("--data", required=True, help="Tick parquet") parser.add_argument("--model", default="models/stage3_tcn_prod.pt") parser.add_argument("--stage2", default="models/stage2_isolation_forest.joblib") parser.add_argument("--threshold", type=float, default=0.5, help="Stage-3 pass gate (operating point: 0.5; or use --thresholds)") parser.add_argument("--thresholds", default=None, help="comma list to sweep, e.g. '0.3,0.4,0.5,0.6,0.7'") parser.add_argument("--stage5-threshold", type=float, default=0.5) parser.add_argument("--cooldown-ms", type=int, default=10000) parser.add_argument("--confirmation-drop-pct", type=float, default=0.0, help="0.0 = OFF. Confirmation gate: an alert counts only when the " "mid is in a REAL sustained decline: " "mids[alert_tick] <= mids[alert_tick - drop_lookback] * (1 - pct/100). " "Rejects transient normal-day spikes, keeps crash-day sustained drops.") parser.add_argument("--drop-lookback", type=int, default=0, help="ticks back to compare mid for the confirmation gate; 0 = W//2 (100). " "Ignored when --confirmation-drop-pct == 0.") parser.add_argument("--min-trailing-vol-bps", type=float, default=0.0, help="Regime gate (bps): only alert if the window's trailing realized " "volatility (mid std/mean over 200 ticks) exceeds this. Suppresses " "calm-day chatter (normal FPs are ~0 bps; crash onsets are 10s of bps). " "0.0 = off.") parser.add_argument("--out-json", default=None, help="Write the threshold-sweep table to a JSON file (for the leaderboard)") parser.add_argument("--crash-pct", type=float, default=2.0) parser.add_argument("--max-ticks", type=int, default=0, help="0 = all") parser.add_argument("--calibrate-on", default=None, help="RSR-11: a SEPARATE parquet day used ONLY to select the " "threshold (max F1). The eval day (--data) is never used " "for threshold selection, which would overfit the metric.") args = parser.parse_args() logging.basicConfig(level=logging.INFO, format="%(asctime)s %(levelname)s %(message)s") dev = device() st = torch.load(args.model, map_location="cpu", weights_only=True) # MLOPS-06 cfg = st["config"] if isinstance(st["config"], TCNConfig) else TCNConfig(**st["config"]) model = TCNDetector(cfg).to(dev) model.load_state_dict(st["model_state"]) model.eval() s2 = Stage2IsolationForest() if Path(args.stage2).exists(): s2.load(args.stage2) s5 = Stage5Bayesian(Stage5Config(alert_threshold=args.stage5_threshold)) scores, times, ticks, mids, crash_ranges = score_day(model, dev, args.data, args.max_ticks, args.crash_pct) logger.info("Found %d ground-truth crash windows on eval day", len(crash_ranges)) thrs = [float(t) for t in args.thresholds.split(",")] if args.thresholds else [args.threshold] # RSR-11: threshold selection must never happen on the evaluation day. if args.calibrate_on: cal = score_day(model, dev, args.calibrate_on, args.max_ticks, args.crash_pct) best_f1, best_thr = -1.0, thrs[0] for thr in thrs: _an, tp, fp, fn, _prec, rec, f1, _med, _rate = evaluate( cal[0], cal[1], cal[2], s2, s5, thr, cal[4], args.cooldown_ms, mids=cal[3], confirmation_drop_pct=args.confirmation_drop_pct, drop_lookback=args.drop_lookback, min_tv_bps=args.min_trailing_vol_bps) if f1 > best_f1: best_f1, best_thr = f1, thr logger.info("RSR-11: selected threshold %.3f on calibration day %s (F1=%.3f); " "eval day %s is untouched by selection", best_thr, Path(args.calibrate_on).name, best_f1, Path(args.data).name) thrs = [best_thr] elif len(thrs) > 1: logger.warning( "RSR-11: sweeping %d thresholds ON the evaluation day overfits the metric — " "use --calibrate-on to select the threshold on a separate day", len(thrs)) print("day=%s ticks=%d crashes=%d gate_drop_pct=%.2f drop_lookback=%d min_tv_bps=%.1f" % (Path(args.data).name, len(ticks), len(crash_ranges), args.confirmation_drop_pct, (args.drop_lookback or W // 2), args.min_trailing_vol_bps)) rows = [] print("threshold,alerts,tp,fp,fn,precision,recall,f1,alert_rate_pct,median_ttd_ms") for thr in thrs: alerts_n, tp, fp, fn, prec, rec, f1, med, rate = evaluate( scores, times, ticks, s2, s5, thr, crash_ranges, args.cooldown_ms, mids=mids, confirmation_drop_pct=args.confirmation_drop_pct, drop_lookback=args.drop_lookback, min_tv_bps=args.min_trailing_vol_bps) print("%.3f,%d,%d,%d,%d,%.3f,%.3f,%.3f,%.3f,%.0f" % (thr, alerts_n, tp, fp, fn, prec, rec, f1, rate, med)) rows.append({"threshold": thr, "alerts": alerts_n, "tp": tp, "fp": fp, "fn": fn, "precision": prec, "recall": rec, "f1": f1, "alert_rate_pct": rate, "median_ttd_ms": med}) if args.out_json: import json as _json payload = { "day": Path(args.data).name, "ticks": len(ticks), "crashes": len(crash_ranges), "gate_bps": args.min_trailing_vol_bps, "cooldown_ms": args.cooldown_ms, "model": Path(args.model).name, "thresholds": rows, } Path(args.out_json).write_text(_json.dumps(payload, indent=2)) logger.info("wrote leaderboard table -> %s", args.out_json) return 0 if __name__ == "__main__": raise SystemExit(main())