#!/usr/bin/env python3
"""score_cue_test.py — the pre-registered test of the p-value cue and the sample-size cue on SCORE
(Tyner et al., Nature 2026; OSF dtzx4, CC0). survives-a-recheck, Stage 5.

SCORE re-ran claims sampled from 54 social- and behavioural-science journals (2009-2018) at a median 99.6% power.
The piece's published falsifiers name exactly this kind of batch. Three tests, each reusing the FReD test's own
estimator (paper-clustered bootstrap, each paper weighted once; `fred_two_cues.diff_with_ci`) and verdict rules:

  A  the falsifier as published: originals in band S (p < .001) against band W (.01-.05), bands by the FReD rule
     (`fred_two_cues.p_band`; a bare threshold such as "p < .05" cannot be placed and is left out)
  B  the card's first sign as a reader applies it: printed exact p if given; a threshold at or below .001 is strong;
     otherwise p from the test statistic (t with df, z, F with one numerator df, chi-squared with one df); with no
     usable statistic the sign counts (weak). Strong against weak.
  C  sample size: the biggest third of original samples against the smallest third (the piece's second falsifier)

Rules, outcomes and what the piece does with each: `PREREG — survives-a-recheck — SCORE cue test — 2026-10-27.md`.

    python3 score_cue_test.py            # real outcomes (repli_outcomes.csv must be present and hash-pinned)
    python3 score_cue_test.py --dummy    # never reads the replication file: random outcomes on every original claim
    python3 score_cue_test.py --plant    # never reads it either: strong-band claims replicate at 65%, the rest at 45%
"""
from __future__ import annotations

import hashlib
import json
import os
import sys

import numpy as np
import pandas as pd
from scipy import stats

HERE = os.path.dirname(os.path.abspath(__file__))
sys.path.insert(0, HERE)
import fred_two_cues as F  # noqa: E402  (same band rule, estimator and verdicts as the FReD test)

DATA = os.path.join(HERE, "data", "score")
ORIG = os.path.join(DATA, "orig_outcomes.csv")
META = os.path.join(DATA, "paper_metadata.csv")
REPLI = os.path.join(DATA, "repli_outcomes.csv")
EXCL = os.path.join(DATA, "repli_exclusions.csv")
ORIG_SHA = "325b445d0a253342eb3d2d318aa4233bcca23c16866b5c1cc2733bab8f26478f"
META_SHA = "0f2beef44537aa2f123f8045281f91cc8d524186aec98e51909e7ea5934c660e"
SEED = 20261027
TYPE_MAP = {"exact": "=", "less-than": "<", "greater-than": ">"}  # 'threshold' and blanks cannot be placed


def sha(p: str) -> str:
    return hashlib.sha256(open(p, "rb").read()).hexdigest()


def truthy(x) -> float:
    s = str(x).strip().lower()
    if s in ("true", "1", "1.0", "yes", "t"):
        return 1.0
    if s in ("false", "0", "0.0", "no", "f"):
        return 0.0
    return np.nan


def p_from_stat(row) -> float:
    """Two-tailed p from the reported statistic, or NaN when the card would not use it."""
    st, v = row["orig_stat_type"], row["orig_stat_value"]
    d1, d2 = row["orig_stat_dof_1"], row["orig_stat_dof_2"]
    if pd.isna(v):
        return np.nan
    v = abs(float(v))
    if st == "t":
        if pd.isna(d1):
            return 0.0 if v > 3.65 else 1.0     # no df: the card's most demanding row (about 30 df)
        return float(2 * stats.t.sf(v, d1))
    if st == "z":
        return float(2 * stats.norm.sf(v))
    if st == "F" and d1 == 1 and not pd.isna(d2):
        return float(2 * stats.t.sf(np.sqrt(v), d2))
    if st == "chi_squared":
        df = d1 if not pd.isna(d1) else d2
        if df == 1:
            return float(2 * stats.norm.sf(np.sqrt(v)))
    return np.nan


def card_strong(row) -> float:
    """The card's first sign: 1 = strong pass, 0 = weak pass (sign counts), NaN = original not significant."""
    t, p = row["orig_p_value_type"], row["orig_p_value"]
    if t == "greater-than" or (t == "exact" and not pd.isna(p) and float(p) > 0.05):
        return np.nan
    if t == "exact" and not pd.isna(p):
        return 1.0 if float(p) < 0.001 else 0.0
    if t in ("less-than", "threshold") and not pd.isna(p) and float(p) <= 0.001:
        return 1.0
    ps = p_from_stat(row)
    if not np.isnan(ps):
        return 1.0 if ps < 0.001 else 0.0
    return 0.0


def load(mode: str) -> pd.DataFrame:
    for p, h in ((ORIG, ORIG_SHA), (META, META_SHA)):
        if sha(p) != h:
            print(f"HASH MISMATCH: {os.path.basename(p)}")
            sys.exit(1)
    o = pd.read_csv(ORIG, low_memory=False)
    meta = pd.read_csv(META, low_memory=False, usecols=["paper_id", "COS_pub_category"])
    o = o.merge(meta, on="paper_id", how="left")
    o["band"] = [F.p_band(TYPE_MAP.get(t), v) for t, v in zip(o["orig_p_value_type"], o["orig_p_value"])]
    o["strong"] = o.apply(card_strong, axis=1)
    o["n_o"] = pd.to_numeric(o["orig_sample_size_value"], errors="coerce")
    o["doi_o"] = o["paper_id"]                       # the FReD estimator clusters on doi_o
    rng = np.random.default_rng(SEED)
    if mode == "dummy":
        o["y"] = rng.integers(0, 2, size=len(o)).astype(float)
        return o
    if mode == "plant":
        o["y"] = (rng.random(len(o)) < np.where(o["band"] == "S", 0.65, 0.45)).astype(float)
        return o
    r = pd.read_csv(REPLI, low_memory=False)
    n_rows = len(r)
    if os.path.exists(EXCL):
        ex = pd.read_csv(EXCL, low_memory=False)
        key = "report_id" if "report_id" in ex.columns else ("claim_id" if "claim_id" in ex.columns else None)
        if key:
            r = r[~r[key].isin(ex[key])]
    if "repli_version_of_record" in r.columns:
        vor = r["repli_version_of_record"].map(truthy)
        r = r[(vor == 1.0) | (~r["claim_id"].isin(r.loc[vor == 1.0, "claim_id"]))]
    r = r.sort_values("report_id").drop_duplicates("claim_id", keep="first")
    r["y"] = r["repli_score_criteria_met"].map(truthy)
    print(f"replication file: {n_rows} rows -> {len(r)} claims after exclusions and one version per claim; "
          f"{int(r['y'].isna().sum())} with no usable outcome")
    m = o.merge(r[["claim_id", "y"]], on="claim_id", how="inner")
    return m.dropna(subset=["y"])


def main() -> int:
    mode = "dummy" if "--dummy" in sys.argv else ("plant" if "--plant" in sys.argv else "real")
    df = load(mode)
    rng = np.random.default_rng(SEED)
    res = {"mode": mode, "claims": int(len(df)), "papers": int(df["paper_id"].nunique()),
           "band_counts": df["band"].value_counts().to_dict(),
           "card_strong_counts": df["strong"].value_counts(dropna=False).to_dict()}
    print(f"SCORE cue test · mode {mode} · {res['claims']} claims from {res['papers']} papers")
    print(f"bands (FReD rule): {res['band_counts']}")

    def block(name, d):
        out = {}
        a = F.run_h1(d, rng)
        out["A_S_vs_W"] = dict(a, verdict=F.h1_verdict(a))
        out["A_middle_band_rate"] = F.diff_with_ci(d[d["band"] == "M"], d[d["band"] == "W"], rng)["rate_a"] if (d["band"] == "M").any() else None
        b = F.diff_with_ci(d[d["strong"] == 1.0], d[d["strong"] == 0.0], rng)
        out["B_card_strong_vs_weak"] = dict(b, verdict=F.h1_verdict(b))
        sig = d[d["band"].isin(F.SIGNIFICANT_BANDS) | d["strong"].notna()]
        c = F.run_h2(sig.dropna(subset=["n_o"]), rng)
        out["C_size_top_vs_bottom"] = dict(c, verdict=F.h2_verdict(c))
        out["A_claim_weighted"] = F.diff_with_ci(d[d["band"] == "S"], d[d["band"] == "W"], rng, weighted=False)
        print(f"\n[{name}]")
        for k in ("A_S_vs_W", "B_card_strong_vs_weak", "C_size_top_vs_bottom"):
            v = out[k]
            print(f"  {k}: {F.pct(v['rate_a'])} ({v['n_a']} claims, {v['papers_a']} papers) vs "
                  f"{F.pct(v['rate_b'])} ({v['n_b']}, {v['papers_b']}) · {F.pp(v['diff'])} "
                  f"[{F.pp(v['ci95'][0])}, {F.pp(v['ci95'][1])}] -> {v['verdict']}")
        cw = out["A_claim_weighted"]
        print(f"  A claim-weighted: {F.pct(cw['rate_a'])} vs {F.pct(cw['rate_b'])} · {F.pp(cw['diff'])} "
              f"[{F.pp(cw['ci95'][0])}, {F.pp(cw['ci95'][1])}]")
        if out["A_middle_band_rate"] is not None:
            print(f"  A middle band (.001-.01) rate: {F.pct(out['A_middle_band_rate'])}")
        return out

    res["all"] = block("all fields (primary)", df)
    ph = df[df["COS_pub_category"] == "psychology and health"]
    if len(ph) >= 20:
        res["psychology_and_health"] = block("psychology and health (secondary)", ph)
    a = res["all"]["A_S_vs_W"]["diff"]
    res["falsifier_p_cue"] = ("FIRED: strong passes failed as often as narrow ones" if a <= 0 else
                              res["all"]["A_S_vs_W"]["verdict"])
    res["falsifier_sample_size"] = ("FIRED" if res["all"]["C_size_top_vs_bottom"]["verdict"].startswith("SIZE HELPS (")
                                    else "not fired")
    print(f"\nfalsifier, p-value cue: {res['falsifier_p_cue']}")
    print(f"falsifier, sample size: {res['falsifier_sample_size']}")
    dst = os.path.join(HERE, f"score-cue-test-results{'' if mode == 'real' else '.' + mode}.json")
    json.dump(res, open(dst, "w"), indent=2, default=float)
    print(f"wrote {os.path.basename(dst)}")
    return 0


if __name__ == "__main__":
    sys.exit(main())