1"""Descriptive figures and tables for WP36 (Phase 4a).23Inputs (read-only): data/processed/listings.parquet, data/processed/price_events.parquet,4data/interim/sync_log.parquet, data/processed/sources_registry.csv,5data/processed/sample_construction.csv.67Outputs: figures/fig01_coverage, fig03_listing_age_and_cuts, fig04_price_paths (PDF + PNG,8copied to paper/figures/), paper/tables/tab01_sources.tex, tab02_sample_construction.tex,9tab03_descriptives.tex, results/descriptives/*.csv|json (every number quoted in captions is10in results/descriptives/summary.json).1112Conventions: day index t = UTC calendar days since 2026-08-12 (day 0, partial: crawl started1312:30 UTC); day 81 (2026-11-01) is partial (snapshot 06:49 UTC). "Exit" = observed14disappearance confirmed by >= 7 days of absence before the snapshot — not a sale.15"""16from __future__ import annotations1718import hashlib19import json20import math21from pathlib import Path2223import matplotlib.pyplot as plt24import numpy as np25import pandas as pd26import polars as pl2728from . import config as C29from . import plotstyle as PS3031OUT = C.RESULTS / "descriptives"32DAY0 = math.floor(C.WINDOW_START_TS / 86400)33T = int(math.floor(C.SNAPSHOT_TS / 86400) - DAY0) + 1 # number of day indices (82)34DAY_LABELS = {0: "2026-08-12", 81: "2026-11-01"}35WITHIN_DAYS = 30 # exit-within-30-days horizon36CLASS_SHORT = {"direct_dealer": "Dealer websites", "full_portal": "Full-inventory portals", "banner": "Dealer groups",37 "sampled_portal": "Sampled portals", "private_marketplace": "Private sellers (Kijiji)", "small_portal": "Small portals"}383940# ----------------------------------------------------------------------------- helpers41def day_of(ts: pl.Expr) -> pl.Expr:42 return ((ts / 86400).floor() - DAY0).cast(pl.Int32)434445def fmt_int(x) -> str:46 """Integer with thin-space thousands separators (LaTeX)."""47 return f"{int(round(float(x))):,}".replace(",", "\\,")484950def fmt1(x) -> str:51 return "--" if x is None or (isinstance(x, float) and np.isnan(x)) else f"{float(x):.1f}"525354def fmt_pct(x) -> str:55 return "--" if x is None or (isinstance(x, float) and np.isnan(x)) else f"{100 * float(x):.1f}"565758def write_tex(path: Path, lines: list[str]) -> None:59 path.parent.mkdir(parents=True, exist_ok=True)60 path.write_text("\n".join(lines) + "\n", encoding="utf-8")616263def load() -> tuple[pl.DataFrame, pl.DataFrame, pl.DataFrame, pd.DataFrame, pd.DataFrame]:64 L = pl.read_parquet(C.PROCESSED / "listings.parquet")65 E = pl.read_parquet(C.PROCESSED / "price_events.parquet")66 S = pl.read_parquet(C.INTERIM / "sync_log.parquet")67 reg = pd.read_csv(C.PROCESSED / "sources_registry.csv")68 sc = pd.read_csv(C.PROCESSED / "sample_construction.csv")69 L = L.with_columns(day_of(pl.col("entry_ts")).alias("e_day"),70 pl.min_horizontal(day_of(pl.col("end_ts")), pl.lit(T - 1)).alias("x_day"))71 return L, E, S, reg, sc727374# ----------------------------------------------------------------------------- fig0175def daily_coverage(L: pl.DataFrame, S: pl.DataFrame) -> tuple[pd.DataFrame, dict]:76 """Daily active stock / entries / confirmed exits in S2L, and crawl cycles per day."""77 U = L.filter(pl.col("s2l_dynamic_with_stock"))78 e = U["e_day"].to_numpy(); x = U["x_day"].to_numpy(); ex = U["exit_confirmed"].to_numpy()79 days = np.arange(T)80 entries = np.bincount(e, minlength=T)[:T]81 exits = np.bincount(x[ex], minlength=T)[:T]82 ends = np.bincount(x, minlength=T)[:T]83 # observed-active on day t: entered on or before t and not yet ended before t (ends are inclusive days)84 active = np.cumsum(entries) - np.concatenate([[0], np.cumsum(ends)[:-1]])85 # crawl cycles per day = sync rows / distinct sources logged that day (one row per source per cycle)86 sd = S.with_columns(day_of(pl.col("ts")).alias("day"))87 agg = (sd.group_by("day").agg(pl.len().alias("rows"), pl.col("source").n_unique().alias("n_sources"),88 (pl.col("ok") == 0).sum().alias("failed_rows"))89 .with_columns((pl.col("rows") / pl.col("n_sources")).alias("cycles")).sort("day").to_pandas())90 cyc = np.zeros(T); cyc[agg["day"].to_numpy()] = agg["cycles"].to_numpy()91 nsrc = np.zeros(T, dtype=int); nsrc[agg["day"].to_numpy()] = agg["n_sources"].to_numpy()92 fail = np.zeros(T); fail[agg["day"].to_numpy()] = (agg["failed_rows"] / agg["rows"]).to_numpy()93 # largest gap between consecutive crawls of the same source (hours)94 gaps = (sd.sort(["source", "ts"]).with_columns(pl.col("ts").diff().over("source").alias("dt"))95 .filter(pl.col("dt").is_not_null()).select((pl.col("dt") / 3600).alias("h")))96 # connector start-up stock: for portal sources first crawled after day 0, the first crawl day's entries are97 # initial stock, not an entry flow (not flagged left_censored because the cut-off is 2026-08-13 00:00 UTC)98 first_day = sd.group_by("source").agg(pl.col("day").min().alias("first")).filter(pl.col("first") > 0)99 st = (U.join(first_day, on="source", how="inner").filter(pl.col("e_day") == pl.col("first"))100 .group_by(["source", "e_day"]).agg(pl.len().alias("n")).filter(pl.col("n") >= 1000).sort("e_day"))101 startup_src = {str(int(r["e_day"])): r["source"] for r in st.iter_rows(named=True)}102 startup_n = {str(int(r["e_day"])): int(entries[int(r["e_day"])]) for r in st.iter_rows(named=True)}103 startup_listings = {r["source"]: int(r["n"]) for r in st.iter_rows(named=True)}104 cov = pd.DataFrame({"day": days, "active_s2l": active, "entries_s2l": entries, "exits_confirmed_s2l": exits,105 "crawl_cycles": cyc, "sources_logged": nsrc, "failed_share": fail})106 cov["date_utc"] = pd.to_datetime((DAY0 + days) * 86400, unit="s").strftime("%Y-%m-%d")107 full = cov[(cov.day >= 1) & (cov.day <= T - 2)]108 info = {"T_days": int(T), "day0_initial_stock_s2l": int(entries[0]),109 "connector_startup_entries": startup_n, "connector_startup_source": startup_src, "connector_startup_listings_s2l": startup_listings,110 "connector_startup_listings_s2_total": int(L.filter(pl.col("s2_dynamic")).join(first_day, on="source", how="inner")111 .filter((pl.col("e_day") == pl.col("first")) & pl.col("source").is_in(list(startup_listings))).height),112 "entries_per_full_day_mean_excl_startup": float(full[~full.day.isin([int(k) for k in startup_n])].entries_s2l.mean()), "entries_day0_left_censored": int(U["left_censored"].sum()),113 "active_s2l_min": int(full.active_s2l.min()), "active_s2l_max": int(full.active_s2l.max()),114 "active_s2l_mean_full_days": float(full.active_s2l.mean()), "active_s2l_at_last_full_day": int(active[T - 2]),115 "entries_per_full_day_mean": float(full.entries_s2l.mean()), "entries_per_full_day_median": float(full.entries_s2l.median()),116 "exits_per_full_day_mean_confirmable": float(cov[(cov.day >= 1) & (cov.day <= T - 1 - C.EXIT_CONFIRM_DAYS)].exits_confirmed_s2l.mean()),117 "crawl_cycles_per_full_day_mean": float(full.crawl_cycles.mean()), "crawl_cycles_per_full_day_min": float(full.crawl_cycles.min()),118 "crawl_cycles_per_full_day_max": float(full.crawl_cycles.max()), "crawl_cycles_total": int(round(cyc.sum())),119 "days_with_fewer_than_4_cycles": [int(d) for d in full.day[full.crawl_cycles < 4]],120 "max_gap_hours_same_source": float(gaps["h"].max()), "median_gap_hours_same_source": float(gaps["h"].median()),121 "sources_logged_min": int(full.sources_logged.min()), "sources_logged_max": int(full.sources_logged.max()),122 "failed_source_cycles_share": float((S["ok"] == 0).mean())}123 return cov, info124125126def region_stock(L: pl.DataFrame) -> pd.DataFrame:127 """Active canonical auto listings by region at the snapshot (all sources)."""128 A = L.filter((pl.col("active") == 1) & pl.col("is_canonical") & (pl.col("kind") == "auto"))129 r = (A.group_by(pl.col("region").fill_null("Region missing")).agg(pl.len().alias("active_canonical"))130 .sort("active_canonical", descending=True).to_pandas())131 r["share"] = r.active_canonical / r.active_canonical.sum()132 return r133134135def fig01(cov: pd.DataFrame, reg: pd.DataFrame, info: dict) -> None:136 fig = plt.figure(figsize=(6.3, 4.6))137 gs = fig.add_gridspec(3, 2, width_ratios=[1.35, 1.0])138 ax_a1 = fig.add_subplot(gs[0, 0]); ax_a2 = fig.add_subplot(gs[1, 0], sharex=ax_a1); ax_b = fig.add_subplot(gs[2, 0], sharex=ax_a1)139 ax_c = fig.add_subplot(gs[:, 1])140 d = cov.day.to_numpy()141 c_act, c_ent, c_ex = PS.CAT[0], PS.CAT[2], PS.CAT[1]142 # (a1) active stock143 ax_a1.plot(d, cov.active_s2l / 1000, color=c_act, lw=1.4)144 ax_a1.set_ylabel("Active listings\n(thousands)")145 ax_a1.set_ylim(0, None)146 ax_a1.set_title("(a) Dynamic-reliable universe (S2L): stock and flows", loc="left")147 # (a2) entries / exits (day 0 = initial left-censored stock, not an entry flow)148 m = d >= 1149 ent = cov.entries_s2l.to_numpy().astype(float)150 spikes = {int(k): int(v) for k, v in info["connector_startup_entries"].items()}151 ycap = float(np.quantile(ent[m & ~np.isin(d, list(spikes))], 1.0)) * 1.5152 ax_a2.plot(d[m], np.minimum(ent[m], ycap), color=c_ent, lw=1.1, label="Entries")153 mx = m & (d <= T - 1 - C.EXIT_CONFIRM_DAYS)154 ax_a2.plot(d[mx], cov.exits_confirmed_s2l[mx], color=c_ex, lw=1.1, label="Confirmed exits")155 ax_a2.plot(list(spikes), [ycap] * len(spikes), ls="none", marker="^", ms=4, color=c_ent, clip_on=False)156 note = "\u25b2 truncated: connector\n start-up stock\n" + "\n".join(157 f" day {k}: {v:,} ({info['connector_startup_source'][str(k)]})".replace(",", "\u2009") for k, v in spikes.items())158 ax_a2.text(0.14, 0.97, note, transform=ax_a2.transAxes, fontsize=6, color=PS.TEXT2, ha="left", va="top", linespacing=1.1)159 ax_a2.axvspan(T - 1 - C.EXIT_CONFIRM_DAYS, T - 1, color=PS.GRID, alpha=0.6, lw=0, zorder=0)160 ax_a2.set_ylabel("Listings per day")161 ax_a2.set_ylim(0, ycap)162 ax_a2.legend(loc="upper right", bbox_to_anchor=(1.0, 1.02), ncol=1, handlelength=1.4, labelspacing=0.3)163 # (b) crawl cycles per day164 ax_b.bar(d, cov.crawl_cycles, width=0.8, color=PS.CAT[0], alpha=0.85)165 ax_b.set_ylabel("Crawl cycles\nper day")166 ax_b.set_xlabel("Day index (0 = 2026-08-12 UTC; days 0 and 81 are partial)")167 ax_b.set_title("(b) Crawl cycles per day, all sources", loc="left")168 ax_b.set_ylim(0, max(12, cov.crawl_cycles.max() * 1.1))169 ax_b.set_xlim(-1, T)170 for ax in (ax_a1, ax_a2):171 plt.setp(ax.get_xticklabels(), visible=False)172 # (c) regions173 rr = reg.sort_values("active_canonical", ascending=True)174 cols = [PS.color_for(PS.REGIONS, k) for k in rr.region]175 ax_c.barh(rr.region, rr.active_canonical / 1000, color=cols, height=0.72)176 ax_c.set_xlabel("Active canonical listings (thousands)")177 ax_c.set_title("(c) Stock by region, 2026-11-01", loc="left")178 ax_c.tick_params(axis="y", labelsize=7)179 ax_c.grid(axis="y", visible=False)180 for y, v in enumerate(rr.active_canonical):181 ax_c.text(v / 1000 + 0.08, y, f"{v:,}".replace(",", " "), va="center", fontsize=6, color=PS.TEXT2)182 ax_c.set_xlim(0, rr.active_canonical.max() / 1000 * 1.22)183 PS.save(fig, "fig01_coverage")184185186# ----------------------------------------------------------------------------- fig03187def duration_and_cuts(L: pl.DataFrame) -> tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame, dict]:188 S2 = L.filter(pl.col("s2_dynamic"))189 done = S2.filter(pl.col("exit_confirmed"))["duration_days"].to_numpy()190 done = np.sort(done)191 ecdf = pd.DataFrame({"duration_days": done, "ecdf": np.arange(1, done.size + 1) / done.size})192 fc = -100 * S2.filter(pl.col("n_cuts") >= 1)["first_cut_pct"].to_numpy()193 cuts = pd.DataFrame({"first_cut_pct": np.sort(fc)})194 bc = (S2.group_by("source_class").agg(pl.len().alias("n"), (pl.col("n_cuts") >= 1).sum().alias("n_with_cut"),195 pl.col("exit_confirmed").sum().alias("n_exit_confirmed"))196 .with_columns((pl.col("n_with_cut") / pl.col("n")).alias("share_with_cut")).to_pandas())197 bc["order"] = bc.source_class.map(PS.CLASS_ORDER.index); bc = bc.sort_values("order").drop(columns="order")198 S2L = L.filter(pl.col("s2l_dynamic_with_stock"))199 q = lambda a, p: float(np.quantile(a, p))200 info = {"S2_n": S2.height, "S2_exit_confirmed": int(S2["exit_confirmed"].sum()),201 "S2_right_censored_share": float(1 - S2["exit_confirmed"].mean()),202 "S2L_n": S2L.height, "S2L_left_censored_n": int(S2L["left_censored"].sum()), "S2L_left_censored_share": float(S2L["left_censored"].mean()),203 "S2_completed_duration_median_days": q(done, .5), "S2_completed_duration_p25_days": q(done, .25), "S2_completed_duration_p75_days": q(done, .75),204 "S2_completed_duration_mean_days": float(done.mean()),205 "S2_completed_share_within_7d": float((done <= 7).mean()), "S2_completed_share_within_30d": float((done <= 30).mean()),206 "S2_share_with_cut": float((S2["n_cuts"] >= 1).mean()), "S2_n_with_cut": int((S2["n_cuts"] >= 1).sum()),207 "S2_share_with_2plus_cuts": float((S2["n_cuts"] >= 2).mean()),208 "S2_first_cut_pct_median": q(fc, .5), "S2_first_cut_pct_mean": float(fc.mean()), "S2_first_cut_pct_p25": q(fc, .25), "S2_first_cut_pct_p75": q(fc, .75),209 "S2_first_cut_share_le_5pct": float((fc <= 5).mean()), "S2_first_cut_share_gt_10pct": float((fc > 10).mean()),210 "S2_share_with_cut_by_class": {r.source_class: float(r.share_with_cut) for r in bc.itertuples()}}211 return ecdf, cuts, bc, info212213214def fig03(ecdf: pd.DataFrame, cuts: pd.DataFrame, bc: pd.DataFrame, info: dict) -> None:215 fig, axes = plt.subplots(1, 3, figsize=(6.3, 2.3), gridspec_kw={"width_ratios": [1.15, 1.15, 1.0]})216 a, b, c = axes217 # (a) ECDF of completed durations218 x = ecdf.duration_days.to_numpy(); y = ecdf.ecdf.to_numpy()219 a.step(x, y, where="post", color=PS.CAT[0], lw=1.2)220 med = info["S2_completed_duration_median_days"]221 a.axvline(med, color=PS.NEUTRAL, lw=0.7, ls=(0, (3, 2)))222 a.text(med + 1.5, 0.08, f"median {med:.1f} d", fontsize=7, color=PS.TEXT2)223 a.set_xlim(0, 75); a.set_ylim(0, 1)224 a.set_xlabel("Completed duration (days)"); a.set_ylabel("Cumulative share of spells")225 a.set_title("(a) Completed spells, S2", loc="left")226 # (b) first markdown size227 v = cuts.first_cut_pct.to_numpy()228 bins = np.arange(1, 25.5, 0.5)229 b.hist(np.clip(v, 1, 25.4), bins=bins, color=PS.CAT[1], alpha=0.9, lw=0, density=False, weights=np.full(v.size, 100 / v.size))230 medc = info["S2_first_cut_pct_median"]231 b.axvline(medc, color=PS.NEUTRAL, lw=0.7, ls=(0, (3, 2)))232 b.text(medc + 0.6, b.get_ylim()[1] * 0.97, f"median {medc:.1f} %", fontsize=7, color=PS.TEXT2, va="top")233 b.set_xlim(1, 25.5)234 b.set_xlabel("First markdown (% of entry price; \u2265 25 pooled)"); b.set_ylabel("Share of first markdowns (%)")235 b.set_title("(b) First markdown size, S2", loc="left")236 # (c) share with >=1 cut by source class237 labels = [CLASS_SHORT[k] for k in bc.source_class]238 cols = [PS.color_for(PS.CLASS_ORDER, k) for k in bc.source_class]239 c.bar(range(len(bc)), 100 * bc.share_with_cut, color=cols, width=0.65)240 c.set_xticks(range(len(bc)), [l.replace(" ", "\n", 1) for l in labels], fontsize=7)241 for i, (s, n) in enumerate(zip(bc.share_with_cut, bc.n)):242 c.text(i, 100 * s + 0.6, f"{100 * s:.1f} %", ha="center", fontsize=7, color=PS.TEXT2)243 c.set_ylabel("Listings with ≥1 markdown (%)")244 c.set_ylim(0, max(100 * bc.share_with_cut) * 1.25)245 c.grid(axis="x", visible=False)246 c.set_title("(c) Markdown incidence", loc="left"); c.set_xlabel("Source class (S2)")247 PS.save(fig, "fig03_listing_age_and_cuts")248249250# ----------------------------------------------------------------------------- fig04251def select_paths(L: pl.DataFrame) -> tuple[pl.DataFrame, dict]:252 """Deterministic selection: S2 listings with >= 3 sustained cuts and 30 <= duration <= 70 days,253 2 per body type among {SUV, Sedan, Pickup} (the three largest body types), chosen by the254 lowest SHA-256 of the anonymised listing_id (hex, ascending)."""255 bodies = ["SUV", "Sedan", "Pickup"]256 cand = (L.filter(pl.col("s2_dynamic") & (pl.col("n_cuts") >= 3) & (pl.col("duration_days") >= 30) & (pl.col("duration_days") <= 70)257 & pl.col("body_filled").is_in(bodies))258 .with_columns(pl.col("listing_id").map_elements(lambda s: hashlib.sha256(s.encode()).hexdigest(), return_dtype=pl.String).alias("h")))259 picks = (cand.sort(["body_filled", "h"]).group_by("body_filled", maintain_order=True).head(2)260 .with_columns(pl.col("body_filled").replace_strict({b: i for i, b in enumerate(bodies)}, return_dtype=pl.Int64).alias("bo"))261 .sort(["bo", "h"]).drop("bo"))262 picks = picks.with_columns(pl.Series("label", [f"Listing {chr(65 + i)}" for i in range(picks.height)]))263 rule = {"rule": "S2 listings with n_cuts >= 3 and 30 <= duration_days <= 70, body in {SUV, Sedan, Pickup}; "264 "within each body type the two listings with the lowest SHA-256(listing_id) hex digest; labelled A-F in body order.",265 "n_candidates": cand.height, "n_candidates_by_body": {b: int(cand.filter(pl.col("body_filled") == b).height) for b in bodies},266 "selected": [{"label": r["label"], "listing_id": r["listing_id"], "sha256_prefix": r["h"][:12], "make": r["make"], "body": r["body_filled"],267 "year": int(r["year"]), "n_cuts": int(r["n_cuts"]), "duration_days": round(float(r["duration_days"]), 1),268 "exit_confirmed": bool(r["exit_confirmed"]), "entry_price": float(r["entry_price"]), "last_price": float(r["last_price"]),269 "total_change_pct": round(100 * float(r["total_change_pct"]), 1)}270 for r in picks.iter_rows(named=True)]}271 return picks, rule272273274def fig04(picks: pl.DataFrame, E: pl.DataFrame) -> pd.DataFrame:275 ev = (E.filter(pl.col("listing_id").is_in(picks["listing_id"].to_list())).join(276 picks.select(["listing_id", "label", "entry_ts", "end_ts", "exit_confirmed", "make", "body_filled", "year"]), on="listing_id")277 .with_columns(((pl.col("ts") - pl.col("entry_ts")) / 86400).alias("days"))278 .sort(["label", "ts"]))279 fig, axes = plt.subplots(2, 3, figsize=(6.3, 3.7), sharex=True)280 rows = []281 for ax, (lab, g) in zip(axes.flat, ev.group_by("label", maintain_order=True)):282 lab = lab[0] if isinstance(lab, tuple) else lab283 d = g["days"].to_numpy(); p = g["price"].to_numpy() / 1000284 d0 = np.concatenate([[0.0], d[1:]]) if d.size else d # entry level starts at day 0285 dur = float((g["end_ts"][0] - g["entry_ts"][0]) / 86400)286 body = g["body_filled"][0]; col = PS.color_for(PS.BODY_ORDER, body)287 ax.step(np.concatenate([d0, [dur]]), np.concatenate([p, [p[-1]]]), where="post", color=col, lw=1.3)288 ax.plot(d[1:], p[1:], ls="none", marker="o", ms=2.6, color=col, mfc=PS.SURFACE, mew=0.9)289 if bool(g["exit_confirmed"][0]):290 ax.plot([dur], [p[-1]], ls="none", marker="x", ms=4.5, color=PS.TEXT2, mew=0.9)291 ax.set_title(f"{lab} — {g['make'][0]} {body} {int(g['year'][0])}", loc="left", fontsize=8)292 ax.set_xlim(0, 72)293 lo, hi = p.min(), p.max(); pad = max(0.6, 0.12 * (hi - lo))294 ax.set_ylim(lo - pad, hi + pad)295 for r in g.iter_rows(named=True):296 rows.append({"label": lab, "days_since_entry": round(r["days"], 3), "price_cad": r["price"], "kind": r["kind"]})297 rows.append({"label": lab, "days_since_entry": round(dur, 3), "price_cad": g["price"][-1], "kind": "exit" if bool(g["exit_confirmed"][0]) else "censored"})298 for ax in axes[1]:299 ax.set_xlabel("Days since entry")300 fig.supylabel("Asking price (CAD, thousands)", fontsize=9)301 from matplotlib.lines import Line2D302 h = [Line2D([], [], color=PS.TEXT2, lw=1.2, label="Sustained price level"),303 Line2D([], [], ls="none", marker="o", ms=2.6, color=PS.TEXT2, mfc=PS.SURFACE, label="Price change"),304 Line2D([], [], ls="none", marker="x", ms=4.5, color=PS.TEXT2, label="Confirmed exit")]305 fig.legend(handles=h, loc="outside lower center", ncol=3, fontsize=7, handlelength=1.6, columnspacing=1.6)306 PS.save(fig, "fig04_price_paths")307 return pd.DataFrame(rows)308309310# ----------------------------------------------------------------------------- tables311def tab01(L: pl.DataFrame, reg: pd.DataFrame) -> pd.DataFrame:312 by = (L.group_by("source_class").agg(pl.col("source").n_unique().alias("n_sources_with_listings"), pl.col("dealer_id").n_unique().alias("n_dealers"),313 pl.len().alias("listings"), pl.col("active").sum().alias("active_at_snapshot"),314 pl.col("s1_hedonic").any().alias("in_s1"), pl.col("s2_dynamic").any().alias("in_s2"),315 pl.col("source").filter(pl.col("s2_dynamic")).n_unique().alias("n_sources_in_s2"),316 pl.col("s1_hedonic").sum().alias("n_s1"), pl.col("s2_dynamic").sum().alias("n_s2"))317 .to_pandas().set_index("source_class"))318 rg = reg.groupby("class").agg(n_sources=("source", "nunique"), found_per_cycle=("mean_found", "sum"), mean_found_per_source=("mean_found", "mean"))319 t = by.join(rg).reindex(PS.CLASS_ORDER)320 t["label"] = [PS.CLASS_LABEL[k] for k in t.index]321 t.index.name = "source_class"322 rows = [r"\begin{tabular}{lrrrrrcc}", r"\toprule",323 r"Source class & Sources & Dealers & Listings & Active & Found/cycle & S1 & S2 \\", r"\midrule"]324 for k, r in t.iterrows():325 s1 = "Yes" if r.in_s1 else "No"326 s2 = ("Yes" if r.n_sources_in_s2 == r.n_sources else f"Yes ({int(r.n_sources_in_s2)} of {int(r.n_sources)})") if r.in_s2 else "No"327 rows.append(f"{CLASS_SHORT[k]} & {fmt_int(r.n_sources)} & {fmt_int(r.n_dealers)} & {fmt_int(r.listings)} & {fmt_int(r.active_at_snapshot)} & {fmt_int(r.found_per_cycle)} & {s1} & {s2} \\\\")328 rows += [r"\midrule",329 f"All sources & {fmt_int(t.n_sources.sum())} & {fmt_int(L['dealer_id'].n_unique())} & {fmt_int(t.listings.sum())} & {fmt_int(t.active_at_snapshot.sum())} & {fmt_int(t.found_per_cycle.sum())} & & \\\\",330 r"\bottomrule", r"\end{tabular}"]331 write_tex(C.TABLES / "tab01_sources.tex", rows)332 return t.reset_index()333334335def tab02(sc: pd.DataFrame) -> pd.DataFrame:336 sc = sc.copy()337 sc["dropped"] = (-sc.listings_remaining.diff()).fillna(0).astype(int)338 rows = [r"\begin{tabular}{lrr}", r"\toprule", r"Step & Remaining & Dropped \\", r"\midrule"]339 marks = {"− make or model missing": None, "− uid identity collision": "S1 (all sources)",340 "− sampled portal (autotrader/cargurus) or private marketplace (kijiji)": "S1 hedonic",341 "− oscillating price path": "S2L", "− left-censored day-1 stock": "S2 dynamic"}342 for _, r in sc.iterrows():343 step = r.step.replace("−", "$-$").replace("[", "[").replace("$", r"\$").replace(r"\$-\$", "$-$")344 tag = marks.get(r.step)345 step = step if tag is None else f"{step} \\hfill $\\Rightarrow$ {tag}"346 rows.append(f"{step} & {fmt_int(r.listings_remaining)} & {fmt_int(r.dropped) if r.dropped else '--'} \\\\")347 rows += [r"\bottomrule", r"\end{tabular}"]348 write_tex(C.TABLES / "tab02_sample_construction.tex", rows)349 return sc350351352def tab03(L: pl.DataFrame) -> tuple[dict, pd.DataFrame]:353 S1 = L.filter(pl.col("s1_hedonic")); S2 = L.filter(pl.col("s2_dynamic"))354355 def block(D: pl.DataFrame, dyn_completed: pl.DataFrame) -> dict:356 n = D.height357 d = {"N": n, "entry_price_mean": D["entry_price"].mean(), "entry_price_median": D["entry_price"].median(),358 "mileage_mean_km": D["mileage_km"].mean(), "mileage_median_km": D["mileage_km"].median(),359 "age_mean_years": D["age_years_entry"].mean(), "age_median_years": D["age_years_entry"].median(),360 "dealers": D["dealer_id"].n_unique(), "active_share": D["active"].mean(),361 "powertrain_share": {k: D.filter(pl.col("powertrain") == k).height / n for k in PS.POWERTRAIN_ORDER},362 "body_share": {k: D.filter(pl.col("body_filled") == k).height / n for k in PS.BODY_ORDER},363 "body_missing_share": D["body_filled"].null_count() / n,364 "duration_median_completed": dyn_completed["duration_days"].median(), "n_completed": dyn_completed.height,365 "share_with_cut": (D["n_cuts"] >= 1).mean(),366 "first_cut_pct_median": -100 * D.filter(pl.col("n_cuts") >= 1)["first_cut_pct"].median()}367 return d368 s1 = block(S1, S1.filter(pl.col("exit_confirmed") & ~pl.col("left_censored")))369 s2 = block(S2, S2.filter(pl.col("exit_confirmed")))370 top_bodies = sorted(PS.BODY_ORDER, key=lambda k: -s1["body_share"][k])[:5]371 # by region (S1), with exit-within-30d among S2 entrants observed long enough to confirm an exit at day 30372 horizon = (WITHIN_DAYS + C.EXIT_CONFIRM_DAYS) * 86400373 S2e = S2.filter(pl.col("entry_ts") <= C.SNAPSHOT_TS - horizon).with_columns(374 (pl.col("exit_confirmed") & (pl.col("duration_days") <= WITHIN_DAYS)).alias("exit30"))375 rs1 = S1.group_by(pl.col("region").fill_null("Region missing")).agg(pl.len().alias("N"), pl.col("entry_price").median().alias("median_price"),376 pl.col("mileage_km").median().alias("median_km"), pl.col("age_years_entry").median().alias("median_age"))377 rs2 = S2e.group_by(pl.col("region").fill_null("Region missing")).agg(pl.len().alias("n_s2_entrants_37d"), pl.col("exit30").mean().alias("share_exit_30d"))378 R = rs1.join(rs2, on="region", how="left").to_pandas().set_index("region")379 order = [r for r in PS.REGIONS if r in R.index] + [r for r in R.index if r not in PS.REGIONS]380 R = R.reindex(order).reset_index()381 tot = {"N": S1.height, "median_price": S1["entry_price"].median(), "median_km": S1["mileage_km"].median(), "median_age": S1["age_years_entry"].median(),382 "n_s2_entrants_37d": S2e.height, "share_exit_30d": S2e["exit30"].mean()}383384 rows = [r"\begin{tabular}{lrr}", r"\toprule", r" & S1 hedonic & S2 dynamic \\", r"\midrule",385 f"Listings (N) & {fmt_int(s1['N'])} & {fmt_int(s2['N'])} \\\\",386 f"Dealers & {fmt_int(s1['dealers'])} & {fmt_int(s2['dealers'])} \\\\",387 f"Entry asking price, mean (\\$) & {fmt_int(s1['entry_price_mean'])} & {fmt_int(s2['entry_price_mean'])} \\\\",388 f"Entry asking price, median (\\$) & {fmt_int(s1['entry_price_median'])} & {fmt_int(s2['entry_price_median'])} \\\\",389 f"Mileage, mean (km) & {fmt_int(s1['mileage_mean_km'])} & {fmt_int(s2['mileage_mean_km'])} \\\\",390 f"Vehicle age at entry, mean (years) & {fmt1(s1['age_mean_years'])} & {fmt1(s2['age_mean_years'])} \\\\",391 r"\addlinespace", r"\multicolumn{3}{l}{\textit{Powertrain (\% of listings)}} \\"]392 for k in PS.POWERTRAIN_ORDER:393 rows.append(f"\\quad {k if k != 'unknown' else 'Unknown'} & {fmt_pct(s1['powertrain_share'][k])} & {fmt_pct(s2['powertrain_share'][k])} \\\\")394 rows += [r"\addlinespace", r"\multicolumn{3}{l}{\textit{Body type, top 5 (\% of listings)}} \\"]395 for k in top_bodies:396 rows.append(f"\\quad {k} & {fmt_pct(s1['body_share'][k])} & {fmt_pct(s2['body_share'][k])} \\\\")397 rows += [r"\addlinespace", r"\multicolumn{3}{l}{\textit{Listing dynamics}} \\",398 f"\\quad Completed spells (confirmed exits) & {fmt_int(s1['n_completed'])} & {fmt_int(s2['n_completed'])} \\\\",399 f"\\quad Duration of completed spells, median (days) & {fmt1(s1['duration_median_completed'])} & {fmt1(s2['duration_median_completed'])} \\\\",400 f"\\quad Listings with $\\geq 1$ sustained markdown (\\%) & {fmt_pct(s1['share_with_cut'])} & {fmt_pct(s2['share_with_cut'])} \\\\",401 f"\\quad First markdown, median (\\% of entry price) & {fmt1(s1['first_cut_pct_median'])} & {fmt1(s2['first_cut_pct_median'])} \\\\",402 r"\bottomrule", r"\end{tabular}", "",403 r"\begin{tabular}{lrrrrr}", r"\toprule",404 r"Region & N (S1) & Median price (\$) & Median km & Median age (y) & Exit $\leq$ 30 d (\%) \\", r"\midrule"]405 for _, r in R.iterrows():406 rows.append(f"{r.region} & {fmt_int(r.N)} & {fmt_int(r.median_price)} & {fmt_int(r.median_km)} & {fmt1(r.median_age)} & {fmt_pct(r.share_exit_30d)} \\\\")407 rows += [r"\midrule",408 f"All regions & {fmt_int(tot['N'])} & {fmt_int(tot['median_price'])} & {fmt_int(tot['median_km'])} & {fmt1(tot['median_age'])} & {fmt_pct(tot['share_exit_30d'])} \\\\",409 r"\bottomrule", r"\end{tabular}"]410 write_tex(C.TABLES / "tab03_descriptives.tex", rows)411 info = {"S1": s1, "S2": s2, "top_bodies_S1": top_bodies, "region_total": tot,412 "exit30_definition": f"S2 listings with entry_ts <= snapshot - {WITHIN_DAYS + C.EXIT_CONFIRM_DAYS} days (so that an exit at day {WITHIN_DAYS} is confirmable); "413 f"exit within {WITHIN_DAYS} d = exit_confirmed and duration_days <= {WITHIN_DAYS}",414 "S1_duration_note": "S1 completed spells exclude left-censored day-1 stock (duration unknown) but include all S1 sources"}415 return info, R416417418# ----------------------------------------------------------------------------- main419def run() -> dict:420 PS.setup()421 OUT.mkdir(parents=True, exist_ok=True)422 L, E, S, reg, sc = load()423 cov, cinfo = daily_coverage(L, S)424 rstock = region_stock(L)425 fig01(cov, rstock, cinfo)426 cov.to_csv(OUT / "fig01_daily_coverage.csv", index=False)427 rstock.to_csv(OUT / "fig01_region_stock.csv", index=False)428429 ecdf, cuts, bc, dinfo = duration_and_cuts(L)430 fig03(ecdf, cuts, bc, dinfo)431 ecdf.iloc[::25].to_csv(OUT / "fig03_duration_ecdf.csv", index=False)432 bc.to_csv(OUT / "fig03_cut_share_by_class.csv", index=False)433434 picks, rule = select_paths(L)435 paths = fig04(picks, E)436 (OUT / "fig04_selection.json").write_text(json.dumps(rule, indent=2, ensure_ascii=False), encoding="utf-8")437 paths.to_csv(OUT / "fig04_price_paths.csv", index=False)438439 t1 = tab01(L, reg); t1.to_csv(OUT / "tab01_sources.csv", index=False)440 t2 = tab02(sc); t2.to_csv(OUT / "tab02_sample_construction.csv", index=False)441 t3info, R = tab03(L); R.to_csv(OUT / "tab03_by_region.csv", index=False)442443 A = L.filter((pl.col("active") == 1) & pl.col("is_canonical") & (pl.col("kind") == "auto"))444 summary = {445 "snapshot_date": C.SNAPSHOT_DATE, "window_start_utc": "2026-08-12 12:30 UTC", "window_end_utc": "2026-11-01 06:49 UTC",446 "window_days": round((C.SNAPSHOT_TS - C.WINDOW_START_TS) / 86400, 1),447 "S0_n": L.height, "S1_n": t3info["S1"]["N"], "S2_n": t3info["S2"]["N"], "S2L_n": dinfo["S2L_n"],448 "active_canonical_auto_at_snapshot": A.height, "active_canonical_region_missing": int(A["region"].null_count()),449 "active_canonical_top_regions": rstock.head(5).to_dict(orient="records"),450 "sources_total": int(reg["source"].nunique()), "dealers_total": int(L["dealer_id"].n_unique()),451 "coverage": cinfo, "dynamics": dinfo, "fig04": rule, "tab01": t1.to_dict(orient="records"), "tab03": t3info,452 "tab03_by_region": R.to_dict(orient="records"),453 "notes": ["exit = observed disappearance from the source confirmed by >= 7 days of absence before the snapshot; it is NOT a sale",454 "day index 0 = 2026-08-12 UTC (partial day, crawl started 12:30 UTC); day 81 = 2026-11-01 (partial, snapshot 06:49 UTC)",455 "crawl cycles per day = sync_log rows / distinct sources logged that day (one row per source per cycle)",456 "S2 contains no left-censored listing by construction; the left-censored share refers to S2L"]}457 (OUT / "summary.json").write_text(json.dumps(summary, indent=2, ensure_ascii=False, default=float), encoding="utf-8")458 return summary459