1"""Replication entry point for F6-conf-02. Run from workspace root."""2import json, os3import numpy as np, pandas as pd4import matplotlib; matplotlib.use("Agg")5import matplotlib.pyplot as plt6from scipy import stats78np.random.seed(0)9os.makedirs("out", exist_ok=True)10f = "data/gc5y_daily.csv"11d = pd.read_csv(f)12n_rows = len(d)13d["date"] = pd.to_datetime(d["date"])14d = d.sort_values("date").reset_index(drop=True)15d["chg_bp"] = 100 * d["yield_5y"].diff()1617issues = []18gaps = d["date"].diff().dt.days19issues.append({"issue": f"Sample covers {d.date.min().date()} to {d.date.max().date()} only, although the treatment date is 2020-03-04; the README's retrieval URL used end_date=2019-12-31.",20 "handling": "No post-cut observations exist; estimand declared not estimable (premise_invalid)."})21issues.append({"issue": f"post_cut is 0 for all {n_rows} rows (n_treated = {int((d.post_cut==1).sum())}).",22 "handling": "Treatment indicator used as provided, not redefined; effect, CI and model RMSE set to null."})23issues.append({"issue": f"Missing values: {int(d.isna().drop(columns='chg_bp').sum().sum())}; duplicate dates: {int(d.date.duplicated().sum())}; calendar gaps up to {int(gaps.max())} days (weekends/holidays, {int((gaps>4).sum())} gaps >4 days).",24 "handling": "No missing/duplicate rows; changes computed between consecutive rows as specified; first row has no change and is dropped from change means."})2526treated = d[d.post_cut == 1]27n_treated = int(len(treated))28baseline = float(d["chg_bp"].mean())2930effect = lo = hi = rmse = None; n_an = None31if n_treated > 0:32 first = treated.index[0]33 pre = d.loc[max(1, first-120):first-1, "chg_bp"]34 post = treated["chg_bp"].dropna()35 effect = float(post.mean() - pre.mean())36 v1, v0 = post.var(ddof=1)/len(post), pre.var(ddof=1)/len(pre)37 se = np.sqrt(v1+v0)38 df = (v1+v0)**2/(v1**2/(len(post)-1)+v0**2/(len(pre)-1))39 t = stats.t.ppf(0.975, df); lo, hi = effect-t*se, effect+t*se40 resid = np.concatenate([post-post.mean(), pre-pre.mean()])41 rmse = float(np.sqrt(np.mean(resid**2))); n_an = int(len(post)+len(pre))4243if effect is None:44 concl, feas = "premise_invalid", "infeasible"45 reason = ("data/gc5y_daily.csv ends on 2019-12-31 and contains no rows with post_cut = 1; "46 "the post-cut mean change (and hence the effect, its CI and the two-window model RMSE) cannot be computed. "47 "No proxy or extrapolated estimate is substituted.")48else:49 concl = "positive_effect" if lo > 0 else "negative_effect" if hi < 0 else "no_detectable_effect"50 feas, reason = "complete", None5152res = {"task_id": "F6-conf-02",53 "sample_start": str(d.date.min().date()), "sample_end": str(d.date.max().date()),54 "n_treated": n_treated, "baseline_metric": round(baseline, 6),55 "model_metric": rmse, "effect_estimate": effect, "effect_ci_low": lo, "effect_ci_high": hi,56 "conclusion": concl, "feasibility": feas, "feasibility_reason": reason,57 "data_issues": issues, "n_input_rows": {"gc5y_daily.csv": n_rows},58 "n_analysed": n_an, "random_seed": 0}59json.dump(res, open("out/results.json", "w"), indent=2)6061t1 = d.assign(year=d.date.dt.year).groupby("year").agg(62 n_days=("yield_5y", "size"), mean_yield=("yield_5y", "mean"), mean_change_bp=("chg_bp", "mean")).reset_index()63t1.round(4).to_csv("out/table1.csv", index=False)6465fig, ax = plt.subplots(figsize=(9, 4))66ax.plot(d.date, d.yield_5y, lw=0.8)67ax.axvline(pd.Timestamp("2020-03-04"), color="r", ls="--", label="2020-03-04 cut (outside sample)")68ax.set_xlim(d.date.min(), pd.Timestamp("2020-06-30"))69ax.set_ylabel("GoC 5-year yield (%)"); ax.set_title("GoC 5-year benchmark yield, file sample 2010-2019")70ax.legend(); fig.tight_layout(); fig.savefig("out/figure1.png", dpi=120)71print(json.dumps(res, indent=2))72