Files
foxhunt/scripts/surfer/equity_vrp2.py
jgrusewski d43fb5e61f feat(surfer): equity-index VRP frontier — real+affordable, measurement needs engineering
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>
2026-06-06 22:13:34 +02:00

88 lines
4.0 KiB
Python
Raw Permalink Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
#!/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()