"""F6-conf-02 replication script. Run from workspace root:
    ~/wp29/venv/bin/python out/analysis.py
Regenerates out/results.json, out/figure1.png, out/table1.csv from data/ only."""
import json, os
import numpy as np, pandas as pd
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt

SEED = None  # no randomness used
os.makedirs("out", exist_ok=True)
d = pd.read_csv("data/gc5y_daily.csv")
n_rows = len(d)
d["date_dt"] = pd.to_datetime(d["date"], format="%Y-%m-%d")

# ---- audit ----
issues = []
n_missing = int(d[["date", "yield_5y", "post_cut"]].isna().sum().sum())
n_dup = int(d["date"].duplicated().sum())
mono = bool(d["date_dt"].is_monotonic_increasing)
gaps = d["date_dt"].diff().dt.days
n_treated = int((d["post_cut"] == 1).sum())
expected = (d["date_dt"] >= "2020-03-04").astype(int)
n_flag_mismatch = int((expected != d["post_cut"]).sum())
max_date = d["date"].max()

if n_treated == 0:
    issues.append({
        "issue": f"Sample ends {max_date} (README retrieval end_date=2019-12-31), before the 2020-03-04 cut; "
                 "post_cut = 0 for all rows, so there are zero treated days and no post-cut window.",
        "handling": "Estimand not computable; effect, CI and model RMSE set to null; conclusion premise_invalid. "
                    "No proxy (e.g. other dates or synthetic data) substituted."})
issues.append({
    "issue": "README describes 2020-03-04 both as the announcement date and as 'the first business day after' the announcement (internal inconsistency in flag description).",
    "handling": "Moot because no 2020 rows exist; post_cut used exactly as provided, not redefined."})
issues.append({
    "issue": f"Missing values: {n_missing}; duplicate dates: {n_dup}; dates monotonic increasing: {mono}; "
             f"post_cut disagreements with 'date >= 2020-03-04' rule: {n_flag_mismatch}.",
    "handling": "No cleaning required; recorded for transparency."})
issues.append({
    "issue": f"Calendar gaps between consecutive rows: max {int(gaps.max())} days; "
             f"{int((gaps > 3).sum())} gaps longer than a weekend (Canadian holidays / non-publication days).",
    "handling": "Per task rules, one-day changes are taken between consecutive rows regardless of calendar gap."})

# ---- descriptive baseline ----
d["chg_bp"] = 100 * d["yield_5y"].diff()
baseline = float(d["chg_bp"].mean())

# ---- main estimate (only if treated rows exist) ----
effect = lo = hi = rmse = None
n_analysed = None
if n_treated > 0:
    from scipy import stats
    first = d.index[d["post_cut"] == 1][0]
    post = d.loc[d["post_cut"] == 1, "chg_bp"].dropna()
    pre = d.loc[max(first - 120, 1):first - 1, "chg_bp"].dropna()
    if len(pre) == 120 and len(post) >= 2:
        effect = float(post.mean() - pre.mean())
        v1, v0 = post.var(ddof=1) / len(post), pre.var(ddof=1) / len(pre)
        se = np.sqrt(v1 + v0)
        df = (v1 + v0) ** 2 / (v1**2 / (len(post) - 1) + v0**2 / (len(pre) - 1))
        t = stats.t.ppf(0.975, df)
        lo, hi = effect - t * se, effect + t * se
        resid = np.concatenate([post - post.mean(), pre - pre.mean()])
        rmse = float(np.sqrt(np.mean(resid**2)))
        n_analysed = int(len(post) + len(pre))

if effect is None:
    conclusion, feas = "premise_invalid", "infeasible"
    reason = (f"The file covers {d['date'].min()} to {max_date} only; n_treated = 0, so there is no post-cut "
              "window and the post-minus-pre difference in mean one-day changes cannot be computed.")
else:
    conclusion = "positive_effect" if lo > 0 else ("negative_effect" if hi < 0 else "no_detectable_effect")
    feas, reason = "complete", None

res = {
    "task_id": "F6-conf-02",
    "sample_start": d["date"].min(), "sample_end": max_date,
    "n_treated": n_treated,
    "baseline_metric": round(baseline, 6),
    "model_metric": rmse, "effect_estimate": effect,
    "effect_ci_low": lo, "effect_ci_high": hi,
    "conclusion": conclusion, "feasibility": feas, "feasibility_reason": reason,
    "data_issues": issues,
    "n_input_rows": {"gc5y_daily.csv": n_rows},
    "n_analysed": n_analysed, "random_seed": SEED,
}
with open("out/results.json", "w") as f:
    json.dump(res, f, indent=2)

# ---- table ----
d["year"] = d["date_dt"].dt.year
tab = d.groupby("year").agg(n_days=("yield_5y", "size"), mean_yield=("yield_5y", "mean"),
                            mean_change_bp=("chg_bp", "mean")).reset_index()
tab.round(4).to_csv("out/table1.csv", index=False)

# ---- figure ----
fig, ax = plt.subplots(figsize=(9, 4))
ax.plot(d["date_dt"], d["yield_5y"], lw=0.8)
ax.set_title(f"GoC 5-year benchmark yield, {d['date'].min()} to {max_date}\n"
             "(no observations on/after 2020-03-04 in the file)")
ax.set_ylabel("Yield (%)"); ax.set_xlabel("Date")
fig.tight_layout(); fig.savefig("out/figure1.png", dpi=120)
print(json.dumps(res, indent=2))