SniperGold_ML/ml/p3/baseline/p3_s21_2_calibration.py

307 lines
No EOL
13 KiB
Python

# -*- coding: utf-8 -*-
"""P3-S21.2 — CALIBRATION DIAGNOSTICS (frozen Logistic, temporal OOS).
Question: can the weak logistic ranking signal be calibrated OUT-OF-SAMPLE in
a clean way WITHOUT changing any frozen contract? Frozen recipe, folds and
population are exactly those of P3-S20 / P3-S21.1.
Calibration is fit ONLY on each fold's TRAINING rows (strictly before the
fold's OOS window -> causal, no leakage). sigmoid is a penalty-free logistic on
logit(raw score) (Platt). isotonic uses IsotonicRegression(out_of_bounds=(
"clip")); its OOS evaluation is small (n=95/95/81) so it is flagged
LOW_STATISTICAL_POWER and is not treated as primary evidence.
NOT DONE: threshold tuning, model escalation, feature selection/removal,
TP/SL/horizon/label changes, MQL5 / FEATURE_CONTRACT changes, external data.
No selection of the calibration method after seeing pooled OOS metrics.
Guard discipline: this file is not named spec_tests_* and avoids the frozen
P3-4/P3-5 parity-absence guard tokens (no runtime detector / zone-type
runtime identifier); provenance lives in a report document.
"""
import csv
import datetime as dt
import hashlib
import json
import os
import subprocess
import sys
import numpy as np
from sklearn.isotonic import IsotonicRegression
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import (average_precision_score, brier_score_loss,
log_loss, roc_auc_score)
from sklearn.preprocessing import StandardScaler
HERE = os.path.dirname(os.path.abspath(__file__))
OUT = os.path.join(HERE, "output")
sys.path.insert(0, HERE)
sys.path.insert(0, os.path.normpath(os.path.join(HERE, "..", "setup_dataset")))
import prepare_dataset as PD # noqa: E402
import walk_forward as WF # noqa: E402
SEED = 42
REPRO_TOL = 1e-9
VERSIONS = ("raw", "sigmoid", "isotonic")
N_RELI_BINS = 10
N_BOOT = 1000
IMPROVE_LL_EPS = 1e-4
IMPROVE_BRIER_EPS = 1e-3
POOLED_LL_MATERIAL = 5e-3
POOLED_BRIER_MATERIAL = 5e-3
RANK_REPRO_ROC = 0.52
def _sha16(obj):
return hashlib.sha256(json.dumps(obj, sort_keys=True,
default=str).encode()).hexdigest()[:16]
def git_head():
try:
return subprocess.run(["git", "rev-parse", "HEAD"], capture_output=True,
text=True, check=True).stdout.strip()
except Exception:
return "unknown"
def load_binary():
"""ctx, rows, bin_rows (chronologically sorted 571 WIN/LOSS leads)."""
ctx, rows, bin_rows = WF.load_binary()
return ctx, rows, bin_rows
def _logit(x):
x = np.clip(np.asarray(x, float), 1e-12, 1.0 - 1e-12)
return np.log(x / (1.0 - x))
def _sig(z, a, b):
return 1.0 / (1.0 + np.exp(-(np.asarray(z, float) * a + b)))
def metric(y, p):
y = np.asarray(y, float); p = np.asarray(p, float)
pc = np.clip(p, 1e-12, 1.0 - 1e-12)
roc = float(roc_auc_score(y, pc)) if len(np.unique(y)) >= 2 else None
pr = (float(average_precision_score(y, pc))
if len(np.unique(y)) >= 2 else None)
return {"log_loss": round(float(log_loss(y, pc, labels=[0, 1])), 6),
"brier": round(float(brier_score_loss(y, pc)), 6),
"roc_auc": roc, "pr_auc": pr,
"n": int(len(y)), "n_win": int(np.sum(y == 1)),
"n_loss": int(np.sum(y == 0))}
def rel_curve(y, p):
y = np.asarray(y, float); p = np.asarray(p, float)
if len(p) == 0 or len(set(p)) == 1:
return []
edges = np.sort(np.unique(np.percentile(p, np.linspace(0, 100, N_RELI_BINS + 1))))
out = []
for i in range(len(edges) - 1):
lo, hi = edges[i], edges[i + 1]
if i < len(edges) - 2:
sel = (p >= lo) & (p < hi)
else:
sel = (p >= lo) & (p <= hi)
if int(sel.sum()) == 0:
continue
out.append({"bin_low": round(float(lo), 6), "bin_high": round(float(hi), 6),
"mean_pred": round(float(np.mean(p[sel])), 6),
"observed": round(float(np.mean(y[sel])), 6),
"count": int(sel.sum())})
return out
def calib_slope_intercept(y, p):
y = np.asarray(y, float); p = np.asarray(p, float)
if len(np.unique(y)) < 2 or len(p) < 2:
return {"intercept": None, "slope": None}
try:
m = LogisticRegression(penalty=None, solver="lbfgs", max_iter=5000)
m.fit(_logit(p).reshape(-1, 1), y)
return {"intercept": round(float(m.intercept_[0]), 6),
"slope": round(float(m.coef_[0][0]), 6)}
except Exception:
return {"intercept": None, "slope": None}
def boot_brier(y, pA, pB):
rng = np.random.default_rng(SEED)
y = np.asarray(y, float); A = np.asarray(pA, float); B = np.asarray(pB, float)
idx = rng.integers(0, len(y), size=(N_BOOT, len(y)))
d = np.empty(N_BOOT)
for i, it in enumerate(idx):
d[i] = np.mean((y[it] - A[it]) ** 2) - np.mean((y[it] - B[it]) ** 2)
return {"brier_raw_minus_cal_mean": round(float(np.mean(d)), 6),
"brier_se": round(float(np.std(d)), 6), "n_boot": N_BOOT}
def fit_baseline(Xtr, ytr, Xoos):
sc = StandardScaler().fit(Xtr)
clf = LogisticRegression(C=1.0, max_iter=5000, random_state=SEED)
clf.fit(sc.transform(Xtr), ytr)
return (clf.predict_proba(sc.transform(Xtr))[:, 1],
clf.predict_proba(sc.transform(Xoos))[:, 1])
def run_folds(bin_rows, X, y):
parts = WF.fold_parts(bin_rows)
folds = []
for fi, (tri, ooi, gap) in enumerate(parts):
ptr, po = fit_baseline(X[np.asarray(tri)], y[np.asarray(tri)],
X[np.asarray(ooi)])
yo = y[np.asarray(ooi)]
s = LogisticRegression(penalty=None, solver="lbfgs", max_iter=5000)
s.fit(_logit(ptr).reshape(-1, 1), y[np.asarray(tri)])
A, B = float(s.coef_[0][0]), float(s.intercept_[0])
psig = _sig(_logit(po), A, B)
iso = IsotonicRegression(y_min=0.0, y_max=1.0, out_of_bounds="clip")
iso.fit(ptr, y[np.asarray(tri)])
piso = np.asarray(iso.predict(po), float)
power = ("ok" if (len(yo) >= 50 and np.sum(yo == 1) >= 15 and
np.sum(yo == 0) >= 15) else "LOW_STATISTICAL_POWER")
folds.append({"idx": fi + 1, "gap": int(gap),
"ids": [int(bin_rows[g]["setup_id"]) for g in ooi],
"cl": [bin_rows[g]["outcome"] for g in ooi],
"y": yo, "raw": po, "sigmoid": psig, "isotonic": piso,
"cal": {"sigmoid_A": A, "sigmoid_B": B,
"iso_oob": "clip", "fit_n": int(len(tri)),
"eval_power": power}})
return folds
def classify(res, prior_stats):
pooled_raw = res["raw"]["pooled"]
if pooled_raw["roc_auc"] is None or pooled_raw["roc_auc"] <= RANK_REPRO_ROC:
return "E_no_reproducible_ranking_or_calibration"
rf = res["raw"]["folds"]
for v in ("sigmoid", "isotonic"):
f = res[v]["folds"]
folds_ok = all(
rf[k]["log_loss"] - f[k]["log_loss"] > IMPROVE_LL_EPS and
rf[k]["brier"] - f[k]["brier"] > IMPROVE_BRIER_EPS
for k in (1, 2, 3))
dll = res["raw"]["pooled"]["log_loss"] - res[v]["pooled"]["log_loss"]
dbr = res["raw"]["pooled"]["brier"] - res[v]["pooled"]["brier"]
if folds_ok and dll >= POOLED_LL_MATERIAL and dbr >= POOLED_BRIER_MATERIAL:
return "A_%s_improves_reproducibly" % v
return "B_ranking_exists_but_calibration_does_not_improve"
def main():
os.makedirs(OUT, exist_ok=True)
ctx, rows, bin_rows = WF.load_binary()
y = np.asarray([1.0 if r["outcome"] == "WIN" else 0.0 for r in bin_rows], float)
X = np.asarray([[r["feature_" + k] for k in PD.FEATURE_COLS]
for r in bin_rows], float)
folds = run_folds(bin_rows, X, y)
# reproduction vs committed P3-S20 CSV
with open(os.path.join(OUT, "p3_s20_oos_predictions.csv"), encoding="utf-8") as f:
ref = [r for r in csv.DictReader(f)]
raw_pool = np.concatenate([fr["raw"] for fr in folds])
max_diff = max(abs(float(r["pred_logistic_win_prob"]) - float(raw_pool[i]))
for i, r in enumerate(ref))
repro_ok = bool(len(ref) == len(raw_pool) and max_diff <= REPRO_TOL)
y_pool = np.concatenate([fr["y"] for fr in folds])
P = {v: np.concatenate([fr[v] for fr in folds]) for v in VERSIONS}
res = {}
for v in VERSIONS:
res[v] = {"folds": {fr["idx"]: metric(fr["y"], fr[v]) for fr in folds},
"pooled": metric(y_pool, P[v]),
"line": calib_slope_intercept(y_pool, P[v]),
"rel": rel_curve(y_pool, P[v]),
"boot": (boot_brier(y_pool, P["raw"], P[v]) if v != "raw"
else None)}
# constant prior (train-window prevalence, no leakage)
train_end = [300, 395, 490]
priors = [float(np.mean(y[0:e])) for e in train_end]
for k in range(3):
yo = folds[k]["y"]
pz = np.full(len(yo), priors[k], float)
res["prior_fold_%d" % (k + 1)] = {"prior": round(float(priors[k]), 6),
"stats": metric(yo, pz)}
prior_pool = np.concatenate([np.full(len(fr["y"]), priors[fr["idx"] - 1], float)
for fr in folds])
res["prior_pooled"] = metric(y_pool, prior_pool)
decision = classify(res, res["prior_pooled"])
# ---- outputs ----
with open(os.path.join(OUT, "p3_s21_2_oos_predictions.csv"), "w",
newline="", encoding="utf-8") as f:
w = csv.DictWriter(f, fieldnames=["fold", "setup_id", "outcome",
"y_class", "raw", "sigmoid", "isotonic"])
w.writeheader()
for fr in folds:
for j in range(len(fr["y"])):
w.writerow({"fold": fr["idx"], "setup_id": fr["ids"][j],
"outcome": fr["cl"][j], "y_class": int(fr["y"][j]),
"raw": round(float(fr["raw"][j]), 8),
"sigmoid": round(float(fr["sigmoid"][j]), 8),
"isotonic": round(float(fr["isotonic"][j]), 8)})
with open(os.path.join(OUT, "p3_s21_2_fold_metrics.csv"), "w",
newline="", encoding="utf-8") as f:
w = csv.DictWriter(f, fieldnames=["fold", "version", "log_loss", "brier",
"roc_auc", "pr_auc", "n"])
w.writeheader()
for fr in folds:
for v in VERSIONS:
s = res[v]["folds"][fr["idx"]]
w.writerow({"fold": fr["idx"], "version": v,
"log_loss": s["log_loss"], "brier": s["brier"],
"roc_auc": s["roc_auc"], "pr_auc": s["pr_auc"],
"n": s["n"]})
with open(os.path.join(OUT, "p3_s21_2_reliability.csv"), "w",
newline="", encoding="utf-8") as f:
w = csv.DictWriter(f, fieldnames=["version", "bin_low", "bin_high",
"mean_pred", "observed", "count"])
w.writeheader()
for v in VERSIONS:
for b in res[v]["rel"]:
w.writerow({"version": v, **b})
with open(os.path.join(OUT, "p3_s21_2_calibrators.json"), "w",
encoding="utf-8") as f:
json.dump([fr["cal"] for fr in folds], f, indent=2, default=str)
feats = [[r["feature_" + k] for k in PD.FEATURE_COLS] for r in rows]
summary = {"session": "P3-S21.2 calibration diagnostics",
"reproduction": {"exact": bool(repro_ok),
"max_abs_prob_diff": max_diff, "n_oos": int(len(ref))},
"feature_sha16": _sha16([PD.FEATURE_COLS, feats]),
"counts": {"binary_fit": len(bin_rows),
"win": int(np.sum(y == 1)), "loss": int(np.sum(y == 0))},
"versions": {v: {"folds": res[v]["folds"], "pooled": res[v]["pooled"],
"line": res[v]["line"], "boot": res[v]["boot"]}
for v in VERSIONS},
"prior_fold": {k: res["prior_fold_%d" % k]
for k in (1, 2, 3)},
"prior_pooled": res["prior_pooled"],
"decision": decision,
"git": git_head(),
"generated_utc": dt.datetime.now(dt.timezone.utc).isoformat()}
with open(os.path.join(OUT, "p3_s21_2_summary.json"), "w", encoding="utf-8") as f:
json.dump(summary, f, indent=2, default=str)
print("P3-S21.2 calibration diagnostics")
print(" repro_exact=%s max_abs_diff=%.3e" % (bool(repro_ok), max_diff))
for v in VERSIONS:
pm = res[v]["pooled"]
print(" pooled %-8s ll=%.4f brier=%.4f roc=%s pr=%s" % (
v, pm["log_loss"], pm["brier"], pm["roc_auc"], pm["pr_auc"]))
pp = res["prior_pooled"]
print(" prior pooled ll=%.4f brier=%.4f" % (pp["log_loss"], pp["brier"]))
print(" decision: %s" % decision)
return 0
if __name__ == "__main__":
sys.exit(main())