spb/wp7_uqo Public

1#!/usr/bin/env python32# =============================================================================3# Author: Simon-Pierre Boucher4# Contact: contact@spboucher.ai5# =============================================================================6"""Step 04 — RQ3: Implied vs realized correlation divergence as a stress7predictor.89Stage 1 (raw-dependent): build the daily implied-correlation index from SPX10and constituent ATM IVs, the 22-day realized correlation from daily returns,11merge with VIX, define stress events, and save the combined series.12Stage 2 (parquet-only): predictive stress regressions and crisis-window13analysis. When the raw stores are absent, stage 1 is skipped and stage 214runs from the shipped ``correlation_divergence.parquet``.1516Inputs : options.duckdb + index_5min.duckdb (stage 1, optional),17         data/processed/realized_vol.parquet (stage 1),18         data/processed/correlation_divergence.parquet (stage 2 fallback)19Outputs: data/processed/correlation_divergence.parquet (stage 1),20         results/rq3_stress_prediction.csv, results/rq3_crisis_analysis.csv21"""2223import warnings2425import numpy as np26import pandas as pd2728import _bootstrap  # noqa: F40129from wp7 import config30from wp7.data_io import (RawDataUnavailableError, load_correlation_divergence,31                         load_realized_vol, open_raw_db)32from wp7.econometrics import add_constant, hc1_tstats, ols, r_squared, standardize3334warnings.filterwarnings('ignore')3536STRESS_HORIZONS = [5, 10, 20]373839# --------------------------------------------------------------------------40# Stage 1 — build the correlation-divergence dataset (needs raw stores)41# --------------------------------------------------------------------------42def build_correlation_dataset() -> pd.DataFrame:43    """Implied correlation (CBOE methodology), realized correlation, VIX."""44    constituents = config.SPX_CONSTITUENTS45    const_str = ",".join([f"'{t}'" for t in constituents])4647    # 3A. SPX and constituent 30-day ATM implied volatilities48    opt_con = open_raw_db("options")49    print("  Extracting SPX implied volatility...")50    spx_iv = opt_con.execute("""51        SELECT trade_date,52               AVG(CASE WHEN bid_iv > 0 AND ask_iv > 0 THEN (bid_iv + ask_iv)/2.053                        WHEN ask_iv > 0 THEN ask_iv ELSE bid_iv END) AS spx_iv_atm54        FROM option_chain55        WHERE ticker = 'SPX'56          AND call_put = 'c'57          AND ABS(delta - 0.5) < 0.158          AND (expiry_date - trade_date) BETWEEN 20 AND 4059          AND ask_price > 060        GROUP BY trade_date61        ORDER BY trade_date62    """).fetchdf()63    print(f"  SPX IV: {len(spx_iv)} days")6465    print("  Extracting constituent implied volatilities...")66    const_iv = opt_con.execute(f"""67        SELECT ticker, trade_date,68               AVG(CASE WHEN bid_iv > 0 AND ask_iv > 0 THEN (bid_iv + ask_iv)/2.069                        WHEN ask_iv > 0 THEN ask_iv ELSE bid_iv END) AS iv_atm70        FROM option_chain71        WHERE ticker IN ({const_str})72          AND call_put = 'c'73          AND ABS(delta - 0.5) < 0.174          AND (expiry_date - trade_date) BETWEEN 20 AND 4075          AND ask_price > 076        GROUP BY ticker, trade_date77        ORDER BY trade_date78    """).fetchdf()79    print(f"  Constituent IVs: {len(const_iv)} rows")80    opt_con.close()8182    # Implied correlation: rho = (sigma²_idx − n⁻¹·avg(sigma²_i)) /83    #                            ((1 − n⁻¹)·avg(sigma_i)²)84    const_iv['iv_sq'] = const_iv['iv_atm'] ** 285    avg_const_iv = const_iv.groupby('trade_date').agg(86        avg_const_iv_sq=('iv_sq', 'mean'),87        avg_const_iv=('iv_atm', 'mean'),88        n_constituents=('ticker', 'nunique'),89    ).reset_index()9091    ic_df = pd.merge(spx_iv, avg_const_iv, on='trade_date', how='inner')92    ic_df['spx_iv_sq'] = ic_df['spx_iv_atm'] ** 293    n = ic_df['n_constituents']94    ic_df['implied_corr'] = (ic_df['spx_iv_sq'] - (1.0 / n) * ic_df['avg_const_iv_sq']) / \95                            ((1.0 - 1.0 / n) * ic_df['avg_const_iv'] ** 2)96    ic_df['implied_corr'] = ic_df['implied_corr'].clip(0, 1)97    print(f"  Implied correlation computed: {len(ic_df)} days")9899    # 3B. Realized correlation: 22-day rolling average pairwise correlation100    #     of constituent daily returns.101    print("\n  Computing realized correlations from daily returns...")102    rv = load_realized_vol()103    rv = rv[rv['ticker'].isin(constituents)]104105    ret_pivot = rv.pivot_table(index='trade_date', columns='ticker',106                               values='daily_return', aggfunc='first')107    ret_pivot = ret_pivot.dropna(axis=0, thresh=10)108109    rolling_corr = []110    dates = ret_pivot.index.tolist()111    for i in range(21, len(dates)):112        window = ret_pivot.iloc[i - 21:i + 1]113        corr_matrix = window.corr()114        mask = np.triu(np.ones(corr_matrix.shape, dtype=bool), k=1)115        rolling_corr.append({116            'trade_date': dates[i],117            'realized_corr': corr_matrix.values[mask].mean(),118            'n_pairs': mask.sum(),119        })120    realized_corr_df = pd.DataFrame(rolling_corr)121    print(f"  Realized correlation computed: {len(realized_corr_df)} days")122123    # 3C. Merge implied and realized correlation124    ic_df['trade_date'] = pd.to_datetime(ic_df['trade_date'])125    corr_merged = pd.merge(ic_df, realized_corr_df, on='trade_date', how='inner')126    corr_merged['corr_divergence'] = corr_merged['implied_corr'] - corr_merged['realized_corr']127    corr_merged['corr_ratio'] = corr_merged['implied_corr'] / \128        corr_merged['realized_corr'].clip(0.01)129130    print(f"  Merged correlation data: {len(corr_merged)} days")131    print(f"  Avg implied corr: {corr_merged['implied_corr'].mean():.4f}")132    print(f"  Avg realized corr: {corr_merged['realized_corr'].mean():.4f}")133    print(f"  Avg divergence: {corr_merged['corr_divergence'].mean():.4f}")134135    # 3D. Stress events from VIX levels and 5-day VIX spikes136    con = open_raw_db("indices_5min")137    vix_daily = con.execute("""138        SELECT CAST(datetime AS DATE) AS trade_date,139               LAST(close) AS vix_close,140               MAX(close) AS vix_high141        FROM ohlcv142        WHERE symbol = 'VIX'143        GROUP BY CAST(datetime AS DATE)144        ORDER BY trade_date145    """).fetchdf()146    con.close()147    vix_daily['trade_date'] = pd.to_datetime(vix_daily['trade_date'])148    corr_merged = pd.merge(corr_merged, vix_daily, on='trade_date', how='left')149150    corr_merged['vix_5d_change'] = corr_merged['vix_close'].pct_change(5)151    corr_merged['stress_event'] = ((corr_merged['vix_close'] > 25) |152                                   (corr_merged['vix_5d_change'] > 0.20)).astype(int)153    for horizon in STRESS_HORIZONS:154        corr_merged[f'stress_fwd_{horizon}d'] = corr_merged['stress_event'].rolling(155            horizon, min_periods=1).max().shift(-horizon)156157    corr_merged.to_parquet(config.CORRELATION_DIVERGENCE_PARQUET, index=False)158    return corr_merged159160161# --------------------------------------------------------------------------162# Stage 2 — stress regressions and crisis windows (parquet-only)163# --------------------------------------------------------------------------164def stress_regressions(corr_merged: pd.DataFrame) -> None:165    """Linear probability models: divergence measures → forward stress."""166    print("\n--- STRESS PREDICTION REGRESSIONS ---")167    features = ['corr_divergence', 'corr_ratio', 'implied_corr', 'realized_corr',168                'spx_iv_atm', 'vix_close']169    rows = []170    for horizon in STRESS_HORIZONS:171        target = f'stress_fwd_{horizon}d'172        sub = corr_merged[features + [target]].dropna()173        if len(sub) < 100:174            continue175176        X = add_constant(standardize(sub[features].values))177        y = sub[target].values178        coefs = ols(X, y)179        r2 = r_squared(y, X @ coefs)180        _, t_stats = hc1_tstats(X, y, coefs, len(features))181182        for i, fname in enumerate(['const'] + features):183            rows.append({184                'horizon': f'{horizon}d',185                'variable': fname,186                'coefficient': coefs[i],187                't_stat': t_stats[i],188                'significant': abs(t_stats[i]) > 1.96,189            })190191        print(f"\n  {horizon}-day stress prediction: R²={r2:.4f}, N={len(sub)}")192        for i, fname in enumerate(['const'] + features):193            sig = "*" if abs(t_stats[i]) > 1.96 else ""194            sig2 = "*" if abs(t_stats[i]) > 2.576 else ""195            print(f"    {fname:20s}: β={coefs[i]:8.4f}, t={t_stats[i]:7.3f} {sig}{sig2}")196197    pd.DataFrame(rows).to_csv(config.RESULTS_DIR / "rq3_stress_prediction.csv", index=False)198199200def crisis_analysis(corr_merged: pd.DataFrame) -> None:201    """Correlation divergence in 30-day pre / 10-day post windows around crises."""202    print("\n--- CORRELATION DIVERGENCE AROUND CRISES ---")203    rows = []204    for name, date_str in config.CRISES.items():205        date = pd.Timestamp(date_str)206        pre = corr_merged[(corr_merged['trade_date'] >= date - pd.Timedelta(days=30)) &207                          (corr_merged['trade_date'] < date)]208        post = corr_merged[(corr_merged['trade_date'] >= date) &209                           (corr_merged['trade_date'] < date + pd.Timedelta(days=10))]210        if len(pre) > 0 and len(post) > 0:211            rows.append({212                'crisis': name,213                'pre_impl_corr': pre['implied_corr'].mean(),214                'pre_real_corr': pre['realized_corr'].mean(),215                'pre_divergence': pre['corr_divergence'].mean(),216                'post_impl_corr': post['implied_corr'].mean(),217                'post_real_corr': post['realized_corr'].mean(),218                'post_divergence': post['corr_divergence'].mean(),219                'divergence_change': post['corr_divergence'].mean() - pre['corr_divergence'].mean(),220                'pre_vix': pre['vix_close'].mean() if 'vix_close' in pre.columns else np.nan,221                'post_vix': post['vix_close'].mean() if 'vix_close' in post.columns else np.nan,222            })223            print(f"  {name}:")224            print(f"    Pre:  IC={pre['implied_corr'].mean():.4f}, "225                  f"RC={pre['realized_corr'].mean():.4f}, "226                  f"Div={pre['corr_divergence'].mean():.4f}")227            print(f"    Post: IC={post['implied_corr'].mean():.4f}, "228                  f"RC={post['realized_corr'].mean():.4f}, "229                  f"Div={post['corr_divergence'].mean():.4f}")230231    pd.DataFrame(rows).to_csv(config.RESULTS_DIR / "rq3_crisis_analysis.csv", index=False)232233234def main():235    print("=" * 70)236    print("RQ3: IMPLIED vs REALIZED CORRELATION DIVERGENCE")237    print("=" * 70)238    config.ensure_output_dirs()239240    try:241        corr_merged = build_correlation_dataset()242    except RawDataUnavailableError as exc:243        print(f"\n[Stage 1 skipped — raw stores unavailable]\n{exc}\n")244        print("Loading shipped correlation_divergence.parquet instead.")245        corr_merged = load_correlation_divergence()246247    stress_regressions(corr_merged)248    crisis_analysis(corr_merged)249    print("\nRQ3 COMPLETE.")250251252if __name__ == "__main__":253    main()254