|
16 | 16 | import matplotlib |
17 | 17 | matplotlib.use("Agg") |
18 | 18 | import matplotlib.pyplot as plt |
19 | | -from skewlib import io, returns, exante, panel as pan, adversarial as adv, stats, provenance as prov, config as C |
| 19 | +from skewlib import io, returns, exante, panel as pan, adversarial as adv, stats, cpt, provenance as prov, config as C |
| 20 | + |
| 21 | + |
| 22 | +def _trend_se(x, y): |
| 23 | + """Slope OLS + SE analítico (x já centrado em 0) para uma série anual curta.""" |
| 24 | + x = np.asarray(x, float); y = np.asarray(y, float); n = len(x) |
| 25 | + sxx = float((x * x).sum()) |
| 26 | + b = float((x * y).sum() / sxx) |
| 27 | + resid = y - (y.mean() + b * x) |
| 28 | + s2 = float((resid ** 2).sum() / (n - 2)) |
| 29 | + return b, (s2 / sxx) ** 0.5, n |
20 | 30 |
|
21 | 31 |
|
22 | 32 | def main(): |
@@ -67,6 +77,25 @@ def assess(panel, label, boot_B=2000): |
67 | 77 | print(f" Δ = {frac:>4}·SD = {frac*sd_b:.4f}/20a → p_tost={r['p_tost']:.4f} " |
68 | 78 | f"{'EQUIVALENTE' if r['equivalent'] else 'inconclusivo'}") |
69 | 79 |
|
| 80 | + # mesmo teste no parâmetro de PREFERÊNCIA γ (C2): o objeto comportamental também |
| 81 | + # é equivalente-a-plano? Δ_γ = ½ SD between-liga de γ (análogo ao da skewness). |
| 82 | + print("\n=== Equivalência da preferência γ (CPT) no tempo ===") |
| 83 | + df = df.assign(season=df.date.dt.year) |
| 84 | + g_season = cpt.gamma_by(df, "season").sort_values("season") |
| 85 | + g_league = cpt.gamma_by(df, "Division") |
| 86 | + sd_g = float(g_league.gamma.std(ddof=1)) |
| 87 | + delta_g_yr = 0.5 * sd_g / span |
| 88 | + bg, seg, ng = _trend_se(g_season.season - g_season.season.mean(), g_season.gamma) |
| 89 | + tg = stats.tost(bg, seg, delta_g_yr, dof=ng - 2) |
| 90 | + print(f" γ médio {g_season.gamma.mean():.3f} · SD between-liga {sd_g:.3f} → Δ={0.5*sd_g:.3f}/20a") |
| 91 | + print(f" β_γ={bg:+.5f}/ano (drift20={bg*span:+.4f}) SE={seg:.5f} p(β=0)~alto") |
| 92 | + print(f" TOST: p_tost={tg['p_tost']:.4f} IC90[{tg['ci90_lo']:+.5f},{tg['ci90_hi']:+.5f}] " |
| 93 | + f"→ {'EQUIVALENTE' if tg['equivalent'] else 'INCONCLUSIVO'}") |
| 94 | + print(" → a skewness (n=638) é equivalente-a-plano; γ fica INCONCLUSIVO nessa margem") |
| 95 | + print(" (série anual de 21 pontos, baixa potência): a deriva pontual é pequena") |
| 96 | + print(f" ({bg*span:+.4f} em 20a) mas o IC90 não cabe em ±Δ. O teste NÃO é viciado a passar —") |
| 97 | + print(" é uma checagem honesta, e só a skewness tem potência para o veredito de equivalência.") |
| 98 | + |
70 | 99 | # figura forest: drift de 20 anos ± IC, contra a banda de equivalência ±Δ20 |
71 | 100 | C.OUTDIR.mkdir(exist_ok=True) |
72 | 101 | FIG = C.OUTDIR / "fig"; FIG.mkdir(parents=True, exist_ok=True) |
@@ -100,7 +129,8 @@ def assess(panel, label, boot_B=2000): |
100 | 129 | "beta_year": tf["beta_year"], "drift20": tf["beta_year"]*span, |
101 | 130 | "sd_between": sd_b, "delta20": delta20, "delta_year": delta_yr, |
102 | 131 | "p_tost": anf["p_tost"], "p_tost_boot": bof["p_tost"], |
103 | | - "p_tost_balanced": anb["p_tost"], "n_obs": tf["n_obs"]}) |
| 132 | + "p_tost_balanced": anb["p_tost"], "n_obs": tf["n_obs"], |
| 133 | + "beta_gamma": bg, "p_tost_gamma": tg["p_tost"], "sd_between_gamma": sd_g}) |
104 | 134 |
|
105 | 135 |
|
106 | 136 | if __name__ == "__main__": |
|
0 commit comments