PaPP v2: cartella Analisi (script, risultati, metodologia incroci)

Co-Authored-By: Devin AI <158243242+devin-ai-integration[bot]@users.noreply.github.com>
This commit is contained in:
Pietro Giacobazzi
2026-06-17 07:40:41 +00:00
co-authored by Devin AI <158243242+devin-ai-integration[bot]@users.noreply.github.com>
parent 2070d44ce6
commit 19c2e0ac72
11 changed files with 639 additions and 0 deletions
+82
View File
@@ -0,0 +1,82 @@
"""Utility condivise per l'analisi PaPP v2.
I dati di input sono i due CSV prodotti dallo script MQL5 `PaPP_CrossExport.mq5`:
- PaPP_crosses_<SYM>_D1.csv (un record per incrocio)
- PaPP_bars_<SYM>_D1.csv (un record per giorno D1 = baseline)
Per default vengono cercati in ../data. Si possono passare percorsi diversi via argv.
"""
import os, sys, numpy as np, pandas as pd
HERE = os.path.dirname(os.path.abspath(__file__))
DATA = os.path.normpath(os.path.join(HERE, "..", "data"))
HZ = [1, 3, 5, 10, 20]
def load(path):
df = pd.read_csv(path)
df["time"] = pd.to_datetime(df["time"], format="%Y.%m.%d", errors="coerce")
if df["time"].isna().all():
df["time"] = pd.to_datetime(df["time"], errors="coerce")
return df.dropna(subset=["time"]).sort_values("time").reset_index(drop=True)
def load_pair(sym="EURUSD", year_min=1999, data_dir=None):
d = data_dir or DATA
cr = load(os.path.join(d, f"PaPP_crosses_{sym}_D1.csv"))
ba = load(os.path.join(d, f"PaPP_bars_{sym}_D1.csv"))
if year_min:
cr = cr[cr.time.dt.year >= year_min].reset_index(drop=True)
ba = ba[ba.time.dt.year >= year_min].reset_index(drop=True)
return cr, ba
def add_regime(df, edges):
"""Regime = trend(sopra/sotto MA365) x tercile(cluster) x tercile(velocita')."""
tr = (df["trend"] > 0).astype(int)
cl = np.digitize(df["cluster_pct"], edges["cl"])
ve = np.digitize(df["vel_med"], edges["ve"])
out = df.copy()
out["regime"] = tr.astype(str) + "_" + cl.astype(str) + "_" + ve.astype(str)
return out
def regime_edges(ba):
return {
"cl": np.quantile(ba["cluster_pct"], [1/3, 2/3]),
"ve": np.quantile(ba["vel_med"], [1/3, 2/3]),
}
def block_bootstrap_p(x, nb=2000, bs=20, seed=42):
"""p-value bilaterale per media != 0 con block bootstrap (gestisce overlap)."""
rng = np.random.default_rng(seed)
x = np.asarray(x, float)
n = len(x)
if n < 10:
return np.nan
nblocks = int(np.ceil(n / bs))
idx = np.arange(n)
means = np.empty(nb)
for k in range(nb):
starts = rng.integers(0, n, size=nblocks)
pick = np.concatenate([
idx[s:s+bs] if s+bs <= n else np.concatenate([idx[s:], idx[:s+bs-n]])
for s in starts])[:n]
means[k] = x[pick].mean()
return 2 * min((means <= 0).mean(), (means >= 0).mean())
def benjamini_hochberg(p):
"""q-value BH per array di p (NaN ammessi)."""
p = np.asarray(p, float)
out = np.full_like(p, np.nan)
idx = np.where(~np.isnan(p))[0]
if len(idx) == 0:
return out
m = len(idx)
order = idx[np.argsort(p[idx])]
q = p[order] * m / np.arange(1, m + 1)
q = np.minimum.accumulate(q[::-1])[::-1]
out[order] = np.clip(q, 0, 1)
return out
@@ -0,0 +1,69 @@
"""Fase 2 raffinata: extra-rendimento vs baseline + significativita'.
Per ogni (coppia, direzione) calcola:
- rendimento medio grezzo a 1/3/5/10/20 g
- EXTRA-rendimento = evento - baseline nello stesso regime (toglie il bias di periodo)
- effetto in sigma del movimento tipico a h giorni
- quota di eventi che battono la baseline
- p-value (block bootstrap) e q-value (Benjamini-Hochberg)
Uso:
python excess_analysis.py [SYMBOL] [YEAR_MIN]
Output:
../results/summary_<SYM>_<YEAR>.csv + stampa dei top a 5/10/20 g
"""
import os, sys, numpy as np, pandas as pd
from _common import load_pair, add_regime, regime_edges, block_bootstrap_p, benjamini_hochberg, HZ, HERE
SYM = sys.argv[1] if len(sys.argv) > 1 else "EURUSD"
YEAR = int(sys.argv[2]) if len(sys.argv) > 2 else 1999
RESULTS = os.path.normpath(os.path.join(HERE, "..", "results"))
cr, ba = load_pair(SYM, YEAR)
print(f"[{SYM} {YEAR}+] incroci={len(cr)} baseline={len(ba)}")
edges = regime_edges(ba)
cr = add_regime(cr, edges)
ba = add_regime(ba, edges)
base_mean = {h: ba.groupby("regime")[f"cret_{h}"].mean() for h in HZ}
base_glob = {h: ba[f"cret_{h}"].mean() for h in HZ}
base_sd = {h: ba[f"cret_{h}"].std() for h in HZ}
rows = []
for (pair, d), g in cr.groupby(["pair", "dir"]):
n = len(g)
if n < 30:
continue
rec = {"pair": pair, "dir": int(d), "n": n}
for h in HZ:
ev = g[f"cret_{h}"].values
bexp = g["regime"].map(base_mean[h]).fillna(base_glob[h]).values
exc = ev - bexp
rec[f"raw_{h}"] = np.nanmean(ev)
rec[f"exc_{h}"] = np.nanmean(exc)
rec[f"eff_{h}"] = np.nanmean(exc) / base_sd[h]
rec[f"pos_{h}"] = (exc > 0).mean()
rec[f"p_{h}"] = block_bootstrap_p(exc)
rows.append(rec)
res = pd.DataFrame(rows)
for h in HZ:
res[f"q_{h}"] = benjamini_hochberg(res[f"p_{h}"].values)
os.makedirs(RESULTS, exist_ok=True)
outp = os.path.join(RESULTS, f"summary_{SYM}_{YEAR}.csv")
res.to_csv(outp, index=False)
print(f"vol giornaliera prezzo (sd cret_1) = {ba['cret_1'].std():.3f}%")
print(f"baseline drift 5/10/20g = {base_glob[5]:.3f}% / {base_glob[10]:.3f}% / {base_glob[20]:.3f}%")
for h in [5, 10, 20]:
sig = res[res[f"q_{h}"] < 0.10].copy()
sig["ae"] = sig[f"eff_{h}"].abs()
sig = sig.sort_values("ae", ascending=False).head(12)
print(f"\n=== TOP {h}g (q<0.10) ===")
for _, r in sig.iterrows():
print(f" {r['pair']:13s} dir={int(r['dir']):+d} n={int(r['n']):5d} "
f"exc={r[f'exc_{h}']:+.3f}% eff={r[f'eff_{h}']:+.2f}sd "
f"pos={r[f'pos_{h}']*100:.0f}% q={r[f'q_{h}']:.3f}")
print(f"\nsalvato: {outp}")
+83
View File
@@ -0,0 +1,83 @@
"""Verifiche di correttezza + tabella completa 36 coppie x 2 direzioni.
Verifiche:
V1 la direzione registrata coincide col segno di (A-B) all'incrocio
V2 cret_1 coincide con ret_1
V3 ricalcolo indipendente della media di un pattern
V4 baseline drift ~0
Tabella: per ogni (coppia, dir) extra-rendimento 5/10/20g, quota che batte il baseline,
q-value a 10g, tasso di ritorno alla Mediana (rev10/rev20), giorni mediani al ritorno,
ed etichetta verbale "dopo 10g".
Uso:
python full_table.py [SYMBOL] [YEAR_MIN]
Output:
../results/tabella_incroci_<SYM>_<YEAR>.csv (+ .md)
"""
import os, sys, numpy as np, pandas as pd
from _common import load_pair, add_regime, regime_edges, block_bootstrap_p, benjamini_hochberg, HZ, HERE
SYM = sys.argv[1] if len(sys.argv) > 1 else "EURUSD"
YEAR = int(sys.argv[2]) if len(sys.argv) > 2 else 1999
RESULTS = os.path.normpath(os.path.join(HERE, "..", "results"))
cr, ba = load_pair(SYM, YEAR)
# ---------------- VERIFICHE ----------------
sub = cr[cr.pair == "PRICExMED"]
chk = ((sub.price - sub.med) > 0).astype(int).replace(0, -1)
print(f"[V1] PRICExMED dir vs segno(price-med): {(np.sign(sub.dir)==np.sign(chk)).mean()*100:.1f}% (atteso ~100%)")
print(f"[V2] max|cret_1 - ret_1| = {(cr.cret_1-cr.ret_1).abs().max():.6f} (atteso 0)")
g0 = cr[(cr.pair=='MA121xMA7') & (cr.dir==1)]
print(f"[V3] MA121xMA7+ n={len(g0)} media cret_10 = {g0.cret_10.mean():.4f}%")
print(f"[V4] baseline cret_10 medio={ba.cret_10.mean():.4f}% sd={ba.cret_10.std():.3f}% (drift ~0)")
# ---------------- TABELLA ----------------
edges = regime_edges(ba)
cr = add_regime(cr, edges); ba = add_regime(ba, edges)
bmean = {h: ba.groupby("regime")[f"cret_{h}"].mean() for h in HZ}
bglob = {h: ba[f"cret_{h}"].mean() for h in HZ}
bsd = {h: ba[f"cret_{h}"].std() for h in HZ}
rows = []
for (pair, d), g in cr.groupby(["pair", "dir"]):
rec = {"pair": pair, "dir": int(d), "n": len(g)}
for h in HZ:
ev = g[f"cret_{h}"].values
bexp = g["regime"].map(bmean[h]).fillna(bglob[h]).values
exc = ev - bexp
rec[f"raw_{h}"] = np.nanmean(ev)
rec[f"exc_{h}"] = np.nanmean(exc)
rec[f"eff_{h}"] = np.nanmean(exc) / bsd[h]
rec[f"pos_{h}"] = (exc > 0).mean() * 100
rec[f"p_{h}"] = block_bootstrap_p(exc)
rec["rev10"] = g["rev_10"].mean() * 100
rec["rev20"] = g["rev_20"].mean() * 100
rec["btr_med"] = g["bars_to_revert"].median()
rows.append(rec)
res = pd.DataFrame(rows)
res["q10"] = benjamini_hochberg(res["p_10"].values)
def verb(r):
if r["q10"] < 0.10 and r["exc_10"] > 0: return "SALE (extra +)"
if r["q10"] < 0.10 and r["exc_10"] < 0: return "SCENDE (extra -)"
return "neutro"
res["dopo_10g"] = res.apply(verb, axis=1)
res = res.sort_values(["q10", "exc_10"]).reset_index(drop=True)
os.makedirs(RESULTS, exist_ok=True)
csvp = os.path.join(RESULTS, f"tabella_incroci_{SYM}_{YEAR}.csv")
res.to_csv(csvp, index=False)
md = [f"# {SYM} (D1, {YEAR}+) — cosa succede dopo ogni incrocio (36 coppie x 2 dir)", "",
"extra5/10/20 = extra-rendimento % vs baseline stesso regime; pos10 = % che batte il baseline;",
"q10 = significativita' (BH, <0.10 robusto); rev10/20 = % ritorno alla Mediana; btr = giorni mediani.", "",
"| coppia | dir | n | extra5 | extra10 | extra20 | pos10 | q10 | rev10 | rev20 | btr | dopo 10g |",
"|---|---|---|---|---|---|---|---|---|---|---|---|"]
for _, x in res.iterrows():
md.append(f"| {x['pair']} | {int(x['dir']):+d} | {int(x['n'])} | {x['exc_5']:+.3f} | "
f"{x['exc_10']:+.3f} | {x['exc_20']:+.3f} | {x['pos_10']:.0f}% | {x['q10']:.3f} | "
f"{x['rev10']:.0f}% | {x['rev20']:.0f}% | {x['btr_med']:.0f} | {x['dopo_10g']} |")
open(os.path.join(RESULTS, f"tabella_incroci_{SYM}_{YEAR}.md"), "w").write("\n".join(md))
print(f"\nsalvato: {csvp} (+ .md) righe={len(res)}")
@@ -0,0 +1,43 @@
"""Grafico delle traiettorie di extra-rendimento (giorno per giorno) dopo l'incrocio.
Disegna l'extra-rendimento cumulato (cret_1..20 - baseline stesso regime) per un
insieme di pattern (default: i piu' consistenti legati alla MA121/MA365).
Uso:
python plot_trajectories.py [SYMBOL] [YEAR_MIN]
Output:
../results/trajectories_<SYM>_<YEAR>.png
"""
import os, sys, numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from _common import load_pair, add_regime, regime_edges, HERE
SYM = sys.argv[1] if len(sys.argv) > 1 else "EURUSD"
YEAR = int(sys.argv[2]) if len(sys.argv) > 2 else 1999
RESULTS = os.path.normpath(os.path.join(HERE, "..", "results"))
HZ = list(range(1, 21))
cr, ba = load_pair(SYM, YEAR)
edges = regime_edges(ba)
cr = add_regime(cr, edges); ba = add_regime(ba, edges)
bmean = {h: ba.groupby("regime")[f"cret_{h}"].mean() for h in HZ}
bglob = {h: ba[f"cret_{h}"].mean() for h in HZ}
sel = [("MA121xMA7", 1), ("MEDxMA121", -1), ("PRICExMA121", -1), ("MA121xMA3", 1), ("MA365xMA7", 1)]
plt.figure(figsize=(11, 6))
for pair, d in sel:
g = cr[(cr.pair == pair) & (cr.dir == d)]
if len(g) < 20:
continue
exc = [(g[f"cret_{h}"] - g["regime"].map(bmean[h]).fillna(bglob[h])).mean() for h in HZ]
plt.plot(HZ, exc, marker="o", ms=3, label=f"{pair} dir={d:+d} (n={len(g)})")
plt.axhline(0, color="k", lw=.8)
plt.title(f"{SYM} ({YEAR}+) - extra-rendimento medio vs baseline (stesso regime)")
plt.xlabel("giorni dopo l'incrocio"); plt.ylabel("extra-rendimento cumulato (%)")
plt.legend(); plt.grid(alpha=.3); plt.tight_layout()
os.makedirs(RESULTS, exist_ok=True)
outp = os.path.join(RESULTS, f"trajectories_{SYM}_{YEAR}.png")
plt.savefig(outp, dpi=130)
print("salvato:", outp)