PythagorasGoal/scripts/research/tp01_dvol_overlay.py

"""TP01 × DVOL — la vol IMPLICITA (forward-looking) migliora il risk-sizing di TP01? (ESEGUIBILE)

A differenza degli sleeve diversificatori (XS01/VRP01/carry = STAT-MODE, non eseguibili a $600),
questo TOCCA il book live: TP01 è BTC/ETH perp su Deribit, già armato. Oggi vol-targeta sulla vol
REALIZZATA 30g (backward-looking). Ipotesi: il DVOL (vol implicita 30g Deribit, forward-looking,
che spesso ANTICIPA i salti di vol realizzata) come denominatore del vol-target → de-risking più
tempestivo prima dei crash → hold-out migliore e/o DD più basso, SENZA peggiorare il FULL.

Onestà: DVOL parte 2021-03 → confronto TUTTE le varianti sulla FINESTRA COMUNE 2021-2026 (perdo
2018-2021, incluso il toro 2021 pre-DVOL). Baseline ricalcolato sulla stessa finestra. Hold-out 2025+.
Tutto causale (vol/segnale ≤ close[i]), fee 0.10% RT, long-flat, leva cap 2x — config CANONICA TP01.

VARIANTI (denominatore del vol-target):
  REALIZED  -> 30g realizzata (baseline canonica)
  DVOL      -> DVOL/100 (implicita)
  BLEND     -> 0.5·realizzata + 0.5·DVOL
  MAX       -> max(realizzata, DVOL)  (sizing più difensivo: la più alta delle due)
  DERISK    -> realizzata, ma posizione ×0.5 quando DVOL > pctl espandente causale (gate crash)

    uv run python scripts/research/tp01_dvol_overlay.py
"""
from __future__ import annotations
import sys
from pathlib import Path

import numpy as np
import pandas as pd

ROOT = Path(__file__).resolve().parents[2]
sys.path.insert(0, str(ROOT))

from src.data.downloader import load_data
from src.strategies.trend_portfolio import (
    resample_1d, simple_returns, realized_vol, tsmom_blend, CANONICAL,
)

RAW = ROOT / "data" / "raw"
SQ = np.sqrt(365.25)
HOLDOUT = pd.Timestamp("2025-01-01", tz="UTC")
TGT_VOL = CANONICAL["target_vol"]; LEV = CANONICAL["leverage"]; FEE = CANONICAL["fee_side"]
HZ = CANONICAL["horizons_days"]; VW = CANONICAL["vol_win_days"]


def _components(asset: str):
    df = resample_1d(load_data(asset, "1h"))
    c = df["close"].values.astype(float)
    idx = pd.to_datetime(df["datetime"])
    if idx.dt.tz is None:
        idx = idx.dt.tz_localize("UTC")
    idx = pd.DatetimeIndex(idx)                            # tz-aware (UTC)
    r = simple_returns(c)
    rv = realized_vol(r, VW, 365.25)                       # 30g realizzata annualizzata (bpd=1)
    direction = np.clip(tsmom_blend(c, HZ), 0, None)       # long-flat
    dv = pd.read_parquet(RAW / f"dvol_{asset.lower()}.parquet")
    dser = pd.Series(dv["close"].values.astype(float) / 100.0,
                     index=pd.to_datetime(dv["timestamp"], unit="ms", utc=True)).sort_index()
    dvol = dser.reindex(idx, method="ffill").values
    return c, r, idx, rv, direction, dvol


def _net_returns(asset: str, mode: str, tvol: float = TGT_VOL) -> pd.Series:
    c, r, idx, rv, direction, dvol = _components(asset)
    derisk = np.ones(len(c))
    if mode == "realized":
        vol = rv
    elif mode == "dvol":
        vol = dvol
    elif mode == "blend":
        vol = 0.5 * rv + 0.5 * dvol
    elif mode == "max":
        vol = np.fmax(rv, dvol)
    elif mode == "derisk":
        vol = rv
        # gate crash causale: DVOL sopra il suo percentile espandente (90%) -> dimezza l'esposizione
        dd = pd.Series(dvol, index=idx)
        rank = dd.expanding(min_periods=60).apply(lambda x: (x[:-1] < x[-1]).mean() if len(x) > 1 else 0.5, raw=True)
        derisk = np.where(rank.values > 0.90, 0.5, 1.0)
    else:
        raise ValueError(mode)
    with np.errstate(divide="ignore", invalid="ignore"):
        scal = np.where((vol > 0) & np.isfinite(vol), tvol / vol, 0.0)
    tgt = np.clip(direction * scal * derisk, -LEV, LEV)
    tgt[~np.isfinite(tgt)] = 0.0
    pos = np.zeros(len(tgt)); pos[1:] = tgt[:-1]           # decisa a close[t-1], tenuta in t
    gross = pos * r
    turn = np.abs(np.diff(pos, prepend=0.0))
    net = np.clip(gross - FEE * turn, -0.99, None); net[0] = 0.0
    return pd.Series(net, index=idx)


def portfolio(mode: str, tvol: float = TGT_VOL) -> pd.Series:
    b = _net_returns("BTC", mode, tvol); e = _net_returns("ETH", mode, tvol)
    J = pd.concat({"B": b, "E": e}, axis=1, join="inner").fillna(0.0)
    return 0.5 * J["B"] + 0.5 * J["E"]


def metrics(daily: pd.Series, lo=None) -> dict:
    if lo is not None:
        daily = daily[daily.index >= lo]
    r = daily.values
    sh = float(np.mean(r) / np.std(r) * SQ) if np.std(r) > 0 else 0.0
    eq = np.cumprod(1.0 + r); pk = np.maximum.accumulate(eq)
    dd = float(np.max((pk - eq) / pk)) if len(eq) else 0.0
    yrs = (daily.index[-1] - daily.index[0]).days / 365.25 if len(daily) > 1 else 1.0
    cagr = eq[-1] ** (1 / yrs) - 1 if yrs > 0 and len(eq) and eq[-1] > 0 else -1.0
    s = pd.Series(eq, index=daily.index); yearly = {}
    for y, g in s.groupby(s.index.year):
        if len(g) > 1:
            yearly[int(y)] = float(g.iloc[-1] / g.iloc[0] - 1)
    return dict(sharpe=sh, dd=dd, cagr=cagr, tot=float(eq[-1] - 1) if len(eq) else 0.0, yearly=yearly)


def main():
    modes = ["realized", "dvol", "blend", "max", "derisk"]
    series = {m: portfolio(m) for m in modes}
    # vero inizio DVOL (dove TUTTE le varianti hanno dati validi) — non il primo indice del prezzo
    dstart = max(pd.read_parquet(RAW / f"dvol_{a.lower()}.parquet")["timestamp"].min() for a in ("BTC", "ETH"))
    dstart = pd.Timestamp(dstart, unit="ms", tz="UTC") + pd.Timedelta(days=VW)   # +warmup vol-win
    series = {m: s[s.index >= dstart] for m, s in series.items()}
    common = series["realized"].index
    base = series["realized"]

    print("=" * 96)
    print(f"  TP01 × DVOL — vol-target con denominatore di vol diverso. Finestra COMUNE "
          f"{common[0].date()} -> {common[-1].date()} ({len(base)}g). Hold-out 2025+.")
    print("=" * 96)
    print(f"  {'variante':<10} {'FULL Sh':>8} {'FULL DD':>8} {'CAGR':>7} | "
          f"{'HOLD Sh':>8} {'HOLD ret':>9} {'HOLD DD':>8} | per-anno PnL")
    for m in modes:
        s = series[m]; f = metrics(s); h = metrics(s, lo=HOLDOUT)
        ys = " ".join(f"{y}:{p*100:+.0f}" for y, p in sorted(f['yearly'].items()))
        tag = "  (baseline)" if m == "realized" else ""
        print(f"  {m:<10} {f['sharpe']:>+8.2f} {f['dd']*100:>7.1f}% {f['cagr']*100:>+6.0f}% | "
              f"{h['sharpe']:>+8.2f} {h['tot']*100:>+8.1f}% {h['dd']*100:>7.1f}% | {ys}{tag}")

    print("\n  DELTA vs baseline (realized) sulla stessa finestra:")
    bf = metrics(base); bh = metrics(base, lo=HOLDOUT)
    for m in modes:
        if m == "realized":
            continue
        f = metrics(series[m]); h = metrics(series[m], lo=HOLDOUT)
        print(f"    {m:<8}: ΔFULL Sh {f['sharpe']-bf['sharpe']:+.2f}  ΔFULL DD {(f['dd']-bf['dd'])*100:+.1f}pp  "
              f"ΔHOLD Sh {h['sharpe']-bh['sharpe']:+.2f}  ΔHOLD ret {(h['tot']-bh['tot'])*100:+.1f}pp")

    print("\n  CONTROLLO DECISIVO — il taglio di DD del DVOL è 'posizioni più piccole' o vero timing?")
    print("  Confronto le varianti DVOL con il realized a target_vol RIDOTTO (stesso de-levering, ma")
    print("  senza DVOL). Se realized-ridotto eguaglia/batte il DVOL a parità di DD → DVOL non aggiunge.")
    for tv in (0.15, 0.13):
        s = portfolio("realized", tvol=tv); s = s[s.index >= dstart]
        f = metrics(s); h = metrics(s, lo=HOLDOUT)
        print(f"    realized @ vol-tgt {tv*100:.0f}%: FULL Sh {f['sharpe']:+.2f}  DD {f['dd']*100:.1f}%  "
              f"CAGR {f['cagr']*100:+.0f}%  |  HOLD Sh {h['sharpe']:+.2f}")
    mx = metrics(series["max"]); mxh = metrics(series["max"], lo=HOLDOUT)
    print(f"    (vs max-DVOL:        FULL Sh {mx['sharpe']:+.2f}  DD {mx['dd']*100:.1f}%  "
          f"CAGR {mx['cagr']*100:+.0f}%  |  HOLD Sh {mxh['sharpe']:+.2f})")


if __name__ == "__main__":
    main()