Buckets:

AsianPlayer
/

TriDrive

Files
xet
AsianPlayer/TriDrive / code /aide /eval /report_hm_subsets.py
AsianPlayer
16 days ago
download
raw
52.9 kB
#!/usr/bin/env python
"""cwm3/aide/eval/report_hm_subsets.py — rule-based alternatives to Driver-WM's official High-Motion (HM) subset.
Pre-registered in docs/decisions.md (row of 2026-09-15, UTC 2026-09-15T03:01:43Z) BEFORE any metric was computed.  Everything here is
DESCRIPTIVE: every subset is reported, none is used for selection or adoption; the official top-10 % HM (n = 60, threshold 72.05 px)
remains the only subset on which Driver-WM's published HM h=5 = 155.82 applies.
Reads ONLY saved test predictions (`aide/results/runs/<run>/pred_test.npz`, read-only), the ground truth (`load_split("test")`, future
block = protocol frames 5..9, all 136 joints, no confidence mask) and the one-label-per-clip AIDE test labels.  No training, no GPU.
Subsets (all rules on the ground truth, all deterministic):
  (A) motion-quantile sweep — the official score (Driver-WM Eq. 16 = `metrics.motion_score`: mean GT displacement over the 4 future
      transitions and all joints), top q of the 609 test clips for q = 5 / 10 / 20 / 30 / 50 / 100 % (n = floor(q N) = 30 / 60 / 121 / 182 / 304 / 609);
  (B) part-specific top 10 % (n = 60): the same score restricted to the 42 hand joints (94–135), the 68 face joints (26–93), the 26 body joints (0–25);
  (C) label-defined subsets: every class of every AIDE task (DBR 7 / DER 5 / TCR 3 / VCR 5) and the pre-specified union "reactive manoeuvre"
      = VCR in {Turning, Backward Moving, Changing Lane} OR DBR in {Looking Around, Body Movement}.
Metrics per subset: MPJPE per horizon and horizon-averaged, d-nMPJPE, PCK@0.05 (96 px); 10k clip-bootstrap CI (seed 0); PAIRED per-clip Δ
(10k, seed 0; the `metrics.paired_diff` construction, verified bit-identical) of every model run vs zero-velocity, the two val-fitted
ZV-shrink rules (α = 0.44; α_h = 0.76·h/5), and — seed-matched — the Driver-WM architecture on CockpitWM tokens and our Driver-WM
reproduction.  Verdict per (model, subset, comparator) = `report_arch_attribution.verdict` (BEATS / WORSE = CI excl. 0 on >= ceil(2n/3) seeds,
same sign on all; NOT_WORSE = all CIs cover 0 and abs(Δ) < 0.25 px; else INCONCLUSIVE).
Integrity: the q = 10 % subset must reproduce every run's HM MPJPE / HM d-n / PCK@0.05 HM / HM h1..h5 of its harness record
`metrics_test.json` and of `master_aide.json`, and the q = 100 % subset the All columns (tolerance --tol, default 1e-6); the script aborts otherwise.
Outputs (new files only): aide/results/tables/hm_subsets.md + .json; aide/results/figs/hm_quantile_curves.{pdf,png}, hm_part_subsets.png.
    python aide/eval/report_hm_subsets.py [--n-boot 10000] [--tol 1e-6] [--out aide/results/tables] [--figs aide/results/figs]
"""
from __future__ import annotations

import os

for _v in ("OMP_NUM_THREADS", "MKL_NUM_THREADS", "OPENBLAS_NUM_THREADS", "NUMEXPR_NUM_THREADS"):
    os.environ.setdefault(_v, "4")

import argparse  # noqa: E402
import json  # noqa: E402
import math  # noqa: E402
import sys  # noqa: E402
import time  # noqa: E402
from pathlib import Path  # noqa: E402

import numpy as np  # noqa: E402

ROOT = Path(__file__).resolve().parents[2]
sys.path.insert(0, str(ROOT))
from aide.data import protocol as P  # noqa: E402
from aide.eval import metrics as M  # noqa: E402
from aide.eval.baselines_analytic import zv_shrink, zv_shrink_h  # noqa: E402
from aide.eval.harness import align, load_split, zero_velocity  # noqa: E402

RUNS = P.RESULTS / "runs"
TABLES = P.RESULTS / "tables"
FIGS = P.RESULTS / "figs"
PREREG_UTC = "2026-09-15T03:01:43Z"
SEEDS3 = [42, 123, 456]
SEEDS5 = [42, 123, 456, 789, 1011]
TOL_NOT_WORSE = 0.25
QUANTILES = [0.05, 0.10, 0.20, 0.30, 0.50, 1.00]
PART_JOINTS = {"hands": np.arange(94, 136), "face": np.arange(26, 94), "body": np.arange(0, 26)}
PART_LABEL = {"hands": "hand motion (42 hand joints, both hands)", "face": "head/face motion (68 face joints)", "body": "body motion (26 body joints)"}
REACTIVE_VCR = ("Turning", "Backward Moving", "Changing Lane")
REACTIVE_DBR = ("Looking Around", "Body Movement")
TASK_NAME = {"dbr": "DBR", "der": "DER", "tcr": "TCR", "vcr": "VCR"}

# key, display label (manuscript names), run template, seeds, master_aide.json group key
MODELS = [
    dict(key="K", label="CockpitWM", tmpl="cwm3d_kp_s{s}_v2", seeds=SEEDS5, master="K5"),
    dict(key="Kpre", label="CockpitWM + BATON pretraining", tmpl="cwm3d_kp_batonpre_trunkp123_s{s}_v2", seeds=SEEDS5, master="K5_pre"),
    dict(key="KS", label="CockpitWM (+ semantic heads)", tmpl="cwm3d_vis_sem_s{s}_v2", seeds=SEEDS3, master="KS3"),
    dict(key="DWMkp", label="Driver-WM architecture on CockpitWM tokens", tmpl="driverwm_kp_s{s}_v2", seeds=SEEDS3, master="DWMkp3"),
    dict(key="DWMrepro", label="Driver-WM reproduction (ours)", tmpl="driverwm_repro_s{s}_v2", seeds=SEEDS3, master="DWMrepro3"),
    dict(key="MB", label="MotionBERT architecture under CockpitWM formulation (supplement-only control)", tmpl="motionbert_kp_s{s}_v2", seeds=SEEDS5, master="MB5"),
]
# key, label, saved run (for the exact-reproduction check), rule
ANALYTIC = [
    dict(key="zv", label="zero-velocity", run="zero_velocity"),
    dict(key="zvs", label="ZV-shrink α=0.44 (val-fitted)", run="zv_shrink_val"),
    dict(key="zvsh", label="ZV-shrink α_h=0.76·h/5 (val-fitted)", run="zv_shrink_h_val"),
]
MASTER_KEY_ANALYTIC = {"zv": "zv", "zvs": "zv_shrink", "zvsh": "zv_shrink_h"}
COMPARATORS = ["zv", "zvs", "zvsh", "DWMkp", "DWMrepro"]          # sections A and B
COMPARATORS_C = ["zv", "DWMkp", "DWMrepro"]                        # section C (pre-registered)
COMP_LABEL = {"zv": "zero-velocity", "zvs": "ZV-shrink α=0.44", "zvsh": "ZV-shrink α_h=0.76·h/5",
              "DWMkp": "Driver-WM architecture on CockpitWM tokens (seed-matched)", "DWMrepro": "Driver-WM reproduction (ours) (seed-matched)"}
COMP_SHORT = {"zv": "ZV", "zvs": "ZV-shrink α", "zvsh": "ZV-shrink α_h", "DWMkp": "Driver-WM arch. on CockpitWM tokens", "DWMrepro": "Driver-WM reproduction (ours)"}


# ------------------------------------------------------------------------------------------------------------ helpers
def verdict(rows, higher_better=False, n_req=None):
    """Verbatim copy of aide/eval/report_arch_attribution.verdict (baton_pretrain.md rule generalised to n seeds)."""
    n = len(rows)
    if n == 0:
        return "INSUFFICIENT"
    n_req = math.ceil(2 * n / 3) if n_req is None else n_req
    d = np.array([r["mean"] for r in rows]); lo = np.array([r["ci_lo"] for r in rows]); hi = np.array([r["ci_hi"] for r in rows])
    improving = (d > 0) if higher_better else (d < 0)
    excl = (lo > 0) | (hi < 0)
    same = np.all(np.sign(d) == np.sign(d.mean())) and d.mean() != 0
    if excl.sum() >= n_req and same and improving.all():
        return "BEATS"
    if excl.sum() >= n_req and same and (~improving).all():
        return "WORSE"
    if np.all((lo <= 0) & (hi >= 0)) and np.all(np.abs(d) < TOL_NOT_WORSE):
        return "NOT_WORSE"
    return "INCONCLUSIVE"


class Boot:
    """Clip-bootstrap indices shared by every statistic of one subset; identical construction to metrics.bootstrap_ci / paired_diff
    (default_rng(seed).integers(0, n, (n_boot, n)) for the CI; default_rng(seed + 1) for p_lt0)."""

    def __init__(self, n: int, n_boot: int, seed: int = 0):
        self.n, self.n_boot, self.seed = n, n_boot, seed
        self.idx = np.random.default_rng(seed).integers(0, n, size=(n_boot, n))
        self.idx1 = np.random.default_rng(seed + 1).integers(0, n, size=(n_boot, n))

    def ci(self, v):
        v = np.asarray(v, dtype=np.float64)
        b = v[self.idx].mean(1)
        return float(np.quantile(b, 0.025)), float(np.quantile(b, 0.975))

    def paired(self, a, b):
        d = np.asarray(a, np.float64) - np.asarray(b, np.float64)
        lo, hi = self.ci(d)
        return {"mean": float(d.mean()), "ci_lo": lo, "ci_hi": hi, "p_lt0": float((d[self.idx1].mean(1) < 0).mean())}


def fmt(x, nd=2, sign=False):
    return (f"{x:+.{nd}f}" if sign else f"{x:.{nd}f}")


def ms(vals, nd=2, sign=False):
    v = np.asarray(vals, dtype=np.float64)
    if len(v) == 1:
        return fmt(v[0], nd, sign)
    return f"{fmt(v.mean(), nd, sign)} ± {v.std(ddof=1):.{nd}f}"


def ci_str(lo, hi, nd=2, sign=False):
    return f"[{fmt(lo, nd, sign)}, {fmt(hi, nd, sign)}]"


def agg(vals):
    v = np.asarray(vals, dtype=np.float64)
    return {"mean": float(v.mean()), "sd": (float(v.std(ddof=1)) if len(v) > 1 else None), "n": int(len(v)), "per_seed": [float(x) for x in v]}


def n_excl(rows):
    return int(sum(1 for r in rows if r["ci_lo"] > 0 or r["ci_hi"] < 0))


def top_mask(score, n_top):
    order = np.argsort(-score, kind="stable")
    m = np.zeros(len(score), bool)
    m[order[:n_top]] = True
    thr = float(score[order[n_top - 1]])
    return m, thr


# ------------------------------------------------------------------------------------------------------------ io
def load_pred(name: str, clip_ids: np.ndarray) -> np.ndarray:
    z = np.load(RUNS / name / "pred_test.npz", allow_pickle=False)
    pred = np.asarray(z["pred"], dtype=np.float64)[align(z, clip_ids)]
    assert pred.shape == (len(clip_ids), P.T_PRED, P.J, 2), (name, pred.shape)
    assert np.isfinite(pred).all(), f"{name}: non-finite predictions"
    return pred


def harness_record(name: str) -> dict:
    return json.loads((RUNS / name / "metrics_test.json").read_text())


# ------------------------------------------------------------------------------------------------------------ main
def main(argv=None) -> int:
    ap = argparse.ArgumentParser()
    ap.add_argument("--n-boot", type=int, default=10000)
    ap.add_argument("--tol", type=float, default=1e-6)
    ap.add_argument("--out", type=Path, default=TABLES)
    ap.add_argument("--figs", type=Path, default=FIGS)
    a = ap.parse_args(argv)
    t0 = time.time()
    stamp = time.strftime("%Y%m%dT%H%M%SZ", time.gmtime())

    # ---- ground truth, labels, analytic rules recomputed from the observed frames
    d = load_split("test")
    N = len(d["clip_id"])
    assert N == P.EXPECTED_SPLIT_SIZES["test"] == 609, N
    S10 = d["S10"]
    gt = S10[:, P.T_OBS:]                                        # [N, 5, 136, 3]
    obs = S10[:, :P.T_OBS]                                       # frames 0..4 only (the rules never see the future)
    labels = {t: np.asarray(d[f"y_{t}"]) for t in P.TASKS}
    fits = {k: json.loads((RUNS / r / "alpha_fit_val.json").read_text())["chosen_alpha"] for k, r in (("zvs", "zv_shrink_val"), ("zvsh", "zv_shrink_h_val"))}
    assert fits == {"zvs": 0.44, "zvsh": 0.76}, fits
    rules = {"zv": zero_velocity(S10), "zvs": zv_shrink(obs, fits["zvs"]), "zvsh": zv_shrink_h(obs, fits["zvsh"])}
    preds: dict[str, np.ndarray] = {}
    rule_check = {}
    for an in ANALYTIC:
        p = np.asarray(rules[an["key"]], np.float64).astype(np.float32).astype(np.float64)    # the saved rule predictions are float32
        saved = load_pred(an["run"], d["clip_id"])
        rule_check[an["key"]] = float(np.abs(p - saved).max())
        assert rule_check[an["key"]] == 0.0, (an["key"], rule_check[an["key"]])
        preds[an["key"]] = p
    print(f"[hm] analytic rules recomputed from GT; max abs diff vs saved predictions: {rule_check}")

    # ---- model runs
    runs = {}                     # run name -> dict(model key, seed)
    groups = {}                   # model key -> list of present runs
    missing = {}
    for mm in MODELS:
        present = []
        for s in mm["seeds"]:
            r = mm["tmpl"].format(s=s)
            if (RUNS / r / "pred_test.npz").exists() and (RUNS / r / "metrics_test.json").exists():
                preds[r] = load_pred(r, d["clip_id"])
                runs[r] = {"model": mm["key"], "seed": s}
                present.append(r)
            else:
                missing.setdefault(mm["key"], []).append(r)
        groups[mm["key"]] = present
        print(f"[hm] {mm['label']}: {len(present)} run(s)" + (f", missing {missing[mm['key']]}" if mm["key"] in missing else ""))
    all_names = [an["key"] for an in ANALYTIC] + [r for mm in MODELS for r in groups[mm["key"]]]
    label_of = {an["key"]: an["label"] for an in ANALYTIC}
    label_of.update({r: next(mm["label"] for mm in MODELS if mm["key"] == runs[r]["model"]) for r in runs})

    # ---- per-clip error tables (computed once)
    thr_pck = P.pck_threshold_px(0.05)
    err = {nme: M.joint_errors(preds[nme], gt) for nme in all_names}                 # [N, 5, 136]
    ce_all = {nme: e.mean(axis=(1, 2)) for nme, e in err.items()}                   # [N]
    ce_h = {nme: e.mean(axis=2) for nme, e in err.items()}                          # [N, 5]
    pck_clip = {nme: 100.0 * (e <= thr_pck).mean(axis=(1, 2)) for nme, e in err.items()}
    ce_part = {p: {nme: e[:, :, idx].mean(axis=(1, 2)) for nme, e in err.items()} for p, idx in PART_JOINTS.items()}
    ce_part_h = {p: {nme: e[:, :, idx].mean(axis=2) for nme, e in err.items()} for p, idx in PART_JOINTS.items()}

    # ---- motion scores and subsets
    score_all = M.motion_score(gt)
    score_part = {p: M.motion_score(gt[:, :, idx]) for p, idx in PART_JOINTS.items()}
    hm_official, thr_official, s_chk = M.hm_mask(gt)
    assert np.array_equal(s_chk, score_all) and int(hm_official.sum()) == 60
    subsets = []                 # dict(section, key, label, mask, n, extra)
    for q in QUANTILES:
        n_top = int(q * N)
        m, thr = top_mask(score_all, n_top)
        if q == 0.10:
            assert np.array_equal(m, hm_official) and abs(thr - thr_official) < 1e-9 and n_top == P.hm_top(N) == 60
        subsets.append(dict(section="A", key=f"q{int(round(q * 100)):03d}", q=q, label=f"top {int(round(q * 100))} %", mask=m, n=n_top, threshold=thr))
    for p in PART_JOINTS:
        m, thr = top_mask(score_part[p], P.hm_top(N))
        subsets.append(dict(section="B", key=f"part_{p}", part=p, label=PART_LABEL[p], mask=m, n=int(m.sum()), threshold=thr,
                            overlap_official=int((m & hm_official).sum())))
    for t in P.TASKS:
        for ci, cname in enumerate(P.CLASSES[t]):
            m = labels[t] == ci
            subsets.append(dict(section="C", key=f"{t}_{ci}", task=t, cls=cname, label=f"{TASK_NAME[t]}: {cname}", mask=m, n=int(m.sum())))
    vcr_idx = [P.CLASSES["vcr"].index(c) for c in REACTIVE_VCR]
    dbr_idx = [P.CLASSES["dbr"].index(c) for c in REACTIVE_DBR]
    m_union = np.isin(labels["vcr"], vcr_idx) | np.isin(labels["dbr"], dbr_idx)
    subsets.append(dict(section="C", key="reactive_union", label="reactive manoeuvre (union)", mask=m_union, n=int(m_union.sum()),
                        definition=f"VCR in {list(REACTIVE_VCR)} OR DBR in {list(REACTIVE_DBR)}",
                        n_vcr_part=int(np.isin(labels['vcr'], vcr_idx).sum()), n_dbr_part=int(np.isin(labels['dbr'], dbr_idx).sum()),
                        n_both=int((np.isin(labels['vcr'], vcr_idx) & np.isin(labels['dbr'], dbr_idx)).sum())))
    for sub in subsets:
        sub["mean_score"] = float(score_all[sub["mask"]].mean())
        sub["n_in_official_hm"] = int((sub["mask"] & hm_official).sum())
    print(f"[hm] {len(subsets)} subsets; official HM threshold {thr_official:.4f} px (paper 72.05)")

    # ---- comparator resolution: analytic -> the rule; model comparator -> the same-seed run of that group (None if absent / self)
    def comp_run(run: str, comp: str):
        if comp in rules:
            return comp if run != comp else None
        if run not in runs or runs[run]["model"] == comp:
            return None
        cand = next((mm["tmpl"] for mm in MODELS if mm["key"] == comp)).format(s=runs[run]["seed"])
        return cand if cand in preds else None

    # ---- evaluation per subset
    def level(mask, boot, nme, part=None):
        c = (ce_all if part is None else ce_part[part])[nme][mask]
        ch = (ce_h if part is None else ce_part_h[part])[nme][mask]
        out = {"mpjpe": float(c.mean()), "ci": boot.ci(c), "dn_mpjpe": 100.0 * float(c.mean()) / P.DIAG,
               "h": [float(ch[:, h].mean()) for h in range(P.T_PRED)], "h_ci": [boot.ci(ch[:, h]) for h in range(P.T_PRED)]}
        if part is None:
            pc = pck_clip[nme][mask]
            out["pck05"] = float(pc.mean()); out["pck05_ci"] = boot.ci(pc)
        return out

    def delta(mask, boot, nme, cr, part=None):
        c, cc = (ce_all if part is None else ce_part[part]), (ce_h if part is None else ce_part_h[part])
        out = boot.paired(c[nme][mask], c[cr][mask])
        out["h"] = [boot.paired(cc[nme][mask][:, h], cc[cr][mask][:, h]) for h in range(P.T_PRED)]
        return out

    results = {}
    boot_check = None
    for sub in subsets:
        ts = time.time()
        mask, n = sub["mask"], sub["n"]
        boot = Boot(n, a.n_boot, 0)
        comps = COMPARATORS if sub["section"] in ("A", "B") else COMPARATORS_C
        parts = [sub["part"]] if sub["section"] == "B" else []
        res = {}
        for nme in all_names:
            r = {"level": level(mask, boot, nme), "delta": {}}
            for p in parts:
                r[f"level_{p}"] = level(mask, boot, nme, p)
            for comp in comps:
                cr = comp_run(nme, comp)
                if cr is None:
                    continue
                r["delta"][comp] = delta(mask, boot, nme, cr)
                r["delta"][comp]["vs_run"] = cr
                for p in parts:
                    r["delta"][comp][f"part_{p}"] = delta(mask, boot, nme, cr, p)
            res[nme] = r
        results[sub["key"]] = res
        if boot_check is None and sub["key"] == "q010":      # bit-identity of the shared-index bootstrap with metrics.paired_diff
            r0 = groups["K"][0]
            ref = M.paired_diff(ce_all[r0][mask], ce_all["zv"][mask], a.n_boot, 0)
            mine = res[r0]["delta"]["zv"]
            boot_check = {k: abs(ref[k] - mine[k]) for k in ("mean", "ci_lo", "ci_hi", "p_lt0")}
            assert max(boot_check.values()) == 0.0, boot_check
        print(f"[hm] subset {sub['key']:18s} n={n:3d} done in {time.time() - ts:5.1f} s")

    # ---- integrity: q = 10 % must reproduce the harness / master HM cells, q = 100 % the All cells
    master = json.loads((TABLES / "master_aide.json").read_text())
    checks, mism = [], []

    def chk(run, key, mine, theirs, src):
        dv = abs(mine - theirs)
        checks.append(dv)
        if dv > a.tol:
            mism.append(f"{run}: {key} recomputed {mine!r} vs {src} {theirs!r} (abs diff {dv:.3e})")

    for nme in all_names:
        run_name = nme if nme in runs else next(an["run"] for an in ANALYTIC if an["key"] == nme)
        rec = harness_record(run_name)["metrics"]
        L10, L100 = results["q010"][nme]["level"], results["q100"][nme]["level"]
        chk(run_name, "mpjpe_hm", L10["mpjpe"], rec["mpjpe_hm"]["value"], "metrics_test.json")
        chk(run_name, "dn_mpjpe_hm", L10["dn_mpjpe"], rec["dn_mpjpe_hm"]["value"], "metrics_test.json")
        chk(run_name, "pck05_hm", L10["pck05"], rec["pck05_hm"]["value"], "metrics_test.json")
        for h in range(P.T_PRED):
            chk(run_name, f"hm_h{h + 1}", L10["h"][h], rec[f"hm_h{h + 1}"]["value"], "metrics_test.json")
            chk(run_name, f"mpjpe_h{h + 1}", L100["h"][h], rec[f"mpjpe_h{h + 1}"]["value"], "metrics_test.json")
        chk(run_name, "mpjpe_all", L100["mpjpe"], rec["mpjpe_all"]["value"], "metrics_test.json")
        chk(run_name, "dn_mpjpe", L100["dn_mpjpe"], rec["dn_mpjpe"]["value"], "metrics_test.json")
        chk(run_name, "pck05", L100["pck05"], rec["pck05"]["value"], "metrics_test.json")
        chk(run_name, "diff_vs_ref_hm", results["q010"][nme]["delta"]["zv"]["mean"] if "zv" in results["q010"][nme]["delta"] else 0.0,
            rec["diff_vs_ref_hm"]["value"], "metrics_test.json")
        mkey = MASTER_KEY_ANALYTIC.get(nme) or next(mm["master"] for mm in MODELS if mm["key"] == runs[nme]["model"])
        mrec = master["groups"][mkey]["per_run"][run_name]["recomputed"]
        chk(run_name, "master mpjpe_hm", L10["mpjpe"], mrec["mpjpe_hm"], "master_aide.json")
        chk(run_name, "master dn_mpjpe_hm", L10["dn_mpjpe"], mrec["dn_mpjpe_hm"], "master_aide.json")
        chk(run_name, "master mpjpe_all", L100["mpjpe"], mrec["mpjpe_all"], "master_aide.json")
        chk(run_name, "master dn_mpjpe", L100["dn_mpjpe"], mrec["dn_mpjpe"], "master_aide.json")
        chk(run_name, "master pck05", L100["pck05"], mrec["pck05"], "master_aide.json")
    n_checked = len(all_names)
    max_chk = max(checks)
    print(f"[hm] integrity: {n_checked} runs × {len(checks) // n_checked} cells, max abs diff {max_chk:.2e}, {len(mism)} mismatch(es)")
    for m_ in mism:
        print("  MISMATCH:", m_)
    assert not mism, "official top-10 % / All cells do not reproduce master_aide.md — aborting (nothing written)"

    # ---- group aggregation: levels (mean ± sd over seeds), Δ per seed -> verdict, seed-mean per-clip CIs for the curves
    group_names = [an["key"] for an in ANALYTIC] + [mm["key"] for mm in MODELS if groups[mm["key"]]]
    glabel = {an["key"]: an["label"] for an in ANALYTIC}
    glabel.update({mm["key"]: mm["label"] for mm in MODELS})
    gruns = {an["key"]: [an["key"]] for an in ANALYTIC}
    gruns.update({mm["key"]: groups[mm["key"]] for mm in MODELS})

    def gsummary(sub, gkey, comps, part=None):
        res = results[sub["key"]]
        rs = gruns[gkey]
        lk = "level" if part is None else f"level_{part}"
        out = {"runs": rs, "n_runs": len(rs),
               "mpjpe": agg([res[r][lk]["mpjpe"] for r in rs]),
               "h": [agg([res[r][lk]["h"][h] for r in rs]) for h in range(P.T_PRED)],
               "dn_mpjpe": agg([res[r][lk]["dn_mpjpe"] for r in rs]),
               "per_run_ci": {r: res[r][lk]["ci"] for r in rs}}
        if part is None:
            out["pck05"] = agg([res[r][lk]["pck05"] for r in rs])
        # CI of the seed-mean level: clip bootstrap of the per-clip error averaged over the group's runs (same indices)
        boot = Boot(sub["n"], a.n_boot, 0)
        cm = np.mean([(ce_all if part is None else ce_part[part])[r][sub["mask"]] for r in rs], axis=0)
        chm = np.mean([(ce_h if part is None else ce_part_h[part])[r][sub["mask"]] for r in rs], axis=0)
        out["seedmean_ci"] = boot.ci(cm)
        out["seedmean_h_ci"] = [boot.ci(chm[:, h]) for h in range(P.T_PRED)]
        out["delta"] = {}
        for comp in comps:
            dk = "" if part is None else f"part_{part}"
            rows = []
            for r in rs:
                if comp in res[r]["delta"]:
                    dd = res[r]["delta"][comp] if part is None else res[r]["delta"][comp][dk]
                    rows.append({**{k: dd[k] for k in ("mean", "ci_lo", "ci_hi", "p_lt0")}, "run": r, "vs_run": res[r]["delta"][comp]["vs_run"],
                                 "h": [{k: hh[k] for k in ("mean", "ci_lo", "ci_hi")} for hh in dd["h"]]})
            if not rows:
                continue
            vd = verdict(rows) if len(rows) >= 3 else ("SINGLE_RULE" if len(rows) == 1 else "INSUFFICIENT")
            dsum = {"per_seed": rows, "n_seeds": len(rows), "n_ci_excl0": n_excl(rows), "verdict": vd,
                    "mean": agg([r_["mean"] for r_ in rows]), "h_mean": [agg([r_["h"][h]["mean"] for r_ in rows]) for h in range(P.T_PRED)],
                    "h_n_ci_excl0": [n_excl([r_["h"][h] for r_ in rows]) for h in range(P.T_PRED)]}
            # seed-mean paired Δ per clip (for the Δ curves): mean over seeds of (model − comparator) per clip, then clip bootstrap
            cr = [rw["vs_run"] for rw in rows]
            dm = np.mean([(ce_h[rw["run"]][sub["mask"]] - ce_h[cr_][sub["mask"]]) for rw, cr_ in zip(rows, cr)], axis=0)
            dsum["seedmean_h"] = [float(dm[:, h].mean()) for h in range(P.T_PRED)]
            dsum["seedmean_h_ci"] = [boot.ci(dm[:, h]) for h in range(P.T_PRED)]
            out["delta"][comp] = dsum
        return out

    summary = {}
    for sub in subsets:
        comps = COMPARATORS if sub["section"] in ("A", "B") else COMPARATORS_C
        summary[sub["key"]] = {g: gsummary(sub, g, comps) for g in group_names}
        if sub["section"] == "B":
            for g in group_names:
                summary[sub["key"]][g][f"part_{sub['part']}"] = gsummary(sub, g, comps, sub["part"])
    print(f"[hm] aggregation done ({time.time() - t0:.1f} s)")

    # ---- figures
    fig_paths = make_figures(a.figs, subsets, summary, group_names, glabel, stamp)

    # ---- markdown
    L = []
    subA = [s for s in subsets if s["section"] == "A"]
    subB = [s for s in subsets if s["section"] == "B"]
    subC = [s for s in subsets if s["section"] == "C"]
    union = next(s for s in subC if s["key"] == "reactive_union")
    classes = [s for s in subC if s["key"] != "reactive_union"]
    n_runs_total = len(runs)
    L += [f"# Rule-based alternatives to the official High-Motion subset — AIDE official test split (n = {N})", "",
          f"Generated {stamp} by `aide/eval/report_hm_subsets.py`; pre-registered in `docs/decisions.md` at {PREREG_UTC} (definitions, model list, "
          f"class names of the reactive union and the reporting rule were fixed before any metric was computed).  Every number is recomputed from the saved "
          f"`pred_test.npz` of every run ({n_runs_total} model runs + 3 analytic rules recomputed from the observed frames) with `aide/eval/metrics.py` on the "
          f"`load_split(\"test\")` ground truth (protocol frames 5..9, all {P.J} joints, no confidence mask).  Locked protocol (`aide/data/protocol.py`, Gate 0): "
          f"driver `{P.DRIVER_RULE}`, grid `{P.SAMPLE_GRID}`, {P.T_OBS} observed → {P.T_PRED} future frames, {P.W}×{P.H} px.  "
          f"**Integrity check: {n_checked} runs, {len(checks)} cells (q = 10 % vs the HM cells and q = 100 % vs the All cells of every `metrics_test.json` and of "
          f"`master_aide.json`), max abs diff {max_chk:.2e}; {len(mism)} mismatch(es)**; analytic rules vs saved predictions max abs diff "
          + ", ".join(f"{k} {v:.1e}" for k, v in rule_check.items()) + f"; shared-index bootstrap vs `metrics.paired_diff` max abs diff {max(boot_check.values()):.1e}.",
          "",
          "## 0. How to read", "",
          "- **All subsets are reported; none was used for selection, adoption or any headline number.**  The official top-10 % HM subset (n = 60, threshold "
          f"{thr_official:.2f} px; the paper's 72.05) is the ONLY subset on which Driver-WM's published HM h=5 = 155.82 applies; every other subset is a rule of ours "
          "and has no published counterpart.",
          "- **Driver-WM reproduction (ours)** has All-MPJPE 130.86 ± 3.11 (master_aide.md) and does NOT reproduce the published geometry (71.47 / 138.03): margins against "
          "it reflect our reproduction, not Driver-WM as published.  The meaningful architecture comparator is **Driver-WM architecture on CockpitWM tokens** "
          "(the published architecture trained on our keypoint part tokens under our formulation; 50.65 ± 0.29 All).",
          "- Levels are mean ± sd over seeds (sample sd, ddof 1) of one test evaluation per run; the bracketed CI is the 10k clip-bootstrap 95 % interval of the "
          "seed-mean per-clip error (the same clip indices for every statistic of a subset).  Δ = paired per-clip difference model − comparator on the subset "
          "(negative = model better), per seed with its own 10k clip-bootstrap CI (seed 0); the table shows mean ± sd over seeds and k/n = seeds whose CI excludes 0.",
          "- Verdict = the existing rule (`report_arch_attribution.verdict`): BEATS / WORSE = CI excl. 0 on >= ceil(2n/3) seeds with the same sign on every seed "
          "(2 of 3, 4 of 5); NOT_WORSE = all CIs cover 0 and abs(Δ) < 0.25 px; else INCONCLUSIVE.  Seed-matched comparators pair the shared seeds 42/123/456 only.  "
          "Analytic comparators (zero-velocity, ZV-shrink) are single deterministic rules: their own Δ rows are labelled SINGLE_RULE and carry only the CI.",
          "- Small subsets (n = 12 … 60) have wide CIs by construction; the per-class subsets of section C are reported in full whatever their size and are "
          "descriptive only.  The MotionBERT row is a supplement-only control.",
          "- PCK@0.05 uses the Gate-0 normaliser (0.05 · max(W, H) = 96 px); d-nMPJPE = 100 · MPJPE / diag = MPJPE / 22.029.", ""]

    # ---------------- section A
    L += ["## 1. Section A — motion-quantile sweep (official motion score, top q of the test clips)", "",
          "Subsets: " + "; ".join(f"top {int(round(s['q'] * 100))} % = n {s['n']} (threshold {s['threshold']:.2f} px, mean score {s['mean_score']:.2f} px)" for s in subA)
          + ".  Motion score = mean GT displacement over the 4 future transitions and all 136 joints (Driver-WM Eq. 16); the top-10 % subset IS the official HM subset "
          "(identical clip set, identical threshold).", "",
          "### 1.1 Owner summary — CockpitWM vs each comparator (horizon-averaged MPJPE on the subset; Δ = CockpitWM − comparator, px)", "",
          "| top q | n clips | comparator | CockpitWM (px) | comparator (px) | Δ mean ± sd over seeds | CI excl. 0 | verdict |", "|---|---|---|---|---|---|---|---|"]
    for s in subA:
        gK = summary[s["key"]]["K"]
        for comp in COMPARATORS:
            if comp not in gK["delta"]:
                continue
            dd = gK["delta"][comp]
            cl = summary[s["key"]][comp]["mpjpe"]
            L.append(f"| {int(round(s['q'] * 100))} % | {s['n']} | {COMP_SHORT[comp]} | {ms(gK['mpjpe']['per_seed'])} | {ms(cl['per_seed'])} | "
                     f"{ms(dd['mean']['per_seed'], sign=True)} | {dd['n_ci_excl0']}/{dd['n_seeds']} | {dd['verdict']} |")
    L += ["", "### 1.2 Levels per model and quantile — horizon-averaged MPJPE (px), mean ± sd over seeds [CI of the seed-mean level]", "",
          "| model | seeds | " + " | ".join(f"top {int(round(s['q'] * 100))} % (n {s['n']})" for s in subA) + " |", "|---|---|" + "---|" * len(subA)]
    for g in group_names:
        cells = []
        for s in subA:
            gg = summary[s["key"]][g]
            cells.append(f"{ms(gg['mpjpe']['per_seed'])} {ci_str(*gg['seedmean_ci'])}")
        L.append(f"| {glabel[g]} | {summary[subA[0]['key']][g]['n_runs']} | " + " | ".join(cells) + " |")
    L += ["", "### 1.3 Per-horizon levels, d-nMPJPE and PCK@0.05 per model and quantile (seed means; sd in the json)", "",
          "| top q | n | model | h=1 | h=2 | h=3 | h=4 | h=5 | avg | d-n | PCK@0.05 |", "|---|---|---|---|---|---|---|---|---|---|---|"]
    for s in subA:
        for g in group_names:
            gg = summary[s["key"]][g]
            L.append(f"| {int(round(s['q'] * 100))} % | {s['n']} | {glabel[g]} | " + " | ".join(fmt(hh["mean"]) for hh in gg["h"])
                     + f" | {fmt(gg['mpjpe']['mean'])} | {fmt(gg['dn_mpjpe']['mean'], 3)} | {fmt(gg['pck05']['mean'])} |")
    L += ["", "### 1.4 Paired Δ per model, quantile and comparator (horizon-averaged and h=5; mean ± sd over seeds; k/n = seeds with CI excl. 0)", "",
          "| top q | n | model | comparator | Δ avg | CI excl. 0 | Δ h=5 | CI excl. 0 (h=5) | verdict (avg) | seed-42 Δ avg [CI] |", "|---|---|---|---|---|---|---|---|---|---|"]
    for s in subA:
        for g in group_names:
            gg = summary[s["key"]][g]
            for comp in COMPARATORS:
                if comp not in gg["delta"]:
                    continue
                dd = gg["delta"][comp]
                r0 = dd["per_seed"][0]
                L.append(f"| {int(round(s['q'] * 100))} % | {s['n']} | {glabel[g]} | {COMP_SHORT[comp]} | {ms(dd['mean']['per_seed'], sign=True)} | "
                         f"{dd['n_ci_excl0']}/{dd['n_seeds']} | {ms(dd['h_mean'][4]['per_seed'], sign=True)} | {dd['h_n_ci_excl0'][4]}/{dd['n_seeds']} | {dd['verdict']} | "
                         f"{fmt(r0['mean'], sign=True)} {ci_str(r0['ci_lo'], r0['ci_hi'], sign=True)} |")
    L += ["", "### 1.5 Curve data — seed-mean MPJPE (px) vs quantile, per horizon (figure `hm_quantile_curves`); CockpitWM rows carry the CI of the seed-mean level", "",
          "| model | horizon | " + " | ".join(f"top {int(round(s['q'] * 100))} %" for s in subA) + " |", "|---|---|" + "---|" * len(subA)]
    for g in group_names:
        for h in range(P.T_PRED):
            cells = []
            for s in subA:
                gg = summary[s["key"]][g]
                c = fmt(gg["h"][h]["mean"])
                if g == "K":
                    c += " " + ci_str(*gg["seedmean_h_ci"][h])
                cells.append(c)
            L.append(f"| {glabel[g]} | h={h + 1} | " + " | ".join(cells) + " |")
    L += ["", "Δ vs zero-velocity curve data (seed-mean paired Δ per clip, px; CockpitWM with CI):", "",
          "| model | horizon | " + " | ".join(f"top {int(round(s['q'] * 100))} %" for s in subA) + " |", "|---|---|" + "---|" * len(subA)]
    for g in group_names:
        if g == "zv":
            continue
        for h in range(P.T_PRED):
            cells = []
            for s in subA:
                dd = summary[s["key"]][g]["delta"]["zv"]
                c = fmt(dd["seedmean_h"][h], sign=True)
                if g == "K":
                    c += " " + ci_str(*dd["seedmean_h_ci"][h], sign=True)
                cells.append(c)
            L.append(f"| {glabel[g]} | h={h + 1} | " + " | ".join(cells) + " |")

    # ---------------- section B
    L += ["", "## 2. Section B — part-specific high motion (top 10 % = n 60 by the motion score of one joint group)", "",
          "| subset | defining joints | n | threshold (px, part score) | mean all-joint score (px) | overlap with the official HM (clips) | pairwise overlap |", "|---|---|---|---|---|---|---|"]
    for s in subB:
        others = [f"{o['part']} {int((s['mask'] & o['mask']).sum())}" for o in subB if o is not s]
        L.append(f"| {s['label']} | {len(PART_JOINTS[s['part']])} | {s['n']} | {s['threshold']:.2f} | {s['mean_score']:.2f} | {s['n_in_official_hm']}/60 | {', '.join(others)} |")
    L += ["", "### 2.1 Owner summary — CockpitWM vs each comparator on each part subset (Δ = CockpitWM − comparator, px; all joints and the defining part)", "",
          "| subset | n | comparator | Δ all joints mean ± sd | CI excl. 0 | verdict (all) | Δ part joints mean ± sd | CI excl. 0 | verdict (part) |", "|---|---|---|---|---|---|---|---|---|"]
    for s in subB:
        gK, gKp = summary[s["key"]]["K"], summary[s["key"]]["K"][f"part_{s['part']}"]
        for comp in COMPARATORS:
            if comp not in gK["delta"]:
                continue
            da, dp = gK["delta"][comp], gKp["delta"][comp]
            L.append(f"| {s['label']} | {s['n']} | {COMP_SHORT[comp]} | {ms(da['mean']['per_seed'], sign=True)} | {da['n_ci_excl0']}/{da['n_seeds']} | {da['verdict']} | "
                     f"{ms(dp['mean']['per_seed'], sign=True)} | {dp['n_ci_excl0']}/{dp['n_seeds']} | {dp['verdict']} |")
    L += ["", "### 2.2 Levels per model and part subset — MPJPE on all joints and on the defining part (px; mean ± sd over seeds [CI of the seed-mean]); per-horizon seed means", "",
          "| subset | model | all joints avg | part avg | all h=1 | all h=2 | all h=3 | all h=4 | all h=5 | part h=1 | part h=2 | part h=3 | part h=4 | part h=5 | PCK@0.05 |",
          "|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|"]
    for s in subB:
        for g in group_names:
            gg, gp = summary[s["key"]][g], summary[s["key"]][g][f"part_{s['part']}"]
            L.append(f"| {s['label']} | {glabel[g]} | {ms(gg['mpjpe']['per_seed'])} {ci_str(*gg['seedmean_ci'])} | {ms(gp['mpjpe']['per_seed'])} {ci_str(*gp['seedmean_ci'])} | "
                     + " | ".join(fmt(hh["mean"]) for hh in gg["h"]) + " | " + " | ".join(fmt(hh["mean"]) for hh in gp["h"]) + f" | {fmt(gg['pck05']['mean'])} |")
    L += ["", "### 2.3 Paired Δ per model, part subset and comparator (all joints and the defining part; horizon-averaged; mean ± sd over seeds)", "",
          "| subset | model | comparator | Δ all avg | CI excl. 0 | verdict | Δ part avg | CI excl. 0 | verdict | Δ all h=5 | Δ part h=5 |", "|---|---|---|---|---|---|---|---|---|---|---|"]
    for s in subB:
        for g in group_names:
            gg, gp = summary[s["key"]][g], summary[s["key"]][g][f"part_{s['part']}"]
            for comp in COMPARATORS:
                if comp not in gg["delta"]:
                    continue
                da, dp = gg["delta"][comp], gp["delta"][comp]
                L.append(f"| {s['label']} | {glabel[g]} | {COMP_SHORT[comp]} | {ms(da['mean']['per_seed'], sign=True)} | {da['n_ci_excl0']}/{da['n_seeds']} | {da['verdict']} | "
                         f"{ms(dp['mean']['per_seed'], sign=True)} | {dp['n_ci_excl0']}/{dp['n_seeds']} | {dp['verdict']} | {ms(da['h_mean'][4]['per_seed'], sign=True)} | {ms(dp['h_mean'][4]['per_seed'], sign=True)} |")

    # ---------------- section C
    L += ["", "## 3. Section C — label-defined subsets (one AIDE label per clip per task; class order = Driver-WM Table 4)", "",
          f"Pre-specified union \"reactive manoeuvre\" = {union['definition']}: n = {union['n']} "
          f"(VCR part {union['n_vcr_part']}, DBR part {union['n_dbr_part']}, both {union['n_both']}); {union['n_in_official_hm']} of the 60 official HM clips fall in it; "
          f"mean motion score {union['mean_score']:.2f} px (all test clips {score_all.mean():.2f} px).  Reason for the class choice (fixed in the pre-registration): the three VCR "
          "classes denote an active vehicle manoeuvre (steering, mirror / shoulder checks) — Forward Moving and Parking are steady states; the two DBR classes denote overt "
          "driver movement — Normal Driving, Talking, Smoking, Making Phone and Dozing Off are sustained postures or non-motor states.", "",
          "### 3.1 Owner summary — CockpitWM vs each comparator on the reactive-manoeuvre union (horizon-averaged MPJPE, px)", "",
          "| subset | n clips | comparator | CockpitWM (px) | comparator (px) | Δ mean ± sd over seeds | CI excl. 0 | verdict |", "|---|---|---|---|---|---|---|---|"]
    gK = summary[union["key"]]["K"]
    for comp in COMPARATORS_C:
        dd = gK["delta"][comp]
        cl = summary[union["key"]][comp]["mpjpe"]
        L.append(f"| reactive manoeuvre (union) | {union['n']} | {COMP_SHORT[comp]} | {ms(gK['mpjpe']['per_seed'])} | {ms(cl['per_seed'])} | "
                 f"{ms(dd['mean']['per_seed'], sign=True)} | {dd['n_ci_excl0']}/{dd['n_seeds']} | {dd['verdict']} |")
    L += ["", "### 3.2 Every class of every task as its own subset — n, GT motion, zero-velocity and CockpitWM levels, CockpitWM paired Δ vs the three comparators (px; verdict in brackets)", "",
          "| task | class | n | in official HM | mean motion score (px) | zero-velocity | CockpitWM | Δ vs ZV | Δ vs Driver-WM arch. on CockpitWM tokens | Δ vs Driver-WM reproduction (ours) |",
          "|---|---|---|---|---|---|---|---|---|---|"]
    for s in classes + [union]:
        gK = summary[s["key"]]["K"]
        cells = []
        for comp in COMPARATORS_C:
            dd = gK["delta"][comp]
            cells.append(f"{ms(dd['mean']['per_seed'], sign=True)} ({dd['n_ci_excl0']}/{dd['n_seeds']}; {dd['verdict']})")
        task = TASK_NAME[s["task"]] if "task" in s else "union"
        cls = s["cls"] if "cls" in s else "reactive manoeuvre"
        L.append(f"| {task} | {cls} | {s['n']} | {s['n_in_official_hm']} | {s['mean_score']:.2f} | {fmt(summary[s['key']]['zv']['mpjpe']['mean'])} | "
                 f"{ms(gK['mpjpe']['per_seed'])} | " + " | ".join(cells) + " |")
    L += ["", "### 3.3 Per-horizon levels per class subset and model (seed means, px)", "",
          "| task | class | n | model | h=1 | h=2 | h=3 | h=4 | h=5 | avg |", "|---|---|---|---|---|---|---|---|---|---|"]
    for s in classes + [union]:
        task = TASK_NAME[s["task"]] if "task" in s else "union"
        cls = s["cls"] if "cls" in s else "reactive manoeuvre"
        for g in group_names:
            gg = summary[s["key"]][g]
            L.append(f"| {task} | {cls} | {s['n']} | {glabel[g]} | " + " | ".join(fmt(hh["mean"]) for hh in gg["h"]) + f" | {ms(gg['mpjpe']['per_seed'])} |")
    L += ["", "### 3.4 Paired Δ per class subset, model and comparator (horizon-averaged; mean ± sd over seeds; k/n = seeds with CI excl. 0)", "",
          "| task | class | n | model | comparator | Δ avg | CI excl. 0 | verdict | Δ h=5 |", "|---|---|---|---|---|---|---|---|---|"]
    for s in classes + [union]:
        task = TASK_NAME[s["task"]] if "task" in s else "union"
        cls = s["cls"] if "cls" in s else "reactive manoeuvre"
        for g in group_names:
            gg = summary[s["key"]][g]
            for comp in COMPARATORS_C:
                if comp not in gg["delta"]:
                    continue
                dd = gg["delta"][comp]
                L.append(f"| {task} | {cls} | {s['n']} | {glabel[g]} | {COMP_SHORT[comp]} | {ms(dd['mean']['per_seed'], sign=True)} | {dd['n_ci_excl0']}/{dd['n_seeds']} | "
                         f"{dd['verdict']} | {ms(dd['h_mean'][4]['per_seed'], sign=True)} |")

    # ---------------- provenance
    L += ["", "## 4. Provenance and definitions", "",
          "| group | runs (aide/results/runs/) | seeds |", "|---|---|---|"]
    for mm in MODELS:
        L.append(f"| {mm['label']} | " + ", ".join(f"`{r}`" for r in groups[mm["key"]]) + " | " + "/".join(str(runs[r]["seed"]) for r in groups[mm["key"]]) + " |")
    for an in ANALYTIC:
        L.append(f"| {an['label']} | recomputed from the observed frames (`baselines_analytic.{ 'zero_velocity' if an['key'] == 'zv' else ('zv_shrink' if an['key'] == 'zvs' else 'zv_shrink_h') }`), "
                 f"bit-identical to `{an['run']}/pred_test.npz` | — |")
    L += ["",
          f"- Motion score (Driver-WM Eq. 16, `metrics.motion_score`): mean over the {P.HM_N_TRANSITIONS} future frame-to-frame transitions and the scored joints of the Euclidean GT "
          f"displacement (px).  Section A / C use all {P.J} joints; section B restricts the score to hands = joints 94–135, face = 26–93, body = 0–25 (`protocol.PARTS`).  "
          "Ranking: stable descending argsort; ties impossible at float64 resolution here.",
          f"- Subset sizes: n = floor(q · {N}) → " + ", ".join(f"{int(round(s['q'] * 100))} % = {s['n']}" for s in subA) + "; part subsets n = floor(0.1 · 609) = 60.",
          "- Class subsets: `y_<task>` of the keypoint cache (`load_split(\"test\")`), class index = position in `protocol.CLASSES` (DBR: Normal Driving, Looking Around, Making Phone, "
          "Body Movement, Talking, Smoking, Dozing Off; DER: Peace, Anxiety, Weariness, Happiness, Anger; TCR: Smooth Traffic, Waiting, Traffic Jam; VCR: Forward Moving, Parking, "
          "Turning, Backward Moving, Changing Lane).",
          f"- ZV-shrink: (1 − α)·last + α·mean of the 5 observed poses with α = {fits['zvs']:.2f} (val-fitted, `zv_shrink_val/alpha_fit_val.json`); per-horizon variant α_h = "
          f"{fits['zvsh']:.2f}·h/5 (`zv_shrink_h_val/alpha_fit_val.json`); both recomputed here from frames 0..4 and cast to float32 as the saved rule predictions were.",
          f"- Bootstrap: {a.n_boot} clip resamples, `numpy.random.default_rng(0)` indices shared by every statistic of a subset (p_lt0 uses seed 1), percentile 2.5 / 97.5; "
          "verified bit-identical to `metrics.paired_diff` on the official HM cell of CockpitWM seed 42 vs zero-velocity.",
          f"- Figures: `{fig_paths['curves_pdf']}` / `.png` (rows: seed-mean MPJPE vs quantile and seed-mean paired Δ vs zero-velocity; one panel per horizon; shaded = CockpitWM "
          f"clip-bootstrap CI) and `{fig_paths['parts']}` (all-joint and defining-part MPJPE per model on the three part subsets, CI of the seed-mean level).",
          f"- Wall time {time.time() - t0:.1f} s, CPU only."]

    a.out.mkdir(parents=True, exist_ok=True)
    out_md, out_js = a.out / "hm_subsets.md", a.out / "hm_subsets.json"
    out_md.write_text("\n".join(L) + "\n", encoding="utf-8")

    def sub_meta(s):
        keep = {k: v for k, v in s.items() if k != "mask"}
        keep["clip_ids"] = d["clip_id"][s["mask"]].astype(str).tolist()
        return keep

    js = {"generated_utc": stamp, "preregistration_utc": PREREG_UTC, "script": str(Path(__file__).resolve().relative_to(ROOT)), "n_test": N, "n_boot": a.n_boot, "tol": a.tol,
          "protocol": {"driver_rule": P.DRIVER_RULE, "sample_grid": P.SAMPLE_GRID, "hm_official_n": int(hm_official.sum()), "hm_official_threshold_px": thr_official,
                       "hm_n_transitions": P.HM_N_TRANSITIONS, "pck05_threshold_px": thr_pck, "diag_px": P.DIAG, "zv_shrink_alpha": fits["zvs"], "zv_shrink_h_a": fits["zvsh"]},
          "integrity": {"n_runs": n_checked, "n_cells": len(checks), "max_abs_diff": max_chk, "mismatches": mism, "analytic_rule_vs_saved_max_abs_diff": rule_check,
                        "bootstrap_vs_metrics_paired_diff_max_abs_diff": max(boot_check.values())},
          "models": [{k: v for k, v in mm.items()} for mm in MODELS], "analytic": ANALYTIC, "groups": {g: gruns[g] for g in group_names}, "labels": glabel,
          "comparators": {"AB": COMPARATORS, "C": COMPARATORS_C, "labels": COMP_LABEL}, "missing_runs": missing,
          "subsets": [sub_meta(s) for s in subsets],
          "summary": summary,          # per subset -> per group: levels (mean/sd/per_seed, seed-mean CI), Δ per comparator (per seed with CI, verdict)
          "per_run": {sk: {nme: {kk: vv for kk, vv in r.items()} for nme, r in res.items()} for sk, res in results.items()},
          "figures": fig_paths, "wall_s": round(time.time() - t0, 1)}
    out_js.write_text(json.dumps(js, indent=1, default=lambda o: o.tolist() if hasattr(o, "tolist") else str(o)), encoding="utf-8")
    print(f"[hm] wrote {out_md}\n[hm] wrote {out_js}  ({time.time() - t0:.1f} s)")
    return 0


# ------------------------------------------------------------------------------------------------------------ figures
# categorical palette (dataviz skill reference palette; validated adjacent-pair CVD/normal-vision floors, light mode): fixed slot order
COLORS = {"K": "#2a78d6", "DWMkp": "#eb6834", "Kpre": "#1baf7a", "KS": "#4a3aa7", "DWMrepro": "#e34948"}
MARKERS = {"K": "o", "DWMkp": "s", "Kpre": "^", "KS": "D", "DWMrepro": "v"}
NEUTRAL = {"zv": ("#52514e", "-"), "zvs": ("#52514e", ":"), "zvsh": ("#52514e", "-.")}
FIG_ORDER = ["K", "Kpre", "KS", "DWMkp", "DWMrepro"]
FIG_LABEL = {"K": "CockpitWM", "Kpre": "CockpitWM + BATON pretraining", "KS": "CockpitWM (+ semantic heads)",
             "DWMkp": "Driver-WM architecture on CockpitWM tokens", "DWMrepro": "Driver-WM reproduction (ours)",
             "zv": "zero-velocity", "zvs": "ZV-shrink α = 0.44", "zvsh": "ZV-shrink α_h = 0.76·h/5"}


def make_figures(figdir: Path, subsets, summary, group_names, glabel, stamp) -> dict:
    import matplotlib
    matplotlib.use("Agg")
    import matplotlib.pyplot as plt

    figdir.mkdir(parents=True, exist_ok=True)
    plt.rcParams.update({"font.size": 8, "axes.spines.top": False, "axes.spines.right": False, "axes.grid": True, "grid.color": "#e6e5e1",
                         "grid.linewidth": 0.6, "axes.edgecolor": "#9a9891", "axes.labelcolor": "#0b0b0b", "xtick.color": "#52514e", "ytick.color": "#52514e",
                         "legend.frameon": False, "pdf.fonttype": 42})
    subA = [s for s in subsets if s["section"] == "A"]
    qs = [int(round(s["q"] * 100)) for s in subA]
    x = np.arange(len(qs))
    fig, axes = plt.subplots(2, P.T_PRED, figsize=(13, 5.6), sharex=True)
    for h in range(P.T_PRED):
        ax = axes[0, h]
        for g in ["zv", "zvs", "zvsh"]:
            y = [summary[s["key"]][g]["h"][h]["mean"] for s in subA]
            c, ls = NEUTRAL[g]
            ax.plot(x, y, ls, color=c, lw=1.6, label=FIG_LABEL[g], zorder=2)
        for g in FIG_ORDER:
            if g not in group_names:
                continue
            y = [summary[s["key"]][g]["h"][h]["mean"] for s in subA]
            ax.plot(x, y, "-", color=COLORS[g], lw=2, marker=MARKERS[g], ms=4.5, mec="white", mew=0.8, label=FIG_LABEL[g], zorder=3)
            if g == "K":
                lo = [summary[s["key"]][g]["seedmean_h_ci"][h][0] for s in subA]
                hi = [summary[s["key"]][g]["seedmean_h_ci"][h][1] for s in subA]
                ax.fill_between(x, lo, hi, color=COLORS[g], alpha=0.15, lw=0, zorder=1)
        ax.set_title(f"horizon h = {h + 1}", fontsize=9)
        if h == 0:
            ax.set_ylabel("MPJPE on the subset (px)")
        ax2 = axes[1, h]
        ax2.axhline(0, color="#9a9891", lw=0.8, zorder=1)
        for g in ["zvs", "zvsh"]:
            y = [summary[s["key"]][g]["delta"]["zv"]["seedmean_h"][h] for s in subA]
            c, ls = NEUTRAL[g]
            ax2.plot(x, y, ls, color=c, lw=1.6, zorder=2)
        for g in FIG_ORDER:
            if g not in group_names:
                continue
            y = [summary[s["key"]][g]["delta"]["zv"]["seedmean_h"][h] for s in subA]
            if g == "DWMrepro":     # far off-scale (+37 … +89 px): stated as text so the 1–20 px structure of the other rows stays legible
                ax2.text(0.98, 0.03, f"reproduction: Δ +{min(y):.0f} … +{max(y):.0f} px (off-scale)", transform=ax2.transAxes, va="bottom", ha="right",
                         fontsize=7, color=COLORS[g])
                continue
            ax2.plot(x, y, "-", color=COLORS[g], lw=2, marker=MARKERS[g], ms=4.5, mec="white", mew=0.8, zorder=3)
            if g == "K":
                lo = [summary[s["key"]][g]["delta"]["zv"]["seedmean_h_ci"][h][0] for s in subA]
                hi = [summary[s["key"]][g]["delta"]["zv"]["seedmean_h_ci"][h][1] for s in subA]
                ax2.fill_between(x, lo, hi, color=COLORS[g], alpha=0.15, lw=0, zorder=1)
        ax2.set_xticks(x)
        ax2.set_xticklabels([f"{q}" for q in qs])
        ax2.set_xlabel("top q % of test clips by GT motion")
        if h == 0:
            ax2.set_ylabel("paired Δ vs zero-velocity (px)")
    for h in range(P.T_PRED):
        axes[0, h].set_ylim(bottom=0)
    handles, labels_ = axes[0, 0].get_legend_handles_labels()
    fig.legend(handles, labels_, loc="lower center", ncol=4, fontsize=8, bbox_to_anchor=(0.5, -0.02))
    fig.suptitle("AIDE test: MPJPE vs. GT-motion quantile (top) and paired Δ vs. zero-velocity (bottom; reproduction off-scale). Shaded: CockpitWM 10k clip-bootstrap CI.\n"
                 f"Top 10 % = the official HM subset (n = 60). Descriptive; pre-registered {PREREG_UTC}.", fontsize=8, y=1.0)
    fig.tight_layout(rect=(0, 0.07, 1, 0.95))
    p_pdf, p_png = figdir / "hm_quantile_curves.pdf", figdir / "hm_quantile_curves.png"
    fig.savefig(p_pdf); fig.savefig(p_png, dpi=200)
    plt.close(fig)

    # ---- part subsets: dot plot per part subset; filled = all joints, hollow = defining part; CI of the seed-mean level
    subB = [s for s in subsets if s["section"] == "B"]
    order = ["zv", "zvs", "zvsh"] + [g for g in FIG_ORDER if g in group_names]
    fig, axes = plt.subplots(1, len(subB), figsize=(13, 3.6), sharey=True)
    for k, s in enumerate(subB):
        ax = axes[k]
        yy = np.arange(len(order))[::-1]
        for yi, g in zip(yy, order):
            gg, gp = summary[s["key"]][g], summary[s["key"]][g][f"part_{s['part']}"]
            col = COLORS.get(g, "#52514e")
            ax.errorbar(gg["mpjpe"]["mean"], yi + 0.16, xerr=[[gg["mpjpe"]["mean"] - gg["seedmean_ci"][0]], [gg["seedmean_ci"][1] - gg["mpjpe"]["mean"]]],
                        fmt="o", color=col, ms=6, mec="white", mew=0.8, capsize=2, lw=1.2, zorder=3)
            ax.errorbar(gp["mpjpe"]["mean"], yi - 0.16, xerr=[[gp["mpjpe"]["mean"] - gp["seedmean_ci"][0]], [gp["seedmean_ci"][1] - gp["mpjpe"]["mean"]]],
                        fmt="D", mfc="white", color=col, ms=5.5, mew=1.4, capsize=2, lw=1.2, zorder=3)
        ax.set_yticks(yy)
        if k == 0:
            ax.set_yticklabels([FIG_LABEL[g] for g in order])
        ax.set_title(f"top 10 % by {s['label']}\n(n = {s['n']}, {s['n_in_official_hm']}/60 shared with the official HM)", fontsize=8.5)
        ax.set_xlabel("MPJPE (px): ● all joints   ◇ defining part")
        ax.axvline(summary[s["key"]]["zv"]["mpjpe"]["mean"], color="#9a9891", lw=0.8, ls="--", zorder=1)
    fig.suptitle(f"AIDE test, part-specific high-motion subsets — seed-mean levels with 10k clip-bootstrap CI; descriptive, pre-registered {PREREG_UTC}", fontsize=8, y=1.02)
    fig.tight_layout()
    p_parts = figdir / "hm_part_subsets.png"
    fig.savefig(p_parts, dpi=200, bbox_inches="tight")
    plt.close(fig)
    rel = lambda p: str(p.relative_to(ROOT))  # noqa: E731
    return {"curves_pdf": rel(p_pdf), "curves_png": rel(p_png), "parts": rel(p_parts)}


if __name__ == "__main__":
    raise SystemExit(main())
Xet Storage Details

Size:: 52.9 kB
Xet hash:: 3c24ef4959667d94ba57a57c45828232e20033c302483210f62761421ad7ed30
Xet efficiently stores files, intelligently splitting them into unique chunks and accelerating uploads and downloads. More info.