"""F5-conf-03: counterexample to 'X, Y normal => X+Y normal'.
Construction: X ~ N(0,1), W independent Rademacher (+/-1 w.p. 1/2), Y = W*X.
Run from workspace root: ~/wp29/venv/bin/python out/analysis.py
"""
import json, os
import numpy as np, pandas as pd
from scipy import stats

props = pd.read_csv("data/propositions.csv")
design = pd.read_csv("data/design.csv")
issues = []
if props["id"].duplicated().any() or set(props["id"]) != {"A", "B"}:
    issues.append({"issue": "unexpected proposition ids", "handling": "used ids A and B as given"})
d = dict(zip(design["parameter"], design["value"]))
seed, n = int(d["seed"]), int(d["n"])

rng = np.random.default_rng(seed)
x = rng.standard_normal(n)
w = rng.choice(np.array([-1.0, 1.0]), size=n)
y = w * x
s = x + y
os.makedirs("out", exist_ok=True)
pd.DataFrame({"x": x, "y": y}).to_csv("out/counterexample.csv", index=False)

m_s, sd_s = s.mean(), s.std(ddof=1)
ks = stats.kstest(s, "norm", args=(m_s, sd_s)).statistic
# marginal checks (illustration only)
ks_x = stats.kstest(x, "norm", args=(x.mean(), x.std(ddof=1))).statistic
ks_y = stats.kstest(y, "norm", args=(y.mean(), y.std(ddof=1))).statistic
frac_zero = float(np.mean(s == 0.0))

res = {
    "task_id": "F5-conf-03",
    "statement_true": False,
    "counterexample_provided": True,
    "counterexample_values": {
        "construction": "X ~ N(0,1); W independent of X with P(W=1)=P(W=-1)=1/2; Y = W*X. "
                        "Then Y ~ N(0,1) by symmetry, but X+Y = (1+W)X equals 0 with probability 1/2 "
                        "and 2X (N(0,4)) with probability 1/2: a non-degenerate law with an atom at 0, hence not normal.",
        "n": int(n), "seed": seed,
        "mean_x": float(x.mean()), "sd_x": float(x.std(ddof=1)),
        "mean_y": float(y.mean()), "sd_y": float(y.std(ddof=1)),
        "mean_sum": float(m_s), "sd_sum": float(sd_s),
        "ks_sum_vs_normal": float(ks),
    },
    "statement_b_true": True,
    "feasibility": "complete",
    "feasibility_reason": None,
    "data_issues": issues,
    "n_input_rows": {"propositions.csv": int(len(props)), "design.csv": int(len(design))},
    "n_analysed": int(n),
    "random_seed": seed,
}
json.dump(res, open("out/results.json", "w"), indent=2)
print(json.dumps(res, indent=2))
print(f"aux: ks_x={ks_x:.5f} ks_y={ks_y:.5f} frac_sum_exactly_zero={frac_zero:.5f} "
      f"KS 5% crit approx {1.358/np.sqrt(n):.5f}")