PythagorasGoal/Old/scripts/waste/W17_vol_premium.py

"""S2-07: Variance Risk Premium harvesting — versione raffinata.
Ottimizzazione del vol selling con:
1. IV/RV ratio dinamico per entry timing
2. Tail risk cutoff (chiudi se move > N sigma)
3. Position sizing proporzionale al premium
4. Combinazione con directional bias (da gap fade)
5. Multi-asset portfolio (ETH + BTC)
"""
from __future__ import annotations
import sys
sys.path.insert(0, ".")

import numpy as np
import pandas as pd
from scipy.stats import norm
from src.data.downloader import load_data

FEE = 0.001
INITIAL = 1000


def realized_vol(close, window=24):
    log_ret = np.diff(np.log(np.where(close == 0, 1e-10, close)))
    result = np.full(len(close), 0.5)
    for i in range(window, len(log_ret)):
        result[i + 1] = np.std(log_ret[i - window : i]) * np.sqrt(24 * 365)
    return result


def run_vrp(asset):
    print(f"\n{'#'*60}")
    print(f"  {asset} 1h — VARIANCE RISK PREMIUM REFINED")
    print(f"{'#'*60}")

    df = load_data(asset, "1h")
    close = df["close"].values
    high = df["high"].values
    low = df["low"].values
    n = len(close)
    split = int(n * 0.7)
    timestamps = pd.to_datetime(df["timestamp"], unit="ms", utc=True)

    rv_24 = realized_vol(close, 24)
    rv_48 = realized_vol(close, 48)
    rv_168 = realized_vol(close, 168)

    configs = [
        # (dte_h, iv_mult, cutoff_sigma, pos_base, entry_hour, dynamic_sizing, name)
        (24, 1.20, 2.5, 0.10, 8, False, "24h_base"),
        (24, 1.25, 2.5, 0.12, 8, False, "24h_highPrem"),
        (24, 1.20, 2.0, 0.10, 8, False, "24h_tightCut"),
        (24, 1.20, 3.0, 0.12, 8, False, "24h_wideCut"),
        (48, 1.20, 2.5, 0.12, 8, False, "48h_base"),
        (48, 1.25, 2.5, 0.15, 8, False, "48h_highPrem"),
        (48, 1.30, 2.5, 0.15, 8, False, "48h_vhighPrem"),
        (48, 1.20, 3.0, 0.15, 8, False, "48h_wideCut"),
        (24, 1.20, 2.5, 0.10, 8, True, "24h_dynSize"),
        (48, 1.20, 2.5, 0.12, 8, True, "48h_dynSize"),
        (24, 1.20, 2.5, 0.10, 0, False, "24h_midnight"),
        (24, 1.20, 2.5, 0.10, 16, False, "24h_afternoon"),
        (36, 1.22, 2.5, 0.12, 8, False, "36h_medium"),
        (24, 1.15, 2.5, 0.08, 8, False, "24h_lowPrem_safe"),
        (48, 1.20, 2.0, 0.10, 8, True, "48h_tight_dyn"),
    ]

    for dte, iv_mult, cutoff, pos_base, entry_h, dyn_size, name in configs:
        capital = float(INITIAL)
        correct = 0
        total = 0
        daily_trades = {}
        peak_capital = capital
        max_dd = 0

        for i in range(max(split, 170), n - dte):
            day = timestamps.iloc[i].strftime("%Y-%m-%d")
            if daily_trades.get(day, 0) >= 1:
                continue

            if timestamps.iloc[i].hour != entry_h:
                continue

            rv_s = rv_24[i]
            rv_l = rv_168[i]
            if rv_s <= 0.05 or rv_l <= 0.05:
                continue

            iv_est = rv_l * iv_mult
            vrp = iv_est - rv_s

            if vrp <= 0:
                continue

            spot = close[i]
            t = dte / (24 * 365)
            daily_std = rv_s / np.sqrt(365)

            # Premium = IV * sqrt(t) * spot * factor
            premium = iv_est * np.sqrt(t) * spot * 0.4
            premium_pct = premium / spot

            # Expected move based on IV
            expected_move = iv_est * np.sqrt(t) * spot

            # Cutoff: close if actual move > cutoff * expected_move
            max_allowed_move = expected_move * cutoff

            # Dynamic sizing: more when VRP is high
            if dyn_size:
                vrp_ratio = vrp / rv_s
                pos_pct = min(pos_base * (1 + vrp_ratio), pos_base * 2)
            else:
                pos_pct = pos_base

            # Check actual path
            exit_idx = min(i + dte, n - 1)
            actual_move = abs(close[exit_idx] - spot)

            # Early exit: check if intra-period move exceeds cutoff
            breached = False
            for j in range(i + 1, exit_idx + 1):
                intra_move = abs(close[j] - spot)
                if intra_move > max_allowed_move:
                    breached = True
                    exit_idx = j
                    actual_move = intra_move
                    break

            if breached:
                loss = min(actual_move / spot, 0.05) * pos_pct
                pnl = -loss
            else:
                profit = premium_pct * pos_pct
                partial_loss = max(0, actual_move / spot - premium_pct) * pos_pct * 0.5
                pnl = profit - partial_loss

            fee_cost = FEE * 2 * pos_pct
            net = pnl - fee_cost

            capital += capital * net
            capital = max(capital, 0)

            if capital > peak_capital:
                peak_capital = capital
            dd = (peak_capital - capital) / peak_capital if peak_capital > 0 else 0
            max_dd = max(max_dd, dd)

            total += 1
            if pnl > 0:
                correct += 1
            daily_trades[day] = daily_trades.get(day, 0) + 1

        if total < 20:
            continue

        acc = correct / total * 100
        ret = (capital - INITIAL) / INITIAL * 100
        test_days = (n - split) / 24
        test_years = test_days / 365.25
        ann = ((capital / INITIAL) ** (1 / test_years) - 1) * 100 if test_years > 0 and capital > 0 else -100
        dpnl = (capital - INITIAL) / test_days if test_days > 0 else 0
        days_active = len(daily_trades)

        tag = "✅✅" if acc >= 70 and ann >= 50 else "✅" if acc >= 65 and ann >= 30 else ""
        print(f"  {name:22s}: trades={total:4d} acc={acc:.1f}% ret={ret:+.1f}% ann={ann:+.1f}% dd={max_dd*100:.1f}% €/day={dpnl:.2f} active={days_active} {tag}")

    return daily_trades


# Run both assets
results = {}
for asset in ["ETH", "BTC"]:
    results[asset] = run_vrp(asset)

# Multi-asset portfolio simulation
print(f"\n{'#'*60}")
print(f"  MULTI-ASSET PORTFOLIO: ETH + BTC")
print(f"{'#'*60}")

df_eth = load_data("ETH", "1h")
df_btc = load_data("BTC", "1h")
close_eth = df_eth["close"].values
close_btc = df_btc["close"].values
n = min(len(close_eth), len(close_btc))
split = int(n * 0.7)
ts = pd.to_datetime(df_eth["timestamp"].values[:n], unit="ms", utc=True)

rv_eth = realized_vol(close_eth[:n], 168)
rv_btc = realized_vol(close_btc[:n], 168)

capital = float(INITIAL)
total = 0
correct = 0
peak = capital
max_dd = 0
daily_trades = {}

for i in range(max(split, 170), n - 48):
    day = ts[i].strftime("%Y-%m-%d")
    if daily_trades.get(day, 0) >= 1:
        continue
    if ts[i].hour != 8:
        continue

    for asset_close, rv_arr, name in [(close_eth[:n], rv_eth, "ETH"), (close_btc[:n], rv_btc, "BTC")]:
        rv = rv_arr[i]
        if rv <= 0.05:
            continue
        iv = rv * 1.22
        spot = asset_close[i]
        t = 48 / (24 * 365)
        premium_pct = iv * np.sqrt(t) * 0.4
        expected_move = iv * np.sqrt(t) * spot
        max_move = expected_move * 2.5

        exit_idx = min(i + 48, n - 1)
        actual_move = abs(asset_close[exit_idx] - spot)

        breached = False
        for j in range(i + 1, exit_idx + 1):
            if abs(asset_close[j] - spot) > max_move:
                breached = True
                actual_move = abs(asset_close[j] - spot)
                break

        pos_pct = 0.07  # 7% per asset = 14% total
        if breached:
            pnl = -min(actual_move / spot, 0.05) * pos_pct
        else:
            profit = premium_pct * pos_pct
            partial = max(0, actual_move / spot - premium_pct) * pos_pct * 0.5
            pnl = profit - partial

        capital += capital * (pnl - FEE * 2 * pos_pct)
        capital = max(capital, 0)
        total += 1
        if pnl > 0:
            correct += 1

    if capital > peak:
        peak = capital
    dd = (peak - capital) / peak if peak > 0 else 0
    max_dd = max(max_dd, dd)
    daily_trades[day] = daily_trades.get(day, 0) + 1

if total > 0:
    acc = correct / total * 100
    ret = (capital - INITIAL) / INITIAL * 100
    test_days = (n - split) / 24
    test_years = test_days / 365.25
    ann = ((capital / INITIAL) ** (1 / test_years) - 1) * 100 if test_years > 0 and capital > 0 else -100
    dpnl = (capital - INITIAL) / test_days if test_days > 0 else 0
    print(f"\n  ETH+BTC 48h portfolio: trades={total:4d} acc={acc:.1f}% ret={ret:+.1f}% ann={ann:+.1f}% dd={max_dd*100:.1f}% €/day={dpnl:.2f} active={len(daily_trades)}")