AASHTO VTC enrichment: per-crash victim-type matrix for 2023-2025

Extends `agg.py`'s VTC scaffolding (25 cols: 5 victim types × 5 conditions)
to AASHTO years by computing per-crash counts from the person supplements.

- `supplement_with_sp.py`: adds `compute_vtc()` helper + person-supplement
  inputs (`-O`/`-P`). Merges 25 VTC cols into `aashto_supplemented_crashes`;
  NJSP-supplement rows land in `uf`. NJSP-residual county list updated to
  match observed lag pattern (Monmouth/Mercer/Essex).
- `to_njdot_persons.py`: ghost-Driver rows now get `pos=0` (routes to 'u'
  bucket via VTC aggregation) instead of being dropped. Preserves the
  fatal count while avoiding driver-fatal inflation. Ped/cyclist remain
  under-reported at person-level (mostly land in 'u'); crash-level `pk`
  stays authoritative.
- DVX wiring: `aashto_supplemented_crashes.parquet.dvc` now declares the
  person supplements as deps; `ymccmc.dvc`/`ymccmcs.dvc` updated md5.
- Regenerated `cmymc.db` + `ys`/`yms`/`yccs`/`ymccs`/`ymccmc(s)/` aggs.

Known: VTC fatal sum is ~5% over crash-level `tk` for AASHTO years due
to fatal-tagged persons appearing in non-fatal-crash records (~140 cases
across 2023-2025). Documented in `.dvc` notes.

Author: ryan-williams

njdot/aashto/supplement_with_sp.py

@@ -3,21 +3,30 @@
 # requires-python = ">=3.11"
 # dependencies = ["click", "pandas", "pyarrow"]
 # ///
-"""Supplement AASHTO crashes with NJSP-only fatals (AASHTO ingestion lag).
+"""Supplement AASHTO crashes with NJSP-only fatals (AASHTO ingestion lag)
+and per-crash VTC (victim-type × condition) counts.
 
 AASHTO 2025 has a fatal-classification lag concentrated in specific
-counties (Middlesex, Hudson, Essex etc.) — fatal status takes time to
+counties (Monmouth, Mercer, Essex etc.) — fatal status takes time to
 propagate (death-cert / 30-day rule), while injury/PDO records arrive
 on schedule. NJSP carries the up-to-date fatal classifications, so we
 back-fill the 111 NJSP-only-2025 fatals into the AASHTO output to keep
 the homepage NJDOT plot honest until AASHTO catches up.
 
+Also computes 25-column VTC matrix (`df`/`ds`/.../`un`) from the AASHTO
+person supplements so `agg.py` can emit per-VT aggregations for AASHTO
+years, matching the legacy master crashes parquet. NJSP-supplement rows
+land in the `uf` bucket (driver/passenger/ped/cyclist breakdown not
+available at residual granularity).
+
 Removable: when a fresh `Crash.csv` shows AASHTO in sync with NJSP,
 delete this step from the pipeline.
 
 Reads:
   - njdot/data/aashto_combined_crashes.parquet (pure AASHTO, schema-mapped)
   - njsp/data/njsp_njdot_residuals.parquet     (NJSP-side residuals from match_njdot)
+  - njdot/data/aashto_supplemented_occupants.parquet  (for VTC, optional)
+  - njdot/data/aashto_supplemented_pedestrians.parquet (for VTC, optional)
 
 Writes:
   - njdot/data/aashto_supplemented_crashes.parquet
@@ -31,15 +40,69 @@
 
 err = partial(print, file=sys.stderr)
 
+# VTC matrix: 5 victim types × 5 conditions = 25 columns
+VICTIM_TYPES = ['d', 'o', 'p', 'b', 'u']  # driver, passenger, pedestrian, bicyclist, unknown
+CONDITIONS = ['f', 's', 'm', 'p', 'n']     # fatal, serious, minor, possible, none
+VTC_COLS = [f'{vt}{c}' for vt in VICTIM_TYPES for c in CONDITIONS]
+CONDITION_MAP = {1: 'f', 2: 's', 3: 'm', 4: 'p', 5: 'n', 0: 'n'}
+CRASH_PK = ['year', 'cc', 'mc', 'case']
+
+
+def _pos_to_vt(pos):
+    if pd.isna(pos) or pos == 0:
+        return 'u'
+    return 'd' if pos == 1 else 'o'
+
+
+def compute_vtc(occupants: pd.DataFrame, pedestrians: pd.DataFrame) -> pd.DataFrame:
+    """Compute 25-col VTC matrix aggregated by (year, cc, mc, case).
+
+    Inputs are AASHTO-supplemented O/P frames with `condition` (Int 1-5),
+    `pos` (occupants, 1=driver / 2-12=passenger), `cyclist` (peds).
+    Returns a DataFrame indexed by CRASH_PK with VTC_COLS columns (int).
+    """
+    o = occupants.copy()
+    o['cond'] = o['condition'].map(CONDITION_MAP).fillna('n')
+    o['vt'] = o['pos'].apply(_pos_to_vt)
+    o['vtc'] = o['vt'] + o['cond']
+
+    p = pedestrians.copy()
+    p['cond'] = p['condition'].map(CONDITION_MAP).fillna('n')
+    p['vt'] = p['cyclist'].apply(lambda b: 'b' if b else 'p')
+    p['vtc'] = p['vt'] + p['cond']
+
+    parts = []
+    for df in (o, p):
+        if not len(df):
+            continue
+        agg = df.groupby(CRASH_PK + ['vtc']).size().unstack(fill_value=0)
+        parts.append(agg)
+    if not parts:
+        return pd.DataFrame(columns=CRASH_PK + VTC_COLS).set_index(CRASH_PK)
+
+    combined = parts[0]
+    for extra in parts[1:]:
+        combined = combined.add(extra, fill_value=0)
+    for col in VTC_COLS:
+        if col not in combined.columns:
+            combined[col] = 0
+    combined = combined[VTC_COLS].fillna(0).astype(int)
+    return combined
+
 
 @click.command()
 @click.option("-a", "--aashto", type=click.Path(path_type=Path),
               default=Path("njdot/data/aashto_combined_crashes.parquet"))
 @click.option("-r", "--residuals", type=click.Path(path_type=Path),
               default=Path("njsp/data/njsp_njdot_residuals.parquet"))
+@click.option("-O", "--occupants-supplement", type=click.Path(path_type=Path),
+              default=Path("njdot/data/aashto_supplemented_occupants.parquet"))
+@click.option("-P", "--pedestrians-supplement", type=click.Path(path_type=Path),
+              default=Path("njdot/data/aashto_supplemented_pedestrians.parquet"))
 @click.option("-o", "--out", type=click.Path(path_type=Path),
               default=Path("njdot/data/aashto_supplemented_crashes.parquet"))
-def main(aashto: Path, residuals: Path, out: Path):
+def main(aashto: Path, residuals: Path, occupants_supplement: Path,
+         pedestrians_supplement: Path, out: Path):
     a = pd.read_parquet(aashto)
     aashto_years = sorted(a["year"].dropna().astype(int).unique())
     err(f"AASHTO: {len(a):,} crashes, years {aashto_years[0]}–{aashto_years[-1]}")
@@ -48,41 +111,64 @@ def main(aashto: Path, residuals: Path, out: Path):
     sp_only = res[(res["side"] == "njsp") & res["year"].isin(aashto_years)].copy()
     err(f"NJSP-side residuals in AASHTO years: {len(sp_only)} "
         f"({sp_only['tk'].sum()} deaths)")
-    if not len(sp_only):
-        a.to_parquet(out, index=False)
-        err(f"No supplement needed; copied AASHTO → {out}")
-        return
-
-    # Build AASHTO-schema rows from NJSP residuals. Only fields that matter
-    # for downstream `agg.py` (year, cc, mc, severity, tk, ti, pk, pi, tv).
-    # `case` carries an `NJSP-{i}` synthetic id for provenance.
-    supp = pd.DataFrame({
-        "year": sp_only["year"].astype("int32").values,
-        "cc": sp_only["cc"].astype("Int8").values,
-        "mc": sp_only["mc"].astype("Int16").values,
-        "case": [f"NJSP-supplement-{i}" for i in range(len(sp_only))],
-        # Use date as datetime (no time-of-day in residuals)
-        "dt": pd.to_datetime(sp_only["date"]).values,
-        "severity": pd.array(["f"] * len(sp_only), dtype="string"),
-        "tk": sp_only["tk"].astype("int8").values,
-        "tk_broad": sp_only["tk"].astype("int8").values,
-        "ti": pd.array([0] * len(sp_only), dtype="int8"),
-        "pk": pd.array([0] * len(sp_only), dtype="int8"),
-        "pi": pd.array([0] * len(sp_only), dtype="int8"),
-        "tv": pd.array([0] * len(sp_only), dtype="int8"),
-        "cc0": sp_only["cc"].astype("Int8").values,
-        "mc0": sp_only["mc"].astype("Int16").values,
-    })
-    # Add any AASHTO columns not yet set, as nulls.
-    for col in a.columns:
-        if col not in supp.columns:
-            supp[col] = pd.NA
-    supp = supp[a.columns]
-
-    out_df = pd.concat([a, supp], ignore_index=True)
+    if len(sp_only):
+        # Build AASHTO-schema rows from NJSP residuals. Only fields that matter
+        # for downstream `agg.py` (year, cc, mc, severity, tk, ti, pk, pi, tv).
+        # `case` carries an `NJSP-{i}` synthetic id for provenance.
+        supp = pd.DataFrame({
+            "year": sp_only["year"].astype("int32").values,
+            "cc": sp_only["cc"].astype("Int8").values,
+            "mc": sp_only["mc"].astype("Int16").values,
+            "case": [f"NJSP-supplement-{i}" for i in range(len(sp_only))],
+            # Use date as datetime (no time-of-day in residuals)
+            "dt": pd.to_datetime(sp_only["date"]).values,
+            "severity": pd.array(["f"] * len(sp_only), dtype="string"),
+            "tk": sp_only["tk"].astype("int8").values,
+            "tk_broad": sp_only["tk"].astype("int8").values,
+            "ti": pd.array([0] * len(sp_only), dtype="int8"),
+            "pk": pd.array([0] * len(sp_only), dtype="int8"),
+            "pi": pd.array([0] * len(sp_only), dtype="int8"),
+            "tv": pd.array([0] * len(sp_only), dtype="int8"),
+            "cc0": sp_only["cc"].astype("Int8").values,
+            "mc0": sp_only["mc"].astype("Int16").values,
+        })
+        # Add any AASHTO columns not yet set, as nulls.
+        for col in a.columns:
+            if col not in supp.columns:
+                supp[col] = pd.NA
+        supp = supp[a.columns]
+        out_df = pd.concat([a, supp], ignore_index=True)
+        n_supplemented = len(supp)
+    else:
+        out_df = a.copy()
+        n_supplemented = 0
+
+    # VTC enrichment from person supplements
+    if occupants_supplement.exists() and pedestrians_supplement.exists():
+        err(f"\nComputing VTC from {occupants_supplement.name} + {pedestrians_supplement.name}…")
+        occ = pd.read_parquet(occupants_supplement)
+        ped = pd.read_parquet(pedestrians_supplement)
+        err(f"  loaded {len(occ):,} occupants + {len(ped):,} pedestrians")
+        vtc = compute_vtc(occ, ped).reset_index()
+        # Align dtypes for merge (out_df has Int8/Int16 for cc/mc; vtc has Int8/Int16 too via supplements)
+        out_df = out_df.merge(vtc, on=CRASH_PK, how='left')
+        for col in VTC_COLS:
+            out_df[col] = out_df[col].fillna(0).astype('int8')
+        # NJSP-supplemented rows: each row's tk into `uf` (we don't know VT
+        # breakdown for residual fatals). cc/mc nullable issues left as-is.
+        is_njsp_supp = out_df['case'].astype(str).str.startswith('NJSP-supplement-')
+        if is_njsp_supp.any():
+            out_df.loc[is_njsp_supp, 'uf'] = out_df.loc[is_njsp_supp, 'tk'].astype('int8')
+            err(f"  marked {is_njsp_supp.sum()} NJSP-supplement rows with uf={out_df.loc[is_njsp_supp, 'uf'].sum()}")
+        err(f"  total VTC fatals (df+of+pf+bf+uf): "
+            f"{int(out_df[['df','of','pf','bf','uf']].sum().sum()):,}")
+    else:
+        err(f"  Person supplements not found; VTC columns omitted "
+            f"({occupants_supplement.name}, {pedestrians_supplement.name})")
+
     out.parent.mkdir(parents=True, exist_ok=True)
     out_df.to_parquet(out, index=False)
-    err(f"Wrote {out} ({len(out_df):,} crashes; +{len(supp)} from NJSP)")
+    err(f"Wrote {out} ({len(out_df):,} crashes; +{n_supplemented} from NJSP)")
 
     # Quick verification: post-supplement fatals per year
     fa = out_df[out_df["severity"] == "f"]

njdot/aashto/to_njdot_persons.py

@@ -31,11 +31,15 @@
     Vehicle` blank. ~45K/yr; ~240/yr carry `Fatal Injury (K)`. These
     are placeholder rows — for fatal-ped crashes the actual victim is
     on a separate row with `Person Type='Pedestrian'` and
-    `Severity Rating='No Apparent Injury (O)'`. We drop ghost-Drivers
-    from the occupants supplement entirely. Net effect: 2024+ ped/
-    cyclist fatal counts are *under*-reported in the person-level
-    breakdown (~19/yr instead of ~228/yr). Crash-level `pk` in
-    `aashto_supplemented_crashes.parquet` remains authoritative.
+    `Severity Rating='No Apparent Injury (O)'`. We reclassify these
+    rows with `pos=0` so VTC aggregation routes them to the 'u'
+    (unknown) bucket — preserves the death count while not inflating
+    driver totals. Net effect: 2024+ ped/cyclist breakdown remains
+    *under*-reported by ~228/yr at the person-level (those land in
+    `u` bucket), but the strict total killed (`tk`/`uf`+`df`+...)
+    matches crash-level. Crash-level `pk` in
+    `aashto_supplemented_crashes.parquet` remains the authoritative
+    pedestrian count.
 """
 import sys
 from functools import partial
@@ -116,15 +120,16 @@ def normalize_age(age_series: pd.Series) -> pd.Series:
 def to_occupants(joined: pd.DataFrame, year: int) -> pd.DataFrame:
     """Convert AASHTO Driver/Passenger rows to DOTr-style occupants frame.
 
-    Drops "ghost-Driver" rows (Driver Person Type + blank Position in
-    Vehicle) — see module docstring."""
+    "Ghost-Driver" rows (Driver Person Type + blank Position in Vehicle)
+    get `pos=0` so they aggregate to the 'u' (unknown-VT) bucket — see
+    module docstring."""
     occ_mask = joined['Person Type'].isin(['Driver', 'Passenger'])
-    pos_blank = joined['Position in Vehicle'].isna() | (joined['Position in Vehicle'] == 'None')
-    ghost_drivers = (joined['Person Type'] == 'Driver') & pos_blank
+    o = joined[occ_mask].copy()
+    pos_blank = o['Position in Vehicle'].isna() | (o['Position in Vehicle'] == 'None')
+    ghost_drivers = (o['Person Type'] == 'Driver') & pos_blank
     n_ghost = ghost_drivers.sum()
     if n_ghost:
-        err(f'        dropping {n_ghost:,} ghost-Driver rows (blank Position)')
-    o = joined[occ_mask & ~ghost_drivers].copy()
+        err(f'        flagging {n_ghost:,} ghost-Driver rows as pos=0 (unknown VT)')
 
     o['year'] = pd.Series(year, index=o.index, dtype='int32')
     # vehicle_index is 0-based in AASHTO; vn is 1-based in DOTr
@@ -135,8 +140,12 @@ def to_occupants(joined: pd.DataFrame, year: int) -> pd.DataFrame:
 
     o['condition'] = o['Severity Rating (Person)'].map(SEVERITY_MAP).astype('Int8')
     o['pos'] = o['Position in Vehicle'].map(POSITION_MAP).astype('Int8')
-    # Driver Person Type forces pos=1 (covers AASHTO 'Driver' position field being blank for some rows)
-    o.loc[o['Person Type'] == 'Driver', 'pos'] = 1
+    # Driver Person Type with a *real* Position field gets pos=1 (handles
+    # cases where AASHTO Position='Driver' string vs Person Type field
+    # disagree). Ghost-Drivers (blank Position) stay at pos=0 → 'u' bucket.
+    real_driver = (o['Person Type'] == 'Driver') & ~ghost_drivers
+    o.loc[real_driver, 'pos'] = 1
+    o.loc[ghost_drivers, 'pos'] = 0
     o['eject'] = o['Ejection'].map(EJECTION_MAP).astype('Int8')
     o['age'] = normalize_age(o['Age'])
     o['sex'] = o['Sex'].map(SEX_MAP).fillna('')

njdot/data/aashto_supplemented_crashes.parquet.dvc

@@ -1,19 +1,27 @@
 outs:
-- md5: 2cd7616a654b357444d08aa05c41b311
-  size: 36111389
+- md5: 92ac7e301eab0f949748c6061317cf9f
+  size: 37286075
   hash: md5
   path: aashto_supplemented_crashes.parquet
 meta:
   computation:
     cmd: njdot/aashto/supplement_with_sp.py
     git_deps:
-      /njdot/aashto/supplement_with_sp.py: 4f25ec341452526bf5181ed876069da555e3a863
+      /njdot/aashto/supplement_with_sp.py: e88a1d8b0d3f43c460f23f76c97df7e557da19f5
     deps:
       /njdot/data/aashto_combined_crashes.parquet: 0ae2595cd33ce71d47ab2d022a22b60e
       /njsp/data/njsp_njdot_residuals.parquet: e50b477c45576b6357ed94060e7b669e
+      /njdot/data/aashto_supplemented_occupants.parquet: 806dd7a46ebafbe426faeb05373359b9
+      /njdot/data/aashto_supplemented_pedestrians.parquet: d3d913d7d0716ee4b5588e8380bb456a
     note: |
-      AASHTO + NJSP-only-fatal supplement. Backfills the NJSP fatals
-      that AASHTO hasn't ingested yet (concentrated in Middlesex /
-      Hudson / Essex 2025 — fatal-classification lag while
-      injury/PDO records are on schedule). Removable when a fresh
-      `Crash.csv` shows AASHTO in sync with NJSP.
+      AASHTO + NJSP-only-fatal supplement + per-crash VTC matrix
+      (5 victim types x 5 conditions). Backfills the NJSP fatals
+      that AASHTO hasn't ingested yet (concentrated in Monmouth /
+      Mercer / Essex 2025 - fatal-classification lag while
+      injury/PDO records are on schedule). NJSP-supplement rows
+      land in the `uf` bucket. Ghost-Driver rows from AASHTO also
+      land in `u*` buckets (Person Type=Driver with blank Position,
+      ~45K/yr, mostly mis-tagged ped/cyclist victims). VTC fatal
+      sum is ~5% over crash-level `tk` due to AASHTO data quality
+      (fatal-tagged persons in non-fatal-crash records).
+      Removable when a fresh `Crash.csv` shows AASHTO in sync with NJSP.

njdot/data/aashto_supplemented_occupants.parquet.dvc

@@ -1,13 +1,13 @@
 outs:
-- md5: 6957204e8f53b2fbde29b0f9775cf33f
-  size: 13286278
+- md5: 806dd7a46ebafbe426faeb05373359b9
+  size: 14335377
   hash: md5
   path: aashto_supplemented_occupants.parquet
 meta:
   computation:
     cmd: python njdot/aashto/to_njdot_persons.py -y 2023,2024,2025
     git_deps:
-      /njdot/aashto/to_njdot_persons.py: fecfc0545896518b3047c8ffb0c1d4c054fbd596
+      /njdot/aashto/to_njdot_persons.py: 2f0d9eb2505089b276cd82eff4d3d2b77f34f615
       /njdot/aashto/to_njdot_schema.py: 0ad39efb3151ef4b58b4fd0d4c9bd8900074fce6
     deps:
       /njdot/data/2023/persons.parquet: 72625f5a685ceb7489789e886f4f3171
@@ -18,10 +18,11 @@ meta:
       /njdot/data/2025/crashes.parquet: 26e6bc7f841bbaa6b5878e57771869dd
       /www/public/njdot/cc2mc2mn.json: 2c12563c5ea3a3fdd443fc14610822e8
     note: |
-      AASHTO persons → DOTr-style occupants supplement for 2023-2025.
-      Drops "ghost-Driver" placeholder rows (Person Type=Driver +
-      blank Position in Vehicle). Pedestrian/cyclist fatals are
-      under-reported here; crash-level pk in
-      `aashto_supplemented_crashes.parquet` remains authoritative.
-      Produced alongside `aashto_supplemented_pedestrians.parquet`
-      from the same adapter run.
+      AASHTO persons -> DOTr-style occupants supplement for 2023-2025.
+      "Ghost-Driver" placeholder rows (Person Type=Driver + blank
+      Position) get `pos=0` so VTC aggregation routes them to the
+      'u' (unknown-VT) bucket. Ped/cyclist fatals remain
+      *under*-reported at person-level (those mostly land in 'u');
+      crash-level `pk` in `aashto_supplemented_crashes.parquet`
+      remains authoritative. Produced alongside
+      `aashto_supplemented_pedestrians.parquet` from the same run.

njdot/data/aashto_supplemented_pedestrians.parquet.dvc

@@ -7,7 +7,7 @@ meta:
   computation:
     cmd: python njdot/aashto/to_njdot_persons.py -y 2023,2024,2025
     git_deps:
-      /njdot/aashto/to_njdot_persons.py: fecfc0545896518b3047c8ffb0c1d4c054fbd596
+      /njdot/aashto/to_njdot_persons.py: 2f0d9eb2505089b276cd82eff4d3d2b77f34f615
       /njdot/aashto/to_njdot_schema.py: 0ad39efb3151ef4b58b4fd0d4c9bd8900074fce6
     deps:
       /njdot/data/2023/persons.parquet: 72625f5a685ceb7489789e886f4f3171

www/public/data/njdot/yccs.parquet

@@ -1,6 +1,6 @@
 outs:
-- md5: 103160954b4356002b8c445dd18402e2.dir
-  size: 810803
+- md5: c9627d129e88f4c9f48b48f692ccd2e8.dir
+  size: 1247206
   hash: md5
   path: ymccmc
 meta:
@@ -10,11 +10,11 @@ meta:
       /njdot/agg.py: ab9ab62177e15b118d5b80bef44e192ffb1bc7c2
     deps:
       /njdot/data/crashes.parquet: cfb770bf205ce85aedb90cc268607540
-      /njdot/data/aashto_supplemented_crashes.parquet: 2cd7616a654b357444d08aa05c41b311
+      /njdot/data/aashto_supplemented_crashes.parquet: 92ac7e301eab0f949748c6061317cf9f
     note: |
       `agg.py` produces multiple outputs in `www/public/data/njdot/`:
       smaller state/county aggregates (`ys`/`yms`/`yccs`/`ymccs`) live
-      in git; the larger muni-level dirs (`ymccmc/`, `ymccmcs/` —
+      in git; the larger muni-level dirs (`ymccmc/`, `ymccmcs/` -
       ~2.4MB across 42 files) are DVX-tracked here. `ymccmcs.dvc` is
       a co-output of the same `agg.py` invocation; dvx runs the cmd
       twice when both are stale (~10s, fine).