Bought XSP (mini-SPX) options ohlcv-1d 2013-2026 = $66.86 (hard-capped $70, get_cost-gated, year-chunked vs 504, 2025 backfilled) + SPY $0.01. CORRECTION: index VRP is NOT $900 (that was SPY root); clean instrument XSP $66.86/13y or SPX $143/13y. THREE measurement attempts all gave spurious NEGATIVE VRP (RV 28-33% vs implied ~16%) = MEASUREMENT ERROR not finding (contradicts decades of SPX VRP evidence). Bugs: noisy parity underlying; ~469/2800 days survive -> multi-day gaps inflate RV; SPY!=XSP divergence corrupts ATM. Root: XSP EOD ohlcv too sparse to reconstruct clean underlying+ATM-IV. Implied (16.7%) reads right. Frontier real+affordable+in-hand but proper extraction is a real options-quant build (denser SPX +$76 or IV-surface w/ quotes), not a gate. Did NOT record -2.7 as a result (artifact). Crypto momentum+VRP remains only deploy-grade edge. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
88 lines
4.0 KiB
Python
88 lines
4.0 KiB
Python
#!/usr/bin/env python3
|
||
"""Equity-index VRP, CORRECTED: clean SPY underlying (not noisy parity) + XSP option-implied IV.
|
||
|
||
Fixes the bug where parity-derived underlying inflated realized vol (spurious -VRP). SPY ~= XSP
|
||
(~SPX/10), so SPY is the clean underlying for ATM reference + realized vol. Overlap 2023-2026.
|
||
VRP = implied (XSP ATM straddle) - realized (clean SPY). Short-vol Sharpe, per-year, tail, tail-managed.
|
||
"""
|
||
import math
|
||
import os
|
||
import sys
|
||
|
||
import numpy as np
|
||
|
||
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
|
||
from equity_vrp import load_xsp, to_epoch_day # noqa: E402
|
||
from signal_sweep import sharpe_t # noqa: E402
|
||
import torch # noqa: E402
|
||
|
||
DEV = "cuda" if torch.cuda.is_available() else "cpu"
|
||
DAY_NS = 86_400 * 10**9
|
||
|
||
|
||
def load_spy():
|
||
import databento as db
|
||
a = db.DBNStore.from_file("data/surfer/spy_ohlcv1d.dbn").to_ndarray()
|
||
day = a["ts_event"].astype(np.int64) // DAY_NS
|
||
close = a["close"].astype(np.float64) / 1e9
|
||
return {int(day[i]): float(close[i]) for i in range(len(day)) if close[i] > 0}
|
||
|
||
|
||
def main():
|
||
spy = load_spy()
|
||
print(f"SPY clean: {len(spy)} days ({min(spy)}..{max(spy)}), level {spy[min(spy)]:.0f}->{spy[max(spy)]:.0f}")
|
||
d = load_xsp()
|
||
d["exp_day"] = d["expiry"].map(to_epoch_day)
|
||
d["ttm"] = d["exp_day"] - d["day"]
|
||
d = d[(d["ttm"] >= 20) & (d["ttm"] <= 45)]
|
||
d = d[d["day"].isin(spy.keys())] # overlap with clean SPY
|
||
|
||
iv = {}
|
||
for day, g in d.groupby("day"):
|
||
S = spy[int(day)]
|
||
exp = g.iloc[(g["ttm"] - 30).abs().argsort()].iloc[0]["exp_day"]
|
||
ge = g[g["exp_day"] == exp]
|
||
calls = ge[ge["right"] == "C"].groupby("strike")["close"].mean()
|
||
puts = ge[ge["right"] == "P"].groupby("strike")["close"].mean()
|
||
common = calls.index.intersection(puts.index)
|
||
if len(common) < 3:
|
||
continue
|
||
katm = common[np.abs(np.array(common) - S).argmin()] # ATM = strike nearest clean SPY level
|
||
straddle = float(calls[katm] + puts[katm])
|
||
T = float(ge["ttm"].iloc[0]) / 365.0
|
||
if straddle <= 0 or T <= 0:
|
||
continue
|
||
iv[int(day)] = straddle / (0.8 * S * math.sqrt(T))
|
||
|
||
days = np.array(sorted(set(iv) & set(spy)))
|
||
S = np.array([spy[x] for x in days]); IV = np.array([iv[x] for x in days])
|
||
r = np.zeros(len(days)); r[1:] = np.log(S[1:] / S[:-1])
|
||
H = 21
|
||
rv_fwd = np.full(len(days), np.nan)
|
||
for t in range(len(days) - H):
|
||
rv_fwd[t] = np.std(r[t + 1:t + 1 + H]) * math.sqrt(252)
|
||
vrp = IV - rv_fwd
|
||
iv_lag = np.concatenate([[np.nan], IV[:-1]])
|
||
svol = (iv_lag ** 2) / 252.0 - r ** 2
|
||
year = np.array([1970 + x / 365.25 for x in days]).astype(int)
|
||
T_ = lambda x: torch.tensor(np.asarray(x)[np.isfinite(np.asarray(x))], device=DEV, dtype=torch.float64)
|
||
fin = np.isfinite(vrp)
|
||
|
||
print(f"\n===== EQUITY-INDEX VRP (CORRECTED, clean SPY) — {len(days)} days =====")
|
||
print(f"mean implied {np.nanmean(IV):.1%} mean fwd-realized(clean) {np.nanmean(rv_fwd):.1%} "
|
||
f"mean VRP {np.nanmean(vrp[fin]):+.1%} frac>0 {np.mean(vrp[fin]>0):.2f}")
|
||
sk = float(((svol[np.isfinite(svol)] - np.nanmean(svol)) ** 3).mean() / np.nanstd(svol) ** 3)
|
||
print(f"short-vol daily P&L: Sharpe {sharpe_t(T_(svol)):+.2f} skew {sk:+.2f} worst {np.nanmin(svol)/np.nanstd(svol):+.1f}σ")
|
||
print("per-year: " + " ".join(f"{y}:{sharpe_t(T_(svol[year==y])):+.1f}" for y in range(2023, 2027) if (year == y).sum() > 40))
|
||
rising = np.zeros(len(days))
|
||
for t in range(6, len(days)):
|
||
rising[t] = 0.0 if IV[t - 1] > IV[t - 6] else 1.0
|
||
tm = rising * svol
|
||
print(f"tail-managed (don't sell into rising IV): Sharpe {sharpe_t(T_(tm)):+.2f} worst {np.nanmin(tm)/np.nanstd(tm):+.1f}σ")
|
||
print("per-year TM: " + " ".join(f"{y}:{sharpe_t(T_(tm[year==y])):+.1f}" for y in range(2023, 2027) if (year == y).sum() > 40))
|
||
print("\nVERDICT: mean VRP>0 + frac>0~0.8 confirms premium exists; short-vol Sharpe (esp tail-managed) = harvestable.")
|
||
|
||
|
||
if __name__ == "__main__":
|
||
main()
|