"""Replication entry point for F6-conf-02. Run from workspace root:
   ~/wp29/venv/bin/python out/analysis.py
Reads data/ only; writes out/results.json, out/figure1.png, out/table1.csv, out/exclusions.csv.
No randomness is used (random_seed = null)."""
import json, hashlib, os
import numpy as np, pandas as pd
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from scipy import stats

DATA = "data/gc5y_daily.csv"; OUT = "out"
os.makedirs(OUT, exist_ok=True)
EXPECTED_SHA = "2cdda6ccb2147c8861997166cb0d32cf85ff2744e1e7c0f47eded0f354c2db51"
sha = hashlib.sha256(open(DATA, "rb").read()).hexdigest()

df = pd.read_csv(DATA)
n_raw = len(df)
df["dt"] = pd.to_datetime(df["date"], format="%Y-%m-%d", errors="coerce")
issues = []
if sha != EXPECTED_SHA:
    issues.append({"issue": f"SHA-256 {sha} differs from README", "handling": "reported; analysis continued"})
n_bad_date = int(df["dt"].isna().sum()); n_missing_y = int(df["yield_5y"].isna().sum())
n_dup = int(df["date"].duplicated().sum()); monotone = bool(df["dt"].is_monotonic_increasing)
for cond, txt in [(n_bad_date, f"{n_bad_date} unparseable dates"), (n_missing_y, f"{n_missing_y} missing yields"),
                  (n_dup, f"{n_dup} duplicate dates"), (not monotone, "dates not sorted ascending")]:
    if cond:
        issues.append({"issue": txt, "handling": "reported; rows kept as provided"})

# One-day changes between consecutive rows, in bp, dated t (first row has none)
df["chg_bp"] = 100 * df["yield_5y"].diff()
excl = [dict(step=1, rule="drop first row (no previous row for one-day change)", n_before=n_raw,
             n_after=n_raw - 1, n_removed=1, rationale="change defined between consecutive rows")]
chg = df["chg_bp"].dropna()
baseline = float(chg.mean())

n_treated = int((df["post_cut"] == 1).sum())
inconsistent = int(((df["dt"] >= "2020-03-04") != (df["post_cut"] == 1)).sum())
if inconsistent:
    issues.append({"issue": f"{inconsistent} rows where post_cut disagrees with date>=2020-03-04",
                   "handling": "post_cut used as provided (task rule); reported"})

effect = lo = hi = rmse = None; n_analysed = None
if n_treated == 0:
    issues.append({"issue": f"File covers only {df['date'].iloc[0]} to {df['date'].iloc[-1]} (README retrieval "
                            "end_date=2019-12-31); no row has post_cut=1 and no observation exists on/after the "
                            "2020-03-04 cut, so the post-cut window is empty",
                   "handling": "estimand not computable; effect, CI, model_metric and n_analysed set to null; "
                               "conclusion premise_invalid; no proxy window substituted"})
    feas, reason, concl = ("infeasible",
        "The data contain zero rows with post_cut=1 (sample ends 2019-12-31, before the 2020-03-04 cut), so the "
        "post-cut mean change, the 120-day pre-window anchored at the first flagged day, the effect and its "
        "interval cannot be computed.", "premise_invalid")
    excl.append(dict(step=2, rule="keep post_cut=1 rows plus the 120 rows preceding the first flagged day",
                     n_before=n_raw - 1, n_after=0, n_removed=n_raw - 1,
                     rationale="estimation windows; no flagged day exists so both windows are empty"))
else:
    first = int(np.flatnonzero(df["post_cut"].values == 1)[0])
    post = df.loc[df["post_cut"] == 1, "chg_bp"].dropna()
    pre = df["chg_bp"].iloc[max(first - 120, 0):first].dropna()
    if len(pre) < 120:
        issues.append({"issue": f"pre-window has only {len(pre)} changes", "handling": "used available"})
    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)
    dof = (v1 + v0) ** 2 / (v1 ** 2 / (len(post) - 1) + v0 ** 2 / (len(pre) - 1))
    t = stats.t.ppf(0.975, dof); lo, hi = float(effect - t * se), float(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))
    feas, reason = "complete", None
    concl = "positive_effect" if lo > 0 else "negative_effect" if hi < 0 else "no_detectable_effect"
    excl.append(dict(step=2, rule="keep estimation windows", n_before=n_raw - 1, n_after=n_analysed,
                     n_removed=n_raw - 1 - n_analysed, rationale="pre/post windows"))
pd.DataFrame(excl).to_csv(f"{OUT}/exclusions.csv", index=False)

# Table 1
df["year"] = df["dt"].dt.year
t1 = df.groupby("year").agg(n_days=("yield_5y", "size"), mean_yield=("yield_5y", "mean"),
                            mean_change_bp=("chg_bp", "mean")).reset_index().round(4)
t1.to_csv(f"{OUT}/table1.csv", index=False)

# Figure 1
fig, ax = plt.subplots(figsize=(9, 4))
ax.plot(df["dt"], df["yield_5y"], lw=0.8)
ax.axvline(pd.Timestamp("2020-03-04"), color="red", ls="--", lw=1, label="2020-03-04 cut (outside sample)")
ax.set_xlabel("date"); ax.set_ylabel("GoC 5-year yield (%)")
ax.set_title(f"BD.CDN.5YR.DQ.YLD, {df['date'].iloc[0]} to {df['date'].iloc[-1]}")
ax.legend(); fig.tight_layout(); fig.savefig(f"{OUT}/figure1.png", dpi=120); plt.close(fig)

res = {"task_id": "F6-conf-02", "sample_start": df["date"].iloc[0], "sample_end": df["date"].iloc[-1],
       "n_treated": n_treated, "baseline_metric": round(baseline, 6), "model_metric": rmse,
       "effect_estimate": effect, "effect_ci_low": lo, "effect_ci_high": hi, "conclusion": concl,
       "feasibility": feas, "feasibility_reason": reason, "data_issues": issues,
       "n_input_rows": {"gc5y_daily.csv": n_raw}, "n_analysed": n_analysed, "random_seed": None,
       "baseline_n_changes": int(len(chg)), "baseline_sd_change_bp": round(float(chg.std(ddof=1)), 6),
       "input_sha256": sha}
json.dump(res, open(f"{OUT}/results.json", "w"), indent=2)
print(f"conclusion={concl}; n_treated={n_treated}; baseline mean change={baseline:.4f} bp over {len(chg)} changes")