foxhunt/scripts/surfer/vrp_diversifier.py

#!/usr/bin/env python3
"""Crypto volatility-risk-premium (VRP) as a diversifier to the momentum book.

VRP = implied vol (Deribit DVOL, BTC+ETH) systematically exceeds realized vol; a vol-seller
harvests it. Short-variance daily P&L proxy:  vrp_t = (DVOL_{t-1}/100)^2/365 - r_t^2
(collect implied daily variance, pay realized squared return). Approximate (no Greeks) but
captures the premium's edge + correlation — enough for a diversifier CHECK.

Then the prize: Sharpe of VRP, per-year (incl. crash years), worst day, CORRELATION to the
crypto momentum book, and the combined two-premium book. A real diversifier needs Sharpe>=~0.4
AND low correlation. Honest: short-vol has negative skew -> Sharpe FLATTERS it; watch the tails.
"""
import json
import math
import os
import sys
import urllib.request

import numpy as np

sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
import pit_sweep  # noqa: E402
from surfer_poc import compute_weights, CFG  # noqa: E402
from signal_sweep import sharpe_t  # noqa: E402
import torch  # noqa: E402

DEV = "cuda" if torch.cuda.is_available() else "cpu"
DAY_MS = 86_400_000
DVOL_DIR = "data/surfer/dvol"


def _get(u):
    return json.load(urllib.request.urlopen(urllib.request.Request(u, headers={"User-Agent": "curl/8"}), timeout=30))


def dvol(ccy):
    os.makedirs(DVOL_DIR, exist_ok=True)
    cache = f"{DVOL_DIR}/{ccy}.json"
    if os.path.exists(cache):
        return {int(k): v for k, v in json.load(open(cache)).items()}
    out, end, start = {}, 1780704000000, 1577836800000
    for _ in range(8):
        u = f"https://www.deribit.com/api/v2/public/get_volatility_index_data?currency={ccy}&start_timestamp={start}&end_timestamp={end}&resolution=86400"
        d = _get(u).get("result", {}).get("data", [])
        if not d:
            break
        for row in d:
            out[int(row[0]) // DAY_MS] = float(row[4])
        end = d[0][0] - 1
        if len(d) < 1000:
            break
    json.dump({str(k): v for k, v in out.items()}, open(cache, "w"))
    return out


def main():
    syms, cdays, cclose, cqv, cfund = pit_sweep.load()
    idx = {s: j for j, s in enumerate(syms)}
    # realized daily returns for BTC/ETH
    def closes(sym):
        j = idx[sym]; c = cclose[:, j]
        return {int(cdays[t]): c[t] for t in range(len(cdays)) if np.isfinite(c[t])}
    btc, eth = closes("BTCUSDT"), closes("ETHUSDT")
    dv_btc, dv_eth = dvol("BTC"), dvol("ETH")
    print(f"DVOL coverage: BTC {len(dv_btc)}d ({min(dv_btc)}..{max(dv_btc)}), ETH {len(dv_eth)}d")

    def vrp_series(px, dv):
        days = sorted(set(px) & set(dv))
        pnl = {}
        for i in range(1, len(days)):
            d0, d1 = days[i - 1], days[i]
            if px[d0] > 0 and px[d1] > 0:
                r = math.log(px[d1] / px[d0])
                iv = dv[d0] / 100.0
                pnl[d1] = (iv * iv) / 365.0 - r * r            # short-variance daily P&L
        return pnl
    vb, ve = vrp_series(btc, dv_btc), vrp_series(eth, dv_eth)
    alld = sorted(set(vb) | set(ve))
    vrp = np.array([np.nanmean([vb.get(d, np.nan), ve.get(d, np.nan)]) for d in alld])
    alld = np.array(alld)
    year = (1970 + alld / 365.25).astype(int)
    sr = sharpe_t(torch.tensor(vrp, device=DEV, dtype=torch.float64))
    print(f"\n===== CRYPTO VRP (short-vol BTC+ETH) =====")
    print(f"days={len(vrp)}  Sharpe {sr:+.2f}  worst-day {np.nanmin(vrp)/np.nanstd(vrp):+.1f}σ  skew {float(((vrp-np.nanmean(vrp))**3).mean()/np.nanstd(vrp)**3):+.2f}")
    print("per-year: " + " ".join(f"{y}:{sharpe_t(torch.tensor(vrp[year==y],device=DEV,dtype=torch.float64)):+.2f}" for y in range(2021, 2027) if (year == y).sum() > 30))

    # crypto momentum book
    cR = np.zeros_like(cclose); cR[1:] = np.log(cclose)[1:] - np.log(cclose)[:-1]; cR = np.where(np.isfinite(cR), cR, 0.0)
    cf = np.where(np.isfinite(cfund), cfund, 0.0)
    cw, _ = compute_weights(cclose, cqv, cdays, CFG)
    mom = np.sum(cw[:-1] * (cR - cf)[1:], axis=1)
    mom_by_day = {int(cdays[1:][i]): mom[i] for i in range(len(mom))}
    common = sorted(set(alld.tolist()) & set(mom_by_day))
    va = np.array([vrp[np.where(alld == d)[0][0]] for d in common])
    ma = np.array([mom_by_day[d] for d in common])
    m = np.isfinite(va) & np.isfinite(ma)
    corr = float(np.corrcoef(va[m], ma[m])[0, 1])
    sv, sm = np.std(va[m]), np.std(ma[m])
    comb = ((1 / sv) * va[m] + (1 / sm) * ma[m]) / (1 / sv + 1 / sm)
    svr = sharpe_t(torch.tensor(va[m], device=DEV, dtype=torch.float64))
    smr = sharpe_t(torch.tensor(ma[m], device=DEV, dtype=torch.float64))
    scr = sharpe_t(torch.tensor(comb, device=DEV, dtype=torch.float64))
    print(f"\n===== DIVERSIFIER CHECK — VRP vs crypto-momentum (overlap {m.sum()} days) =====")
    print(f"VRP Sharpe {svr:+.2f}  momentum Sharpe {smr:+.2f}  CORRELATION {corr:+.2f}  combined {scr:+.2f}")
    opt = math.sqrt(max(svr, 0)**2 + max(smr, 0)**2)
    print(f"optimal-weight combined (uncorrelated approx) ~ {opt:.2f}  (vs momentum alone {smr:+.2f})")
    print("\nVERDICT: VRP Sharpe>=~0.4 + low |corr| + combined>momentum-alone = genuine diversifier (mind the short-vol tail).")


if __name__ == "__main__":
    main()