Buckets:
| """cwm3/eval/road_dynamics_eval.py — F5c four-arm road-dynamics evaluation -> results/tables/road_dynamics.md (+ .json; rows appended to results/paired_compare.csv). | |
| WP4 batch 1, pre-registration docs/decisions.md 2026-09-16T06:46:35Z and the road-side row 2026-09-16T14:10:00Z. CPU, numpy only. | |
| Question: can road DEMAND be predicted from the road's own dynamics rather than only recognised from the current signal (kill_report_final K7)? | |
| Targets = results/targets_road_dyn/road_dyn_{regime}_test.npz (future lead_dRel / lead_vRel / lead_vLead / lane_center_offset at t+0.2k, k = 1..25, validity). | |
| Arms (per-anchor forecasts mu[N,25,4]): | |
| persistence y_k = y_0 (train_road_dynamics.persistence_np) | |
| cv constant velocity / linear lane offset (train_road_dynamics.cv_forecast_np) | |
| dyn_s<k> GRU residual net on the 8 s CAN road history + 40 CAN anchor features (results/runs/RD_dyn_s<k>/pred_rd_{regime}_test.npz) | |
| dynvis_s<k> dyn + the frozen per-anchor vf256 (results/runs/RD_dynvis_s<k>/...) | |
| Metrics per regime (test split): | |
| (a) MAE per channel at 0.2 / 1 / 2 / 3 / 5 s on rows valid at the anchor and at k (learned arms: mean +- sd over seeds); paired device-cluster bootstrap CI of | |
| MAE(arm) - MAE(cv) at 3 s (shared resamples, n_boot 300). | |
| (b) Gaussian NLL per channel at 1 / 3 / 5 s: learned arms from (mu, sigma); rule arms with a sigma per (k, channel) fitted on the VAL residuals ('+ val-fitted sigma'). | |
| (c) Derived demand-onset PR-AUC@3 s per family (ttc, lane_excursion, lead_braking) on EXACTLY the rows of eval/pressure_eval.populations (inactive & valid at the | |
| anchor, valid over the first 15 steps, moving) — the population of the FINAL model's P chains. PRIMARY score of a forecast = max over k <= 15 of the | |
| threshold-normalised margin ratio of the forecast MEAN (ttc: 4 / ttc_k with ttc_k = dRel_k / max(-vRel_k, 0.1); lateral: 0.9 / (w_0/2 - |off_k|) with the current | |
| lane width; lead braking: aLead_k / -2 with aLead_k = the 1 s backward difference of the forecast vLead, history-filled for k < 5); persistence reproduces | |
| the current-margin rule S0 for ttc / lateral (its lead-braking score IS the current lead deceleration margin). SECONDARY (learned arms) = Monte-Carlo onset | |
| probability from 32 independent Gaussian draws per (k, channel). Comparators: the margin rule S0 (pressure pack; = current_signal_S0 for the road | |
| families) and the FINAL model's P chains `I1_all_minus_road2drv_s<k>_final` (p_onset_p at 3 s; seed-matched for the learned arms). Paired device-cluster | |
| bootstrap (shared resamples, n_boot 500, seed 0); verdict over the 3 seeds = kill_checks.beats (report_joint.verdict) where 3 paired rows exist. | |
| Reading (fixed in the pre-registration row): 'road demand is PREDICTABLE beyond recognition' iff a learned arm BEATS S0 on >= 1 family in >= 1 regime; otherwise | |
| K7 stands as recognition-only. Everything is reported regardless of direction. | |
| python eval/road_dynamics_eval.py [--regimes manual,engaged] [--seeds 0,1,2] [--n-boot 500] [--no-append] | |
| No openpilot DMS field is read (inputs: canspec CAN road channels of the road_dyn pack, the RD predictions, the FINAL layouts' P block). | |
| """ | |
| from __future__ import annotations | |
| import argparse | |
| import csv | |
| import json | |
| import sys | |
| import time | |
| from pathlib import Path | |
| import numpy as np | |
| import pandas as pd | |
| CWM3 = Path(__file__).resolve().parents[1] | |
| EXPLORE = CWM3.parent | |
| for p in (EXPLORE / "cwm2" / "eval", EXPLORE / "s12_wm", CWM3 / "common", CWM3 / "eval", CWM3 / "train"): | |
| sys.path.insert(0, str(p)) | |
| import importlib.util as _ilu # noqa: E402 | |
| _spec = _ilu.spec_from_file_location("cwm3_runmeta", CWM3 / "common" / "runmeta.py") | |
| _rm = _ilu.module_from_spec(_spec); _spec.loader.exec_module(_rm) | |
| write_run_json = _rm.write_run_json | |
| from pred_layout import load_pred, col # noqa: E402 | |
| from paired_compare import ap_weighted, COLS, OUT as PC_OUT # noqa: E402 | |
| import pressure_eval as PE # noqa: E402 | |
| import report_joint as RJ # noqa: E402 | |
| import train_road_dynamics as TRD # noqa: E402 | |
| RD = CWM3 / "results" / "targets_road_dyn" | |
| RUNS = CWM3 / "results" / "runs" | |
| TAB = CWM3 / "results" / "tables" | |
| FINAL = "I1_all_minus_road2drv_s{k}_final" | |
| CH = ["lead_dRel", "lead_vRel", "lead_vLead", "lane_center_offset"] | |
| CH_SHORT = ["dRel", "vRel", "vLead", "lane_off"] | |
| FAMS = ["ttc", "lane_excursion", "lead_braking"] | |
| K_MAE = [1, 5, 10, 15, 25] | |
| K_NLL = [5, 15, 25] | |
| TAU3 = 15 | |
| THR = dict(ttc=4.0, veh_half_w=0.9, lead_brake_a=-2.0, ttc_none=30.0, lane_w_min=2.5, lane_w_max=5.0) | |
| N_MC = 32 | |
| # ----------------------------------------------------------------------------- forecasts | |
| def load_rd(regime, split="test"): | |
| with np.load(RD / f"road_dyn_{regime}_{split}.npz", allow_pickle=True) as z: | |
| d = {k: z[k] for k in ("H", "Hvalid", "Yf", "Yvalid", "M0", "M0valid", "uid", "t", "dongle")} | |
| d["H32"] = d["H"].astype(np.float32) | |
| d["av"] = TRD.anchor_valid_np(d["H32"]) | |
| return d | |
| def rule_arms(d): | |
| return {"persistence": TRD.persistence_np(d["H32"]), "cv": TRD.cv_forecast_np(d["H32"], d["Hvalid"])} | |
| def learned_arm(arm, seed, regime, split="test"): | |
| p = RUNS / f"RD_{arm}_s{seed}" / f"pred_rd_{regime}_{split}.npz" | |
| if not p.exists() or not (RUNS / f"RD_{arm}_s{seed}" / "PRED_DONE").exists(): | |
| return None | |
| with np.load(p, allow_pickle=True) as z: | |
| return {"mu": z["mu"].astype(np.float32), "sigma": z["sigma"].astype(np.float32), "row_ok": z["row_ok"].astype(bool), "uid": z["uid"], "t": z["t"]} | |
| def event_scores(mu, d, persistence=False): | |
| """PRIMARY per-family onset score of a deterministic forecast mu[N,25,4] (max over k <= 15 of the threshold-normalised margin ratio); 0 where the family's | |
| signal is absent at the anchor. persistence=True: the lead-braking score is the current 1 s lead deceleration margin (M0 aLead / -2 = the S0 rule).""" | |
| N = len(mu) | |
| H32, av, M0 = d["H32"], d["av"], d["M0"] | |
| m = mu[:, :TAU3] | |
| # ttc | |
| approaching = m[:, :, 1] < 0 | |
| ttc = np.where(approaching, m[:, :, 0] / np.maximum(-m[:, :, 1], 0.1), THR["ttc_none"]) | |
| s_ttc = np.where(av[:, 0], (THR["ttc"] / np.maximum(ttc, 0.1)).max(1), 0.0) | |
| # lane excursion (current lane width, lane usable at the anchor) | |
| w0 = M0[:, 6] | |
| lane_ok = av[:, 3] & (w0 >= THR["lane_w_min"]) & (w0 <= THR["lane_w_max"]) | |
| lat = w0[:, None] / 2.0 - np.abs(m[:, :, 3]) | |
| s_lane = np.where(lane_ok, (THR["veh_half_w"] / np.maximum(lat, 0.05)).max(1), 0.0) | |
| # lead braking: 1 s backward difference of the forecast vLead at each k <= 15 (history-filled for k < 5) | |
| if persistence: | |
| s_brk = np.where(av[:, 0] & d["M0valid"][:, 4], M0[:, 4] / THR["lead_brake_a"], 0.0) | |
| else: | |
| vl = m[:, :, 2] | |
| ratios = np.zeros((N, TAU3), np.float32) | |
| for k in range(1, TAU3 + 1): # k = 1..15 (index k-1); vLead at t + 0.2 (k-5) | |
| j = k - 5 | |
| if j >= 1: | |
| prev, pv = vl[:, j - 1], np.ones(N, bool) | |
| elif j == 0: | |
| prev, pv = H32[:, 39, 2], d["Hvalid"][:, 39] & (H32[:, 39, 6] > 0.5) | |
| else: | |
| prev, pv = H32[:, 39 + j, 2], d["Hvalid"][:, 39 + j] & (H32[:, 39 + j, 6] > 0.5) | |
| ratios[:, k - 1] = np.where(pv, (vl[:, k - 1] - prev) / 1.0 / THR["lead_brake_a"], 0.0) | |
| s_brk = np.where(av[:, 0], ratios.max(1), 0.0) | |
| return {"ttc": s_ttc.astype(np.float64), "lane_excursion": s_lane.astype(np.float64), "lead_braking": s_brk.astype(np.float64)} | |
| def event_probs(mu, sigma, d, n_mc=N_MC, seed=0, chunk=40_000): | |
| """SECONDARY score: Monte-Carlo probability that the derived event is active at some k <= 15 (independent Gaussian draws per (k, channel)).""" | |
| rng = np.random.default_rng(seed) | |
| N = len(mu) | |
| out = {f: np.zeros(N, np.float64) for f in FAMS} | |
| av, M0, H32 = d["av"], d["M0"], d["H32"] | |
| w0 = M0[:, 6] | |
| lane_ok = av[:, 3] & (w0 >= THR["lane_w_min"]) & (w0 <= THR["lane_w_max"]) | |
| for i in range(0, N, chunk): | |
| sl = slice(i, min(i + chunk, N)) | |
| m, s = mu[sl, :TAU3], sigma[sl, :TAU3] | |
| n = m.shape[0] | |
| hit = {f: np.zeros(n, np.float64) for f in FAMS} | |
| for _ in range(n_mc): | |
| x = m + s * rng.standard_normal(m.shape).astype(np.float32) | |
| ttc = np.where(x[:, :, 1] < 0, x[:, :, 0] / np.maximum(-x[:, :, 1], 0.1), THR["ttc_none"]) | |
| hit["ttc"] += (ttc < THR["ttc"]).any(1) | |
| lat = w0[sl, None] / 2.0 - np.abs(x[:, :, 3]) | |
| hit["lane_excursion"] += (lat < THR["veh_half_w"]).any(1) | |
| vl = x[:, :, 2] | |
| a = np.zeros((n, TAU3), np.float32) | |
| for k in range(1, TAU3 + 1): | |
| j = k - 5 | |
| if j >= 1: | |
| a[:, k - 1] = vl[:, k - 1] - vl[:, j - 1] | |
| else: | |
| hv = d["Hvalid"][sl, 39 + j] & (H32[sl, 39 + j, 6] > 0.5) | |
| a[:, k - 1] = np.where(hv, vl[:, k - 1] - H32[sl, 39 + j, 2], 0.0) | |
| hit["lead_braking"] += (a < THR["lead_brake_a"]).any(1) | |
| out["ttc"][sl] = np.where(av[sl, 0], hit["ttc"] / n_mc, 0.0) | |
| out["lane_excursion"][sl] = np.where(lane_ok[sl], hit["lane_excursion"] / n_mc, 0.0) | |
| out["lead_braking"][sl] = np.where(av[sl, 0], hit["lead_braking"] / n_mc, 0.0) | |
| return out | |
| # ----------------------------------------------------------------------------- bootstrap | |
| def ap_sorted(y_sorted, w_sorted): | |
| """weighted average precision with rows pre-sorted by decreasing score (ties in that order) — O(N) per weight vector.""" | |
| tp = np.cumsum(w_sorted * y_sorted) | |
| fp = np.cumsum(w_sorted * (~y_sorted)) | |
| pos = tp[-1] | |
| if pos <= 0: | |
| return np.nan | |
| prec = tp / np.maximum(tp + fp, 1e-12) | |
| return float((prec * w_sorted * y_sorted).sum() / pos) | |
| def boot_ap_pair(y, sa, sb, dongle, n_boot, seed=0): | |
| """paired device-cluster bootstrap of PR-AUC(A) - PR-AUC(B) with shared resamples (paired_compare.boot_pair 'ap' semantics, pre-sorted for speed).""" | |
| y = np.asarray(y, bool) | |
| oa, ob = np.argsort(-sa, kind="mergesort"), np.argsort(-sb, kind="mergesort") | |
| ya, yb = y[oa], y[ob] | |
| uniq, dev = np.unique(dongle, return_inverse=True) | |
| ones = np.ones(len(y)) | |
| pa, pb = ap_sorted(ya, ones), ap_sorted(yb, ones) | |
| rng = np.random.default_rng(seed) | |
| da, db = np.full(n_boot, np.nan), np.full(n_boot, np.nan) | |
| for i in range(n_boot): | |
| cnt = np.bincount(rng.integers(0, len(uniq), len(uniq)), minlength=len(uniq)).astype(np.float64) | |
| w = cnt[dev] | |
| if (w * y).sum() <= 0 or (w * y).sum() >= w.sum(): | |
| continue | |
| da[i], db[i] = ap_sorted(ya, w[oa]), ap_sorted(yb, w[ob]) | |
| dd = da - db | |
| ok = np.isfinite(dd) | |
| pct = lambda v: (float(np.nanpercentile(v, 2.5)), float(np.nanpercentile(v, 97.5))) # noqa: E731 | |
| return {"value_a": pa, "value_b": pb, "ci_a": pct(da), "ci_b": pct(db), "delta": pa - pb, "ci_d": pct(dd), "frac_gt0": float((dd[ok] > 0).mean()) if ok.any() else np.nan, | |
| "n_dev": int(len(uniq)), "n": int(len(y)), "n_pos": int(y.sum())} | |
| def boot_mean_pair(va, vb, dongle, n_boot, seed=0): | |
| uniq, dev = np.unique(dongle, return_inverse=True) | |
| rng = np.random.default_rng(seed) | |
| d = np.full(n_boot, np.nan) | |
| for i in range(n_boot): | |
| cnt = np.bincount(rng.integers(0, len(uniq), len(uniq)), minlength=len(uniq)).astype(np.float64) | |
| w = cnt[dev] | |
| d[i] = ((w * va).sum() - (w * vb).sum()) / w.sum() | |
| return {"delta": float(va.mean() - vb.mean()), "ci_d": (float(np.nanpercentile(d, 2.5)), float(np.nanpercentile(d, 97.5)))} | |
| def final_scores(regime, seed, tg_uid, tg_t): | |
| rd = RUNS / FINAL.format(k=seed) | |
| p = rd / f"pred_layout_{regime}_test.npz" | |
| if not p.exists(): | |
| return None | |
| pa = load_pred(p, {"uid": tg_uid, "t": tg_t}, keys=("p_onset_p", "p_names", "onset_taus")) | |
| names = [str(x) for x in pa["p_names"]] | |
| ti = list(map(int, pa["onset_taus"])).index(TAU3) | |
| return {f: pa["p_onset_p"][:, ti, col(names, f)].astype(np.float64) for f in FAMS} | |
| def verdict_rows(rows): | |
| """rows: list of dicts with seed / delta / delta_ci_lo / delta_ci_hi / value_a / value_b -> kill_checks.beats verdict (3 rows) or a CI-only flag.""" | |
| if len(rows) >= 3: | |
| g = pd.DataFrame(rows) | |
| return RJ.verdict(g)["verdict"] | |
| if len(rows) == 1: | |
| lo, hi = rows[0]["delta_ci_lo"], rows[0]["delta_ci_hi"] | |
| return "CI>0" if lo > 0 else ("CI<0" if hi < 0 else "CI incl. 0") | |
| return "INSUFFICIENT" | |
| def fmt_ms(v): | |
| v = [x for x in v if np.isfinite(x)] | |
| return "—" if not v else (f"{np.mean(v):.3f} ± {np.std(v, ddof=1):.3f}" if len(v) > 1 else f"{v[0]:.3f}") | |
| # ----------------------------------------------------------------------------- main | |
| def main(): | |
| ap = argparse.ArgumentParser() | |
| ap.add_argument("--regimes", default="manual,engaged") | |
| ap.add_argument("--seeds", default="0,1,2") | |
| ap.add_argument("--n-boot", type=int, default=500) | |
| ap.add_argument("--n-boot-mae", type=int, default=300) | |
| ap.add_argument("--no-append", action="store_true") | |
| ap.add_argument("--out", default=str(TAB / "road_dynamics.md")) | |
| a = ap.parse_args() | |
| t0 = time.time() | |
| stamp = time.strftime("%Y%m%dT%H%M%SZ", time.gmtime()) | |
| seeds = [int(s) for s in a.seeds.split(",")] | |
| J = {"stamp": stamp, "regimes": {}, "arms": {}, "n_boot": a.n_boot} | |
| L = ["# Road dynamics — four arms (F5c): can road demand be PREDICTED rather than only recognised (K7)?", "", | |
| f"Pre-registration docs/decisions.md 2026-09-16T06:46:35Z (WP4 batch 1) and 2026-09-16T14:10:00Z (road-side definitions and reading). Test split, both regimes; " | |
| f"generated {stamp} by eval/road_dynamics_eval.py. Targets `results/targets_road_dyn` (future lead dRel / vRel / vLead and lane-centre offset at k = 1..25, canspec " | |
| f"whitelist CAN through grid10 columns ≥ 12 only). Arms: persistence → constant velocity (CV) → `RD_dyn` (GRU-128 on the 8 s CAN road history + 40 CAN anchor " | |
| f"features, Gaussian residual over CV, 149 k params) → `RD_dynvis` (+ frozen vf256, 183 k params); 3 seeds each, VAL-NLL selection. Comparators for the derived " | |
| f"demand onsets: the current-margin rule S0 (recognition) and the FINAL model's P chains `I1_all_minus_road2drv_s<k>_final` (seed-matched). Paired device-cluster " | |
| f"bootstrap with shared resamples (PR-AUC n_boot {a.n_boot}; MAE n_boot {a.n_boot_mae}); verdicts over 3 seeds = kill_checks.beats (≥ 2/3 CIs excl. 0, same sign, " | |
| f"|mean Δ| > seed sd).", "", | |
| "**Reading rule (pre-registered):** 'road demand is predictable beyond recognition' iff a learned arm BEATS S0 on ≥ 1 family in ≥ 1 regime (primary score); otherwise K7 stands as recognition-only.", ""] | |
| csv_rows = [] | |
| headline = [] | |
| for regime in a.regimes.split(","): | |
| d = load_rd(regime) | |
| dval = load_rd(regime, "val") | |
| N = len(d["t"]) | |
| arms = rule_arms(d) | |
| arms_val = rule_arms(dval) | |
| sig_rule = {} | |
| for name in ("persistence", "cv"): # val-fitted sigma per (k, channel) for the rule arms' NLL | |
| r = dval["Yf"] - arms_val[name] | |
| m = dval["Yvalid"] & dval["av"][:, None, :] | |
| sig_rule[name] = np.array([[r[:, k, c][m[:, k, c]].std() + 1e-3 for c in range(4)] for k in range(25)], np.float32) | |
| learned = {} | |
| for arm in ("dyn", "dynvis", "dynfull"): # dynfull = POST-HOC full-pass sensitivity arm (labelled; not part of the pre-registered reading) | |
| for s in seeds: | |
| p = learned_arm(arm, s, regime) | |
| if p is not None: | |
| assert np.array_equal(p["uid"], d["uid"]) and np.allclose(p["t"], d["t"]), (arm, s, "row alignment") | |
| learned[f"{arm}_s{s}"] = p | |
| J["regimes"][regime] = {"N": int(N), "learned_arms": sorted(learned)} | |
| L += [f"## {regime} (N = {N:,} test anchors; learned arms available: {', '.join(sorted(learned)) or 'none yet'})", ""] | |
| # ---------------- (a) MAE | |
| mask = d["Yvalid"] & d["av"][:, None, :] | |
| def mae_table(mu): | |
| return np.array([[np.abs(mu[:, k - 1, c] - d["Yf"][:, k - 1, c])[mask[:, k - 1, c]].mean() if mask[:, k - 1, c].any() else np.nan for c in range(4)] for k in K_MAE]) | |
| mae = {n: mae_table(m) for n, m in arms.items()} | |
| for n, p in learned.items(): | |
| mae[n] = mae_table(p["mu"]) | |
| L += ["### (a) MAE of the forecast (rows valid at the anchor and at k; learned arms mean ± sd over seeds)", "", | |
| "| arm | channel | " + " | ".join(f"{k * 0.2:.1f} s" for k in K_MAE) + " | n rows @3 s |", "|---|---|" + "---|" * len(K_MAE) + "---|"] | |
| for arm in ("persistence", "cv", "dyn", "dynvis", "dynfull"): | |
| for c in range(4): | |
| if arm in mae: | |
| vals = [f"{mae[arm][i, c]:.3f}" for i in range(len(K_MAE))] | |
| else: | |
| ks = [k for k in mae if k.startswith(arm + "_s")] | |
| if not ks: | |
| continue | |
| vals = [fmt_ms([mae[k][i, c] for k in ks]) for i in range(len(K_MAE))] | |
| L.append(f"| {arm} | {CH_SHORT[c]} | " + " | ".join(vals) + f" | {int(mask[:, 14, c].sum()):,} |") | |
| L.append("") | |
| # paired MAE delta vs cv at 3 s | |
| L += ["#### ΔMAE at 3 s vs constant velocity (paired device-cluster bootstrap; negative = better than CV)", "", "| arm | channel | Δ per seed [CI] |", "|---|---|---|"] | |
| mae_json = {} | |
| for arm in ("persistence", "dyn", "dynvis", "dynfull"): | |
| ks = [arm] if arm in arms else [k for k in learned if k.startswith(arm + "_s")] | |
| for c in range(4): | |
| cells = [] | |
| mk = mask[:, 14, c] | |
| eb = np.abs(arms["cv"][:, 14, c] - d["Yf"][:, 14, c])[mk] | |
| for k in ks: | |
| mu = arms[k] if k in arms else learned[k]["mu"] | |
| ea = np.abs(mu[:, 14, c] - d["Yf"][:, 14, c])[mk] | |
| r = boot_mean_pair(ea, eb, d["dongle"][mk], a.n_boot_mae) | |
| cells.append(f"{r['delta']:+.3f} [{r['ci_d'][0]:+.3f},{r['ci_d'][1]:+.3f}]") | |
| mae_json[f"{k}|{CH_SHORT[c]}"] = r | |
| if cells: | |
| L.append(f"| {arm} | {CH_SHORT[c]} | " + " / ".join(cells) + " |") | |
| L.append("") | |
| # ---------------- (b) NLL | |
| def nll_table(mu, sigma): | |
| out = np.full((len(K_NLL), 4), np.nan) | |
| for i, k in enumerate(K_NLL): | |
| for c in range(4): | |
| m = mask[:, k - 1, c] | |
| if m.any(): | |
| s = sigma[:, k - 1, c] if sigma.ndim == 3 else np.full(N, sigma[k - 1, c]) | |
| r = (d["Yf"][:, k - 1, c] - mu[:, k - 1, c])[m] | |
| out[i, c] = np.mean(0.5 * (r / s[m]) ** 2 + np.log(s[m]) + 0.5 * np.log(2 * np.pi)) | |
| return out | |
| nll = {n: nll_table(arms[n], sig_rule[n]) for n in arms} | |
| for n, p in learned.items(): | |
| nll[n] = nll_table(p["mu"], p["sigma"]) | |
| L += ["### (b) Gaussian NLL per channel (nats; rule arms with a per-(k, channel) sigma fitted on VAL; learned arms from their predicted sigma)", "", | |
| "| arm | channel | " + " | ".join(f"{k * 0.2:.0f} s" for k in K_NLL) + " |", "|---|---|" + "---|" * len(K_NLL)] | |
| for arm in ("persistence", "cv", "dyn", "dynvis", "dynfull"): | |
| for c in range(4): | |
| if arm in nll: | |
| vals = [f"{nll[arm][i, c]:.3f}" for i in range(len(K_NLL))] | |
| else: | |
| ks = [k for k in nll if k.startswith(arm + "_s")] | |
| if not ks: | |
| continue | |
| vals = [fmt_ms([nll[k][i, c] for k in ks]) for i in range(len(K_NLL))] | |
| L.append(f"| {arm}{' + val-fitted σ' if arm in arms else ''} | {CH_SHORT[c]} | " + " | ".join(vals) + " |") | |
| L.append("") | |
| # ---------------- (c) derived onset PR-AUC@3 s | |
| tg, pops = PE.populations(regime, "test") | |
| assert np.array_equal(tg["uid"], d["uid"]) and np.allclose(tg["t"], d["t"]) | |
| fin = {s: final_scores(regime, s, d["uid"], d["t"]) for s in seeds} | |
| scores = {"persistence": event_scores(arms["persistence"], d, persistence=True), "cv": event_scores(arms["cv"], d)} | |
| for n, p in learned.items(): | |
| scores[n] = event_scores(p["mu"], d) | |
| scores[n + "_prob"] = event_probs(p["mu"], p["sigma"], d) | |
| L += ["### (c) Derived demand-onset PR-AUC@3 s per family (population = pressure_eval rows; A = arm score, B = comparator; Δ = A − B)", "", | |
| "Primary score = max over k ≤ 15 of the forecast mean's margin ratio; `_prob` = Monte-Carlo onset probability (secondary). Persistence ≡ the margin rule S0 (Δ vs S0 = 0 by construction for ttc / lateral / lead braking).", "", | |
| "| family | arm | PR-AUC A (per seed) | vs S0: Δ per seed [CI] | verdict vs S0 | vs FINAL P chain: Δ per seed [CI] | verdict vs FINAL | n / pos |", "|---|---|---|---|---|---|---|---|"] | |
| cj = {} | |
| for fam in FAMS: | |
| pop = pops[fam] | |
| sel, y, dg, s0 = pop["sel"], pop["y"], pop["dongle"], pop["margin"].astype(np.float64) | |
| base = float(y.mean()) | |
| L.append(f"| {fam} | *base rate* | {base:.4f} | | | | | {len(sel):,} / {int(y.sum()):,} |") | |
| arm_groups = [("persistence", ["persistence"]), ("cv", ["cv"])] | |
| for arm in ("dyn", "dynvis", "dynfull"): | |
| ks = [k for k in learned if k.startswith(arm + "_s")] | |
| if ks: | |
| arm_groups += [(arm, ks), (arm + "_prob", [k + "_prob" for k in ks])] | |
| arm_groups.append(("FINAL P chain", [f"final_s{s}" for s in seeds if fin[s] is not None])) | |
| for gname, members in arm_groups: | |
| rows_s0, rows_fin, vals = [], [], [] | |
| for mem in members: | |
| if mem.startswith("final_s"): | |
| s = int(mem[7:]) | |
| sa = fin[s][fam][sel] | |
| else: | |
| sa = scores[mem][fam][sel] | |
| keep = np.isfinite(sa) | |
| if keep.sum() < 1000 or y[keep].sum() < 20: | |
| continue | |
| seed = int(mem.split("_s")[1].split("_")[0]) if "_s" in mem else 0 | |
| r = boot_ap_pair(y[keep], sa[keep], s0[keep], dg[keep], a.n_boot) | |
| vals.append(r["value_a"]) | |
| rows_s0.append({"seed": seed, "delta": r["delta"], "delta_ci_lo": r["ci_d"][0], "delta_ci_hi": r["ci_d"][1], "value_a": r["value_a"], "value_b": r["value_b"], "status": "secondary"}) | |
| csv_rows.append(dict(run_stamp=stamp, tag=f"roaddyn_s{seed}_{mem.replace('final_s' + str(seed), 'final')}_vs_S0", regime=regime, split="test", cell=f"rd_onset|{fam}|3s", metric="pr_auc", | |
| a=mem, b="margin_rule_S0", value_a=r["value_a"], ci_a_lo=r["ci_a"][0], ci_a_hi=r["ci_a"][1], value_b=r["value_b"], ci_b_lo=r["ci_b"][0], ci_b_hi=r["ci_b"][1], | |
| delta=r["delta"], delta_ci_lo=r["ci_d"][0], delta_ci_hi=r["ci_d"][1], frac_boot_d_gt0=r["frac_gt0"], n=r["n"], n_pos=r["n_pos"], n_dev=r["n_dev"], n_boot=a.n_boot)) | |
| if not mem.startswith("final_s"): | |
| fseeds = [seed] if (mem.startswith("dyn") and fin.get(seed) is not None) else [s for s in seeds if fin[s] is not None] | |
| for fs in fseeds: | |
| sb = fin[fs][fam][sel] | |
| kf = keep & np.isfinite(sb) | |
| r2 = boot_ap_pair(y[kf], sa[kf], sb[kf], dg[kf], a.n_boot) | |
| rows_fin.append({"seed": fs, "delta": r2["delta"], "delta_ci_lo": r2["ci_d"][0], "delta_ci_hi": r2["ci_d"][1], "value_a": r2["value_a"], "value_b": r2["value_b"], "status": "secondary"}) | |
| csv_rows.append(dict(run_stamp=stamp, tag=f"roaddyn_s{fs}_{mem}_vs_final", regime=regime, split="test", cell=f"rd_onset|{fam}|3s", metric="pr_auc", | |
| a=mem, b=FINAL.format(k=fs), value_a=r2["value_a"], ci_a_lo=r2["ci_a"][0], ci_a_hi=r2["ci_a"][1], value_b=r2["value_b"], ci_b_lo=r2["ci_b"][0], ci_b_hi=r2["ci_b"][1], | |
| delta=r2["delta"], delta_ci_lo=r2["ci_d"][0], delta_ci_hi=r2["ci_d"][1], frac_boot_d_gt0=r2["frac_gt0"], n=r2["n"], n_pos=r2["n_pos"], n_dev=r2["n_dev"], n_boot=a.n_boot)) | |
| if not vals: | |
| continue | |
| v_s0, v_fin = verdict_rows(rows_s0), (verdict_rows(rows_fin) if rows_fin else "—") | |
| cj[f"{fam}|{gname}"] = {"pr_auc": vals, "vs_S0": rows_s0, "verdict_vs_S0": v_s0, "vs_final": rows_fin, "verdict_vs_final": v_fin, "n": int(len(sel)), "base_rate": base} | |
| f_s0 = " / ".join(f"{r['delta']:+.3f} [{r['delta_ci_lo']:+.3f},{r['delta_ci_hi']:+.3f}]" for r in rows_s0) | |
| f_fin = " / ".join(f"{r['delta']:+.3f} [{r['delta_ci_lo']:+.3f},{r['delta_ci_hi']:+.3f}]" for r in rows_fin) or "—" | |
| L.append(f"| {fam} | {gname} | {' / '.join(f'{v:.3f}' for v in vals)} | {f_s0} | {v_s0} | {f_fin} | {v_fin} | {len(sel):,} / {int(y.sum()):,} |") | |
| if gname in ("dyn", "dynvis") and v_s0 == "BEATS": # the pre-registered reading uses the two pre-registered arms only (dynfull is descriptive) | |
| headline.append(f"{regime}|{fam}|{gname}") | |
| L.append("") | |
| J["regimes"][regime].update({"mae": {k: v.tolist() for k, v in mae.items()}, "mae_delta_vs_cv_3s": mae_json, "nll": {k: v.tolist() for k, v in nll.items()}, "onset": cj, | |
| "k_mae": K_MAE, "k_nll": K_NLL, "channels": CH}) | |
| del d, dval, arms, learned, scores | |
| # ---------------- verdict | |
| n_learned = sum(len(v["learned_arms"]) for v in J["regimes"].values()) | |
| if n_learned == 0: | |
| verdict, detail = "NOT DECIDABLE", "no learned arm finished yet (rule arms only)" | |
| elif headline: | |
| verdict, detail = "PREDICTABLE BEYOND RECOGNITION (pre-registered criterion met)", "learned arm BEATS the margin rule S0 on: " + ", ".join(headline) | |
| else: | |
| verdict, detail = "RECOGNITION-ONLY (K7 stands)", "no learned arm BEATS the margin rule S0 on any family in either regime (primary score)" | |
| J.update(verdict=verdict, detail=detail, headline_cells=headline) | |
| L.insert(6, f"## VERDICT: **{verdict}** — {detail}") | |
| L.insert(7, "") | |
| out = Path(a.out) | |
| out.parent.mkdir(parents=True, exist_ok=True) | |
| out.write_text("\n".join(L) + "\n") | |
| (out.parent / (out.stem + ".json")).write_text(json.dumps(J, indent=1, default=float)) # json name follows --out (post-hoc re-runs never overwrite the pre-registered files) | |
| if not a.no_append and csv_rows: | |
| new = not PC_OUT.exists() | |
| with open(PC_OUT, "a", newline="") as f: | |
| w = csv.DictWriter(f, fieldnames=COLS) | |
| if new: | |
| w.writeheader() | |
| for r in csv_rows: | |
| w.writerow(r) | |
| print(f"appended {len(csv_rows)} rows to {PC_OUT}") | |
| rd_eval = RUNS / "RD_eval" | |
| write_run_json(rd_eval, f"road_dynamics_eval_{stamp}", [Path(__file__), Path(TRD.__file__)], [RD / f"road_dyn_{r}_test.npz" for r in a.regimes.split(",")], 0, t0, | |
| extra={"verdict": verdict, "detail": detail, "table": str(out), "n_boot": a.n_boot, "dms_fields_used": []}, input_listing=False) | |
| print(f"VERDICT {verdict}: {detail}\nwritten {out} ({time.time() - t0:.0f}s)") | |
| if __name__ == "__main__": | |
| main() | |
Xet Storage Details
- Size:
- 28.2 kB
- Xet hash:
- 7da23a06251fa6a17ec56bf4e660a0c357ec838509f482291de2e841eca38c85
·
Xet efficiently stores files, intelligently splitting them into unique chunks and accelerating uploads and downloads. More info.