1"""Replication entry point for F6-conf-02. Run from workspace root."""2import json, os3import numpy as np, pandas as pd4import matplotlib5matplotlib.use("Agg")6import matplotlib.pyplot as plt7from scipy import stats89SEED = None # no randomness used10os.makedirs("out", exist_ok=True)11d = pd.read_csv("data/gc5y_daily.csv")12n_rows = len(d)13d["date_dt"] = pd.to_datetime(d["date"])14d["chg_bp"] = 100 * d["yield_5y"].diff()1516# ---- audit ----17issues = []18gaps = d["date_dt"].diff().dt.days19n_missing = int(d[["date", "yield_5y", "post_cut"]].isna().sum().sum())20n_dup = int(d["date"].duplicated().sum())21mono = bool(d["date_dt"].is_monotonic_increasing)22mismatch = int(((d["date"] >= "2020-03-04") != (d["post_cut"] == 1)).sum())23n_treated = int((d["post_cut"] == 1).sum())24audit = dict(n_missing=n_missing, n_dup=n_dup, monotonic=mono, flag_mismatch=mismatch,25 max_gap_days=int(gaps.max()), n_gaps_gt4=int((gaps > 4).sum()))26print("audit:", audit)27if n_treated == 0:28 issues.append({"issue": f"Sample covers {d['date'].iloc[0]} to {d['date'].iloc[-1]} only (the Valet request used end_date=2019-12-31); "29 "no row has post_cut = 1, so there are zero post-cut trading days and the event (2020-03-04) lies outside the file.",30 "handling": "Effect estimand, its CI and the two-window model RMSE are not estimable; reported as null and conclusion set to premise_invalid. "31 "No proxy or substitute window was used."})32 issues.append({"issue": "The 120-trading-day pre-window is defined relative to the first flagged day, which does not exist in the file.",33 "handling": "Pre-window not constructed; not substituted by the last 120 rows of 2019."})34if n_missing: issues.append({"issue": f"{n_missing} missing values", "handling": "rows with missing yield would be dropped"})35if n_dup: issues.append({"issue": f"{n_dup} duplicate dates", "handling": "flagged"})36if mismatch: issues.append({"issue": f"{mismatch} rows where post_cut disagrees with date>=2020-03-04", "handling": "post_cut used as provided"})37issues.append({"issue": f"Calendar gaps between consecutive rows up to {int(gaps.max())} days ({int((gaps>4).sum())} gaps >4 days), consistent with weekends and Canadian holidays; no missing-value rows.",38 "handling": "One-day changes computed between consecutive rows as specified by TASK.md, regardless of calendar gap."})39issues.append({"issue": "First row has no previous observation, so no one-day change.", "handling": "Excluded from change-based statistics (2497 changes)."})4041baseline = float(d["chg_bp"].mean())4243effect = lo = hi = rmse = None; n_an = None44if n_treated > 0:45 first = d.index[d["post_cut"] == 1][0]46 pre = d.loc[first - 120:first - 1, "chg_bp"].dropna()47 post = d.loc[d["post_cut"] == 1, "chg_bp"].dropna()48 effect = float(post.mean() - pre.mean())49 v1, v2, n1, n2 = post.var(ddof=1), pre.var(ddof=1), len(post), len(pre)50 se = np.sqrt(v1 / n1 + v2 / n2)51 df = (v1 / n1 + v2 / n2) ** 2 / ((v1 / n1) ** 2 / (n1 - 1) + (v2 / n2) ** 2 / (n2 - 1))52 t = stats.t.ppf(0.975, df); lo, hi = effect - t * se, effect + t * se53 resid = np.concatenate([post - post.mean(), pre - pre.mean()])54 rmse = float(np.sqrt(np.mean(resid ** 2))); n_an = int(n1 + n2)5556if effect is None:57 concl, feas = "premise_invalid", "infeasible"58 reason = ("data/gc5y_daily.csv ends on 2019-12-31 and contains zero rows with post_cut = 1; "59 "the March-2020 post-cut window and the 120-day window preceding it are absent, so the estimand cannot be computed.")60else:61 concl = "positive_effect" if lo > 0 else "negative_effect" if hi < 0 else "no_detectable_effect"62 feas, reason = "complete", None6364res = {"task_id": "F6-conf-02", "sample_start": d["date"].iloc[0], "sample_end": d["date"].iloc[-1],65 "n_treated": n_treated, "baseline_metric": round(baseline, 6), "model_metric": rmse,66 "effect_estimate": effect, "effect_ci_low": lo, "effect_ci_high": hi, "conclusion": concl,67 "feasibility": feas, "feasibility_reason": reason, "data_issues": issues,68 "n_input_rows": {"gc5y_daily.csv": n_rows}, "n_analysed": n_an, "random_seed": SEED}69json.dump(res, open("out/results.json", "w"), indent=2)7071# table72d["year"] = d["date_dt"].dt.year73t1 = d.groupby("year").agg(n_days=("yield_5y", "size"), mean_yield=("yield_5y", "mean"),74 mean_change_bp=("chg_bp", "mean")).reset_index().round(4)75t1.to_csv("out/table1.csv", index=False)7677# figure78fig, ax = plt.subplots(figsize=(9, 4))79ax.plot(d["date_dt"], d["yield_5y"], lw=0.8)80ax.axvline(pd.Timestamp("2020-03-04"), color="r", ls="--", label="2020-03-04 cut (outside sample)")81ax.set_xlim(d["date_dt"].iloc[0], pd.Timestamp("2020-06-30"))82ax.set_ylabel("5-year GoC yield (%)"); ax.set_title("GoC 5-year benchmark yield, BD.CDN.5YR.DQ.YLD (file sample)")83ax.legend(); fig.tight_layout(); fig.savefig("out/figure1.png", dpi=120)84print(json.dumps({k: v for k, v in res.items() if k != "data_issues"}, indent=1))85