Commit
Stage E exhibits (Agent 9): nine manuscript figures from the v2 results with float fragments — fig_coverage, fig_rq2_loss_ratios, fig_rq2_subperiod, fig_rq3_headline, fig_rq3_series, fig_irf_bands (core); fig_spec_curve, fig_magnet, fig_rq5_clipping (appendix); scripts/33 owns them, scripts/12 and 17 no longer overwrite (WP07-ECON-18, ECON-13, THEO-10, EDIT-15, ADVE-17, cluster H 'no figures')
Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
added
figures/fig_coverage.pdf
+0 −0
Binary file not shown.
added
figures/fig_coverage.png
+0 −0
Binary file not shown.
modified
figures/fig_irf_bands.pdf
+0 −0
Binary file not shown.
modified
figures/fig_irf_bands.png
+0 −0
Binary file not shown.
added
figures/fig_magnet.pdf
+0 −0
Binary file not shown.
added
figures/fig_magnet.png
+0 −0
Binary file not shown.
added
figures/fig_rq2_loss_ratios.pdf
+0 −0
Binary file not shown.
added
figures/fig_rq2_loss_ratios.png
+0 −0
Binary file not shown.
added
figures/fig_rq2_subperiod.pdf
+0 −0
Binary file not shown.
added
figures/fig_rq2_subperiod.png
+0 −0
Binary file not shown.
added
figures/fig_rq3_headline.pdf
+0 −0
Binary file not shown.
added
figures/fig_rq3_headline.png
+0 −0
Binary file not shown.
added
figures/fig_rq3_series.pdf
+0 −0
Binary file not shown.
added
figures/fig_rq3_series.png
+0 −0
Binary file not shown.
added
figures/fig_rq5_clipping.pdf
+0 −0
Binary file not shown.
added
figures/fig_rq5_clipping.png
+0 −0
Binary file not shown.
modified
figures/fig_spec_curve.pdf
+0 −0
Binary file not shown.
added
figures/fig_spec_curve.png
+0 −0
Binary file not shown.
added
paper/figures/fig_coverage.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_coverage from results/descriptive_coverage_complete.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_coverage.pdf} | |
| 5 | +\caption{Coverage of complete-predictor ticker-days by year and asset group, 2010--2025. A ticker-day is complete when the ten option-implied predictors, $RV_t$, $RV^{(w)}_t$ and the next-day return are all observed; the 90-day tenor (needed for the term slope) is the binding constraint. Lines: share of complete ticker-days within each asset group (markers and line styles distinguish the groups in grayscale). Over the full panel 119{,}081 of 264{,}383 ticker-days (45.0\%) are complete; the share rises from 15\% in 2010 to 63\% in 2025. Source: \texttt{results/descriptive\_coverage\_complete.csv} (package L1-E3).} | |
| 6 | +\label{fig:coverage} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_irf_bands.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_irf_bands from results/irf_bands.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_irf_bands.pdf} | |
| 5 | +\caption{Orthogonalized impulse responses of the bivariate VAR(5) in ATM implied volatility and realized variance (standardized levels), averaged across the 68 tickers with at least 1{,}000 complete observations. Lines: cross-ticker mean response to a one-standard-deviation shock under the two recursive orderings (solid, circles: RV ordered first, so that RV shocks move IV contemporaneously but not the reverse; dashed, squares: IV first); shaded bands (plain / hatched): cross-ticker average of the per-ticker residual-bootstrap 5th and 95th percentiles (500 draws, seed 42). Under the RV-first ordering the RV response to an IV shock peaks at 0.150 SD at horizon 1 (band [0.126, 0.173]). Variance decompositions and stationarity diagnostics are in Table~\ref{tab:fevd_v2}. Source: \texttt{results/irf\_bands.csv} (package L1-E8).} | |
| 6 | +\label{fig:irf_bands} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_magnet.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_magnet from results/rq4_price_magnet_v2.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_magnet.pdf} | |
| 5 | +\caption{Price-magnet test on the shipped panel (10{,}480 ticker-days of eight stocks, 2010--2025, open within 0.1\%--10\% of the maximum-open-interest strike among options expiring within 30 days). Circles: observed share of days on which the close is relatively closer to that strike than the open, with $\pm1.96$ date-clustered standard errors; squares: the permutation null that re-pairs each (day, distance, strike side) with the open-to-close move of a random other day of the same ticker (1000 draws, seed 42), mean and 95\% interval; crosses: the mirror-strike placebo. The dotted line is the 50\% benchmark of version~1.2, which is not the correct null because a driftless move reaches a strike at distance $d$ only if it has the right sign and does not overshoot. Over all ticker-days the observed share is 47.01\% against a null of 46.45\% [45.6, 47.3] (permutation $p$ = 0.20; minimum detectable difference 2.0 pp at 80\% power). Source: \texttt{results/rq4\_price\_magnet\_v2.csv} (Table~\ref{tab:magnet_v2}).} | |
| 6 | +\label{fig:magnet} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_rq2_loss_ratios.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_rq2_loss_ratios from results/rq2_loss_ratio_ci.csv, results/rq2_diebold_mariano_v2.csv, results/rq2_oos_summary_v2.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_rq2_loss_ratios.pdf} | |
| 5 | +\caption{Out-of-sample loss of HAR-RV $+$ IV surface relative to HAR-RV, ticker by ticker (67 tickers, rolling 500/250-day windows, training-window winsorization, 2010--2025). Each point is the ratio of the summed out-of-sample loss of the augmented model to that of HAR-RV (below one: the IV surface lowers the loss); vertical segments are 95\% stationary-bootstrap intervals of the ratio (1{,}000 draws, mean block 10 days); marker shape and colour give the sign of the per-ticker Diebold--Mariano statistic (Newey--West variance, $h$ lags). For $RV_{t+1}$ the median MSE ratio is 1.001 (6 tickers significantly better, 19 worse) and the median QLIKE ratio is 1.284 (0 better, 38 worse); the panel DM statistics (pooled loss differential, date-aggregated HAC variance) are -0.94 and -6.49. The Diebold--Mariano statistic is biased toward ``worse'' for the nested augmented model (Table~\ref{tab:rq2_null_mc}); Clark--West statistics are in Table~\ref{tab:rq2_cw}. Sources: \texttt{results/rq2\_loss\_ratio\_ci.csv}, \texttt{rq2\_diebold\_mariano\_v2.csv}, \texttt{rq2\_oos\_summary\_v2.csv}.} | |
| 6 | +\label{fig:rq2_loss_ratios} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_rq2_subperiod.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_rq2_subperiod from results/subperiod_results_v2.csv, results/rq2_subperiod_oos.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_rq2_subperiod.pdf} | |
| 5 | +\caption{Stability of the realized-variance results across the nine subperiods of the sample. Panel (a): adjusted in-sample $R^2$ of HAR-RV and of HAR-RV $+$ IV surface for $RV_{t+1}$, re-estimated within each subperiod on the common sample (per-ticker winsorization); the relative gain ranges from 4.3\% to 114.1\% and, because the models are nested, cannot be negative in sample. Panel (b): median across tickers (at least 30 forecasts in the subperiod, at least five tickers) of the ratio of the out-of-sample loss of the augmented model to that of HAR-RV, from the rolling 500/250 forecasts with training-window winsorization; the 2010--12 subperiod has forecasts for one ticker only and is not shown. For the IV-surface model the median MSE ratio ranges from 0.958 to 1.116 and the median QLIKE ratio from 0.956 to 1.638. Sources: \texttt{results/subperiod\_results\_v2.csv}, \texttt{rq2\_subperiod\_oos.csv} (package L1-B3).} | |
| 6 | +\label{fig:rq2_subperiod} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_rq3_headline.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_rq3_headline from results/rq3_headline.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_rq3_headline.pdf} | |
| 5 | +\caption{Implied--realized correlation ratio as a predictor of market stress: the nine headline coefficients with episode-bootstrap intervals. Linear probability models of the indicator that a post-origin stress event (definition S2: VIX$_{t+j}>25$ or VIX$_{t+j}/$VIX$_t-1>20\%$ for some $j\le h$) occurs within $h$ trading days, on the standardized ratio and the extended controls, rebuilt series 2010--2025 (3{,}668 days at $h=5$). Markers: coefficient in percentage points of probability per standard deviation of the transformed ratio; bars: 95\% interval from resampling the maximal runs of equal target value (stress episodes and calm spells) with replacement, 1{,}000 draws; labels: Holm step-down $p$-values over the nine tests. For the log ratio the coefficients are 1.77, 4.20 and 3.78 pp (Newey--West($h$) $t$ of 1.55, 2.50, 1.73); no interval excludes zero and the smallest Holm $p$-value in the figure is 0.005 (raw ratio, $h=10$; Romano--Wolf $p$ = 0.12). Source: \texttt{results/rq3\_headline.csv} (Table~\ref{tab:rq3_headline}).} | |
| 6 | +\label{fig:rq3_headline} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_rq3_series.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_rq3_series from results/correlation_divergence_v2.parquet — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_rq3_series.pdf} | |
| 5 | +\caption{Implied and realized correlation among S\&P~500 constituents, 2010--2025 (3{,}678 trading days, series rebuilt from the shipped parquets by stage~1b of \texttt{scripts/04}, 29 constituents, days with fewer than 20 constituents dropped). Panel (a): the CBOE-style implied correlation from the ATM implied variances of SPX and of the constituents (eq.~(\ref{eq:implcorr}), clipped to $[0,1]$) and the 22-day pairwise realized correlation; panel (b): the log ratio that enters the stress regressions (mean $-$0.252, standard deviation 0.444); panel (c): the VIX close from the complete daily series of the FinData OS lake. Shaded: days on which the VIX closes above 25 (13.3\% of days), the level component of the post-origin stress definition. Source: \texttt{data/processed/correlation\_divergence\_v2.parquet} (regenerable; diagnostics in Table~\ref{tab:rq3_diag}).} | |
| 6 | +\label{fig:rq3_series} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_rq5_clipping.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_rq5_clipping from results/rq5_clipping_variant.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_rq5_clipping.pdf} | |
| 5 | +\caption{Dependence of the RQ5 benchmark comparison on the treatment of the test-year regressors (SPX, next-day realized variance, rolling annual re-estimation 2016--2025, FinData OS lake data, \texttt{scripts/32}). Lines: out-of-sample $R^2$ of the HAR-RV model by test year (vertical axis truncated below the lowest HAR-RV value; off-scale values of the VIX benchmark are marked with their value) when the test-year regressors are clipped at the training window's 1\%/99\% quantiles (the scheme of the RQ2 design), left unclipped, or when the model is estimated in logs; crosses: the VIX benchmark. Clipping caps $RV_t$ in 2020 far below its realized level: the HAR-RV $R^2_{OOS}$ for 2020 is 0.058 clipped against 0.449 unclipped, and over 2016--2025 it is 0.158 clipped against 0.465 unclipped, versus 0.347 for the VIX benchmark (Diebold--Mariano $t$ of the unclipped HAR-RV against the VIX: 1.28). Source: \texttt{results/rq5\_clipping\_variant.csv} (Table~\ref{tab:rq5_v2}, Panel C).} | |
| 6 | +\label{fig:rq5_clipping} | |
| 7 | +\end{figure} | |
added
paper/figures/fig_spec_curve.tex
+7 −0
@@ -0,0 +1,7 @@ | ||
| 1 | +% Generated by scripts/33_make_exhibit_figures.py::fig_spec_curve from results/rq1_spec_curve.csv — do not edit; regenerate with `make tables` (figure). | |
| 2 | +\begin{figure}[tbp] | |
| 3 | +\centering | |
| 4 | +\includegraphics[width=\textwidth]{fig_spec_curve.pdf} | |
| 5 | +\caption{Specification curve of the 25-delta skew coefficient in the return regressions (400 specifications: winsorization cutoff $\times$ standard-error estimator $\times$ sample $\times$ horizon $\times$ controls $\times$ return basis; Table~\ref{tab:spec_curve}). Panel (a): $t$-statistics sorted in increasing order, marker and colour by standard-error estimator; dashed lines at $\pm1.96$, dotted at $\pm3.02$. Panel (b): the specification behind each point. The median $t$ is 1.93; $|t|>1.96$ in 50.0\% and $|t|>3.02$ in 18.5\% of specifications; among the Driscoll--Kraay specifications the share with $|t|>3.02$ is 5.0\%. Source: \texttt{results/rq1\_spec\_curve.csv} (package L1-C5).} | |
| 6 | +\label{fig:spec_curve} | |
| 7 | +\end{figure} | |
modified
scripts/12_make_figures.py
+2 −0
@@ -118,6 +118,8 @@ def fig_irf_bands() -> None: | ||
| 118 | 118 | to an RV shock under both recursive orderings (levels specification), from |
| 119 | 119 | results/irf_bands.csv written by scripts/11. New output fig_irf_bands; the |
| 120 | 120 | legacy fig_irf (reduced-form IRFs, 20 tickers) is unchanged.""" |
| 121 | + print(" [fig_irf_bands: the manuscript version is drawn by scripts/33_make_exhibit_figures.py (Stage E design system); skipped here]") | |
| 122 | + return | |
| 121 | 123 | path = config.RESULTS_DIR / "irf_bands.csv" |
| 122 | 124 | if not path.exists(): |
| 123 | 125 | print(" [fig_irf_bands skipped: results/irf_bands.csv missing — run scripts/11]") |
modified
scripts/17_specification_curve.py
+6 −1
@@ -131,7 +131,12 @@ def main(): | ||
| 131 | 131 | share_196=('abs_t_gt_196', 'mean'), share_302=('abs_t_gt_302', 'mean')) |
| 132 | 132 | print(f"\n by {key}:\n" + g.round(3).to_string()) |
| 133 | 133 | |
| 134 | − # Figure | |
| 134 | + # Figure: the manuscript version (figures/fig_spec_curve.pdf + paper/figures/fig_spec_curve.tex) is drawn by | |
| 135 | + # scripts/33_make_exhibit_figures.py from results/rq1_spec_curve.csv (Stage E design system); the legacy drawing is kept below but not written. | |
| 136 | + if True: | |
| 137 | + print("\n [fig_spec_curve: drawn by scripts/33_make_exhibit_figures.py]") | |
| 138 | + print("\nSPECIFICATION CURVE COMPLETE.") | |
| 139 | + return | |
| 135 | 140 | import matplotlib |
| 136 | 141 | matplotlib.use('Agg') |
| 137 | 142 | import matplotlib.pyplot as plt |
added
scripts/33_make_exhibit_figures.py
+442 −0
@@ -0,0 +1,442 @@ | ||
| 1 | +#!/usr/bin/env python3 | |
| 2 | +# ============================================================================= | |
| 3 | +# Author: Simon-Pierre Boucher | |
| 4 | +# Contact: contact@spboucher.ai | |
| 5 | +# ============================================================================= | |
| 6 | +"""Step 33 — Manuscript figures (Stage E, 2026-11-04; Part X.2 standard). | |
| 7 | + | |
| 8 | +Every figure is drawn from a validated result file of the Stage D revision | |
| 9 | +(results/*.csv, or the regenerable processed parquets for the RQ3 series), | |
| 10 | +in one design system (src/wp7/figstyle.py: Okabe-Ito palette, markers and | |
| 11 | +line styles for grayscale, vector PDF at \\textwidth, embedded serif face), | |
| 12 | +and comes with a float fragment paper/figures/<name>.tex whose caption | |
| 13 | +quotes only numbers read from the same files. | |
| 14 | + | |
| 15 | +Core figures | |
| 16 | + fig_coverage share of complete ticker-days by year and asset group | |
| 17 | + fig_rq2_loss_ratios per-ticker OOS loss ratios with bootstrap intervals | |
| 18 | + fig_rq2_subperiod in-sample adj. R^2 and OOS loss ratios by subperiod | |
| 19 | + fig_rq3_headline RQ3 coefficients with episode-bootstrap intervals | |
| 20 | + fig_rq3_series implied vs realized correlation, log ratio, VIX | |
| 21 | + fig_irf_bands orthogonalized IRFs with residual-bootstrap bands | |
| 22 | +Appendix figures | |
| 23 | + fig_spec_curve specification curve of the 25-delta skew | |
| 24 | + fig_magnet price-magnet shares against the permutation null | |
| 25 | + fig_rq5_clipping RQ5 R^2_OOS by year under the three clipping rules | |
| 26 | + | |
| 27 | +Inputs : results/descriptive_coverage_complete.csv, rq2_loss_ratio_ci.csv, | |
| 28 | + rq2_diebold_mariano_v2.csv, rq2_oos_summary_v2.csv, subperiod_results_v2.csv, | |
| 29 | + rq2_subperiod_oos.csv, rq3_headline.csv, irf_bands.csv, rq1_spec_curve.csv, | |
| 30 | + rq4_price_magnet_v2.csv, rq5_clipping_variant.csv, | |
| 31 | + data/processed/correlation_divergence_v2.parquet (regenerable, scripts/04) | |
| 32 | +Outputs: figures/<name>.pdf (+ .png), paper/figures/<name>.tex | |
| 33 | +""" | |
| 34 | + | |
| 35 | +import warnings | |
| 36 | + | |
| 37 | +import numpy as np | |
| 38 | +import pandas as pd | |
| 39 | +from matplotlib.lines import Line2D | |
| 40 | +from matplotlib.patches import Patch | |
| 41 | + | |
| 42 | +import _bootstrap # noqa: F401 | |
| 43 | +from wp7 import config | |
| 44 | +from wp7 import figstyle as fs | |
| 45 | +from wp7.tables import SUBPERIOD_LABEL, f, figure, num | |
| 46 | + | |
| 47 | +warnings.filterwarnings('ignore') | |
| 48 | +R = config.RESULTS_DIR | |
| 49 | +P = config.DATA_PROCESSED | |
| 50 | + | |
| 51 | + | |
| 52 | +def _txt(x, d=2): | |
| 53 | + """Plain-text number for captions (LaTeX minus).""" | |
| 54 | + return f(x, d) | |
| 55 | + | |
| 56 | + | |
| 57 | +# ---------------------------------------------------------------------------- | |
| 58 | +# core: data coverage | |
| 59 | +# ---------------------------------------------------------------------------- | |
| 60 | +def fig_coverage(): | |
| 61 | + c = pd.read_csv(R / "descriptive_coverage_complete.csv") | |
| 62 | + g = c.groupby(['group', 'year']).agg(n=('n_rows', 'sum'), k1=('n_complete_rq1', 'sum'), k2=('n_complete_rq2', 'sum')).reset_index() | |
| 63 | + a = c.groupby('year').agg(n=('n_rows', 'sum'), k1=('n_complete_rq1', 'sum'), k2=('n_complete_rq2', 'sum')).reset_index() | |
| 64 | + fig, ax = fs.new(width=6.5, height=2.9) | |
| 65 | + for i, grp in enumerate(['Stocks', 'ETFs', 'Indices']): | |
| 66 | + s = g[g['group'] == grp] | |
| 67 | + ax.plot(s['year'], s['k1'] / s['n'] * 100, color=fs.SERIES[i], marker=fs.MARKERS[i], ls=fs.LINESTYLES[i], label=f"{grp} (n = {c[c['group']==grp]['ticker'].nunique()})") | |
| 68 | + ax.plot(a['year'], a['k1'] / a['n'] * 100, color=fs.INK, marker='x', ls='-', lw=1.0, label='All 69 underlyings') | |
| 69 | + ax.set_ylim(0, 100) | |
| 70 | + ax.set_xlabel('Year') | |
| 71 | + ax.set_ylabel('Complete ticker-days (% of ticker-days)') | |
| 72 | + ax.set_xticks(sorted(a['year'].unique())) | |
| 73 | + ax.tick_params(axis='x', labelrotation=45) | |
| 74 | + ax.legend(loc='upper left', ncol=2) | |
| 75 | + fs.save(fig, 'fig_coverage') | |
| 76 | + tot_n, tot_k = int(a['n'].sum()), int(a['k1'].sum()) | |
| 77 | + first, last = a.iloc[0], a.iloc[-1] | |
| 78 | + cap = (rf"Coverage of complete-predictor ticker-days by year and asset group, 2010--2025. A ticker-day is complete when the ten option-implied " | |
| 79 | + rf"predictors, $RV_t$, $RV^{{(w)}}_t$ and the next-day return are all observed; the 90-day tenor (needed for the term slope) is the binding " | |
| 80 | + rf"constraint. Lines: share of complete ticker-days within each asset group (markers and line styles distinguish the groups in grayscale). " | |
| 81 | + rf"Over the full panel {num(tot_k)} of {num(tot_n)} ticker-days ({_txt(tot_k/tot_n*100,1)}\%) are complete; the share rises from " | |
| 82 | + rf"{_txt(first['k1']/first['n']*100,0)}\% in {int(first['year'])} to {_txt(last['k1']/last['n']*100,0)}\% in {int(last['year'])}. " | |
| 83 | + rf"Source: \texttt{{results/descriptive\_coverage\_complete.csv}} (package L1-E3).") | |
| 84 | + figure('fig_coverage', cap, 'fig:coverage') | |
| 85 | + | |
| 86 | + | |
| 87 | +# ---------------------------------------------------------------------------- | |
| 88 | +# core: RQ2 per-ticker loss ratios | |
| 89 | +# ---------------------------------------------------------------------------- | |
| 90 | +def fig_rq2_loss_ratios(): | |
| 91 | + lr = pd.read_csv(R / "rq2_loss_ratio_ci.csv") | |
| 92 | + dm = pd.read_csv(R / "rq2_diebold_mariano_v2.csv") | |
| 93 | + summ = pd.read_csv(R / "rq2_oos_summary_v2.csv") | |
| 94 | + summ = summ[(summ['scheme'] == 'training-window winsorization') & (summ['subset'] == 'all')] | |
| 95 | + fig, axes = fs.new(width=6.5, height=4.6, nrows=2, ncols=2, sharex=False) | |
| 96 | + panels = [('1-Day RV', 'mse', '(a) MSE, target $RV_{t+1}$'), ('1-Day RV', 'qlike', '(b) QLIKE, target $RV_{t+1}$'), | |
| 97 | + ('5-Day RV', 'mse', r'(c) MSE, target $\sum_{j=1}^{5} RV_{t+j}$'), ('5-Day RV', 'qlike', r'(d) QLIKE, target $\sum_{j=1}^{5} RV_{t+j}$')] | |
| 98 | + col_better, col_worse, col_ns = fs.OI['blue'], fs.OI['vermillion'], fs.OI['grey'] | |
| 99 | + for ax, (tgt, loss, title) in zip(axes.ravel(), panels): | |
| 100 | + s = lr[lr['target'] == tgt].merge(dm[dm['target'] == tgt][['ticker', 'dm_mse', 'dm_qlike']], on='ticker') | |
| 101 | + s = s.sort_values(f'{loss}_ratio').reset_index(drop=True) | |
| 102 | + x = np.arange(len(s)) | |
| 103 | + d = s[f'dm_{loss}'] | |
| 104 | + colors = np.where(d > 1.96, col_better, np.where(d < -1.96, col_worse, col_ns)) | |
| 105 | + markers = np.where(d > 1.96, 'o', np.where(d < -1.96, 's', 'o')) | |
| 106 | + ax.vlines(x, s[f'{loss}_ratio_p025'], s[f'{loss}_ratio_p975'], color=colors, lw=0.8, alpha=0.8) | |
| 107 | + for mk in ('o', 's'): | |
| 108 | + sel = markers == mk | |
| 109 | + ax.scatter(x[sel], s.loc[sel, f'{loss}_ratio'], s=12, marker=mk, c=colors[sel], edgecolor=fs.INK, linewidth=0.3, zorder=3) | |
| 110 | + ax.axhline(1.0, color=fs.INK, lw=0.8, ls='--') | |
| 111 | + fs.panel_label(ax, title) | |
| 112 | + ax.set_xlabel('Tickers ordered by loss ratio') | |
| 113 | + ax.set_ylabel('Loss ratio, augmented / HAR-RV') | |
| 114 | + ax.set_xticks([]) | |
| 115 | + q = summ[(summ['target'] == tgt) & (summ['loss'] == loss)].iloc[0] | |
| 116 | + ax.text(0.02, 0.97, f"median {q['median_ratio']:.3f}; DM $>1.96$: {int(q['n_dm_gt_1.96'])}, $<-1.96$: {int(q['n_dm_lt_-1.96'])} of {int(q['n_tickers'])}\n" | |
| 117 | + f"Holm 5%: {int(q['holm_pos'])} / {int(q['holm_neg'])}; panel DM $t$ = {q['panel_dm_t']:.2f}", | |
| 118 | + transform=ax.transAxes, va='top', ha='left', fontsize=7.5) | |
| 119 | + handles = [Line2D([], [], marker='o', color=col_better, ls='', markeredgecolor=fs.INK, markeredgewidth=0.3, label='DM $> 1.96$ (augmented model better)'), | |
| 120 | + Line2D([], [], marker='s', color=col_worse, ls='', markeredgecolor=fs.INK, markeredgewidth=0.3, label='DM $< -1.96$ (augmented model worse)'), | |
| 121 | + Line2D([], [], marker='o', color=col_ns, ls='', markeredgecolor=fs.INK, markeredgewidth=0.3, label='$|$DM$| \\leq 1.96$'), | |
| 122 | + Line2D([], [], color=fs.INK2, lw=0.8, label='95% stationary-bootstrap interval')] | |
| 123 | + fig.legend(handles=handles, loc='outside lower center', ncol=2) | |
| 124 | + fs.save(fig, 'fig_rq2_loss_ratios') | |
| 125 | + q1m = summ[(summ['target'] == '1-Day RV') & (summ['loss'] == 'mse')].iloc[0] | |
| 126 | + q1q = summ[(summ['target'] == '1-Day RV') & (summ['loss'] == 'qlike')].iloc[0] | |
| 127 | + n_t = int(q1m['n_tickers']) | |
| 128 | + cap = (rf"Out-of-sample loss of HAR-RV $+$ IV surface relative to HAR-RV, ticker by ticker ({n_t} tickers, rolling 500/250-day windows, " | |
| 129 | + rf"training-window winsorization, 2010--2025). Each point is the ratio of the summed out-of-sample loss of the augmented model to that of HAR-RV " | |
| 130 | + rf"(below one: the IV surface lowers the loss); vertical segments are 95\% stationary-bootstrap intervals of the ratio (1{{,}}000 draws, mean block 10 days); " | |
| 131 | + rf"marker shape and colour give the sign of the per-ticker Diebold--Mariano statistic (Newey--West variance, $h$ lags). " | |
| 132 | + rf"For $RV_{{t+1}}$ the median MSE ratio is {q1m['median_ratio']:.3f} ({int(q1m['n_dm_gt_1.96'])} tickers significantly better, {int(q1m['n_dm_lt_-1.96'])} worse) " | |
| 133 | + rf"and the median QLIKE ratio is {q1q['median_ratio']:.3f} ({int(q1q['n_dm_gt_1.96'])} better, {int(q1q['n_dm_lt_-1.96'])} worse); the panel DM statistics " | |
| 134 | + rf"(pooled loss differential, date-aggregated HAC variance) are {q1m['panel_dm_t']:.2f} and {q1q['panel_dm_t']:.2f}. The Diebold--Mariano statistic " | |
| 135 | + rf"is biased toward ``worse'' for the nested augmented model (Table~\ref{{tab:rq2_null_mc}}); Clark--West statistics are in Table~\ref{{tab:rq2_cw}}. " | |
| 136 | + rf"Sources: \texttt{{results/rq2\_loss\_ratio\_ci.csv}}, \texttt{{rq2\_diebold\_mariano\_v2.csv}}, \texttt{{rq2\_oos\_summary\_v2.csv}}.") | |
| 137 | + figure('fig_rq2_loss_ratios', cap, 'fig:rq2_loss_ratios') | |
| 138 | + | |
| 139 | + | |
| 140 | +# ---------------------------------------------------------------------------- | |
| 141 | +# core: RQ2 by subperiod | |
| 142 | +# ---------------------------------------------------------------------------- | |
| 143 | +def fig_rq2_subperiod(): | |
| 144 | + d = pd.read_csv(R / "subperiod_results_v2.csv") | |
| 145 | + o = pd.read_csv(R / "rq2_subperiod_oos.csv") | |
| 146 | + o = o[(o['target'] == '1-Day RV')] | |
| 147 | + sps = list(config.SUBPERIODS) | |
| 148 | + labels = [SUBPERIOD_LABEL[s] for s in sps] | |
| 149 | + fig, axes = fs.new(width=6.5, height=3.1, nrows=1, ncols=2) | |
| 150 | + ax = axes[0] | |
| 151 | + har = [d[(d['subperiod'] == s) & (d['model'] == 'HAR-RV')]['adj_r2_model'].iloc[0] for s in sps] | |
| 152 | + hiv = [d[(d['subperiod'] == s) & (d['model'] == 'HAR+IV')]['adj_r2_model'].iloc[0] for s in sps] | |
| 153 | + x = np.arange(len(sps)) | |
| 154 | + ax.plot(x, har, color=fs.SERIES[0], marker=fs.MARKERS[0], ls=fs.LINESTYLES[0], label='HAR-RV') | |
| 155 | + ax.plot(x, hiv, color=fs.SERIES[1], marker=fs.MARKERS[1], ls=fs.LINESTYLES[1], label='HAR-RV $+$ IV surface') | |
| 156 | + ax.set_xticks(x); ax.set_xticklabels(labels, rotation=45, ha='right') | |
| 157 | + ax.set_ylabel('Adjusted $R^2$, in sample ($RV_{t+1}$)') | |
| 158 | + ax.set_ylim(0, max(hiv) * 1.15) | |
| 159 | + ax.legend(loc='upper left') | |
| 160 | + fs.panel_label(ax, '(a) In-sample fit by subperiod') | |
| 161 | + ax = axes[1] | |
| 162 | + for i, (comp, lab) in enumerate([('HAR-RV + IV surface / HAR-RV', 'HAR-RV $+$ IV surface'), ('HAR-RV + IV level / HAR-RV', r'HAR-RV $+$ $IV_{ATM,30d}$')]): | |
| 163 | + for j, (loss, ll) in enumerate([('mse', 'MSE'), ('qlike', 'QLIKE')]): | |
| 164 | + ys = [] | |
| 165 | + for s in sps: | |
| 166 | + r = o[(o['subperiod'] == s) & (o['comparison'] == comp)] | |
| 167 | + ys.append(r[f'{loss}_ratio_median'].iloc[0] if len(r) and r['n_tickers'].iloc[0] >= 5 else np.nan) | |
| 168 | + k = i * 2 + j | |
| 169 | + ax.plot(x, ys, color=fs.SERIES[k], marker=fs.MARKERS[k], ls=fs.LINESTYLES[k], label=f'{lab}, {ll}') | |
| 170 | + ax.axhline(1.0, color=fs.INK, lw=0.8, ls='--') | |
| 171 | + ax.set_xticks(x); ax.set_xticklabels(labels, rotation=45, ha='right') | |
| 172 | + ax.set_ylabel('Median OOS loss ratio vs HAR-RV ($RV_{t+1}$)') | |
| 173 | + ax.legend(loc='upper left', fontsize=7) | |
| 174 | + fs.panel_label(ax, '(b) Out-of-sample loss ratios by subperiod') | |
| 175 | + fs.save(fig, 'fig_rq2_subperiod') | |
| 176 | + srf = o[o['comparison'] == 'HAR-RV + IV surface / HAR-RV'] | |
| 177 | + srf = srf[(srf['n_tickers'] >= 5) & (srf['subperiod'] != 'Full sample')] | |
| 178 | + dr = [(h2 - h1) / h1 * 100 for h1, h2 in zip(har, hiv)] | |
| 179 | + cap = (rf"Stability of the realized-variance results across the nine subperiods of the sample. Panel (a): adjusted in-sample $R^2$ of HAR-RV and of " | |
| 180 | + rf"HAR-RV $+$ IV surface for $RV_{{t+1}}$, re-estimated within each subperiod on the common sample (per-ticker winsorization); the relative gain ranges from " | |
| 181 | + rf"{_txt(min(dr),1)}\% to {_txt(max(dr),1)}\% and, because the models are nested, cannot be negative in sample. Panel (b): median across tickers " | |
| 182 | + rf"(at least 30 forecasts in the subperiod, at least five tickers) of the ratio of the out-of-sample loss of the augmented model to that of HAR-RV, from the rolling " | |
| 183 | + rf"500/250 forecasts with training-window winsorization; the 2010--12 subperiod has forecasts for one ticker only and is not shown. For the IV-surface model " | |
| 184 | + rf"the median MSE ratio ranges from {srf['mse_ratio_median'].min():.3f} to {srf['mse_ratio_median'].max():.3f} and the median QLIKE ratio from " | |
| 185 | + rf"{srf['qlike_ratio_median'].min():.3f} to {srf['qlike_ratio_median'].max():.3f}. " | |
| 186 | + rf"Sources: \texttt{{results/subperiod\_results\_v2.csv}}, \texttt{{rq2\_subperiod\_oos.csv}} (package L1-B3).") | |
| 187 | + figure('fig_rq2_subperiod', cap, 'fig:rq2_subperiod') | |
| 188 | + | |
| 189 | + | |
| 190 | +# ---------------------------------------------------------------------------- | |
| 191 | +# core: RQ3 headline coefficients | |
| 192 | +# ---------------------------------------------------------------------------- | |
| 193 | +def fig_rq3_headline(): | |
| 194 | + h = pd.read_csv(R / "rq3_headline.csv") | |
| 195 | + trs = [('raw', 'Raw ratio'), ('log', 'Log ratio'), ('rank', 'Rank')] | |
| 196 | + hz = ['5d', '10d', '20d'] | |
| 197 | + fig, ax = fs.new(width=6.5, height=3.0) | |
| 198 | + for i, (tr, lab) in enumerate(trs): | |
| 199 | + xs, ys, lo, hi, ph = [], [], [], [], [] | |
| 200 | + for j, z in enumerate(hz): | |
| 201 | + r = h[(h['transform'] == tr) & (h['horizon'] == z)].iloc[0] | |
| 202 | + xs.append(j + (i - 1) * 0.24); ys.append(r['coefficient_pp']); lo.append(r['boot_ci_lo_pp']); hi.append(r['boot_ci_hi_pp']); ph.append(r['p_holm']) | |
| 203 | + xs, ys, lo, hi = map(np.array, (xs, ys, lo, hi)) | |
| 204 | + ax.errorbar(xs, ys, yerr=[ys - lo, hi - ys], fmt=fs.MARKERS[i], color=fs.SERIES[i], capsize=3, lw=1.1, label=lab, markeredgecolor=fs.INK, markeredgewidth=0.3) | |
| 205 | + for xx, ll, p in zip(xs, lo, ph): | |
| 206 | + ax.annotate(f"Holm\n$p$ = {p:.3f}", xy=(xx, ll), xytext=(0, -4), textcoords='offset points', ha='center', va='top', fontsize=6.8, color=fs.INK2) | |
| 207 | + ax.axhline(0, color=fs.INK, lw=0.8) | |
| 208 | + ax.set_xticks(range(3)); ax.set_xticklabels([f'$h = {z[:-1]}$ days' for z in hz]) | |
| 209 | + ax.set_ylabel('Stress probability, pp per SD of the ratio') | |
| 210 | + ax.set_xlabel('Forecast horizon') | |
| 211 | + fig.legend(loc='outside lower center', ncol=3, title='Transform of $\\hat\\rho_{impl}/\\hat\\rho_{real}$', title_fontsize=8) | |
| 212 | + ax.margins(x=0.12, y=0.18) | |
| 213 | + fs.save(fig, 'fig_rq3_headline') | |
| 214 | + n = int(h['n_obs'].max()); lg = h[h['transform'] == 'log'].set_index('horizon') | |
| 215 | + cap = (rf"Implied--realized correlation ratio as a predictor of market stress: the nine headline coefficients with episode-bootstrap intervals. " | |
| 216 | + rf"Linear probability models of the indicator that a post-origin stress event (definition S2: VIX$_{{t+j}}>25$ or VIX$_{{t+j}}/$VIX$_t-1>20\%$ for some $j\le h$) " | |
| 217 | + rf"occurs within $h$ trading days, on the standardized ratio and the extended controls, rebuilt series 2010--2025 ({num(n)} days at $h=5$). Markers: coefficient in " | |
| 218 | + rf"percentage points of probability per standard deviation of the transformed ratio; bars: 95\% interval from resampling the maximal runs of equal target value " | |
| 219 | + rf"(stress episodes and calm spells) with replacement, 1{{,}}000 draws; labels: Holm step-down $p$-values over the nine tests. For the log ratio the coefficients are " | |
| 220 | + rf"{_txt(lg.loc['5d','coefficient_pp'])}, {_txt(lg.loc['10d','coefficient_pp'])} and {_txt(lg.loc['20d','coefficient_pp'])} pp " | |
| 221 | + rf"(Newey--West($h$) $t$ of {_txt(lg.loc['5d','t_nw'])}, {_txt(lg.loc['10d','t_nw'])}, {_txt(lg.loc['20d','t_nw'])}); no interval excludes zero and the smallest " | |
| 222 | + rf"Holm $p$-value in the figure is {h['p_holm'].min():.3f} (raw ratio, $h=10$; Romano--Wolf $p$ = {h.loc[h['p_holm'].idxmin(),'p_romano_wolf']:.2f}). " | |
| 223 | + rf"Source: \texttt{{results/rq3\_headline.csv}} (Table~\ref{{tab:rq3_headline}}).") | |
| 224 | + figure('fig_rq3_headline', cap, 'fig:rq3_headline') | |
| 225 | + | |
| 226 | + | |
| 227 | +# ---------------------------------------------------------------------------- | |
| 228 | +# core: RQ3 series | |
| 229 | +# ---------------------------------------------------------------------------- | |
| 230 | +def fig_rq3_series(): | |
| 231 | + path = P / "correlation_divergence_v2.parquet" | |
| 232 | + if not path.exists(): | |
| 233 | + print(" [fig_rq3_series skipped: data/processed/correlation_divergence_v2.parquet missing — run scripts/04]") | |
| 234 | + return | |
| 235 | + d = pd.read_parquet(path).sort_values('trade_date') | |
| 236 | + d['log_ratio'] = np.log(d['corr_ratio'].where(d['corr_ratio'] > 0)) | |
| 237 | + fig, axes = fs.new(width=6.5, height=5.6, nrows=3, ncols=1, sharex=True) | |
| 238 | + hi_vix = d['vix_close'] > 25 | |
| 239 | + for ax in axes: | |
| 240 | + ax.fill_between(d['trade_date'], 0, 1, where=hi_vix, transform=ax.get_xaxis_transform(), color=fs.OI['grey'], alpha=0.25, lw=0) | |
| 241 | + ax = axes[0] | |
| 242 | + ax.plot(d['trade_date'], d['implied_corr'], color=fs.SERIES[0], lw=0.9, ls='-', label='Implied correlation (29 constituents, CBOE-style)') | |
| 243 | + ax.plot(d['trade_date'], d['realized_corr'], color=fs.SERIES[1], lw=0.9, ls='--', label='Realized correlation (22-day pairwise)') | |
| 244 | + ax.set_ylabel('Correlation') | |
| 245 | + ax.set_ylim(0, 1) | |
| 246 | + ax.legend(loc='upper right', ncol=2) | |
| 247 | + fs.panel_label(ax, '(a) Implied and realized correlation among S&P 500 constituents') | |
| 248 | + ax = axes[1] | |
| 249 | + ax.plot(d['trade_date'], d['log_ratio'], color=fs.SERIES[2], lw=0.8) | |
| 250 | + ax.axhline(0, color=fs.INK, lw=0.7) | |
| 251 | + ax.set_ylabel(r'$\ln(\hat\rho_{impl}/\hat\rho_{real})$') | |
| 252 | + fs.panel_label(ax, '(b) Log ratio of implied to realized correlation') | |
| 253 | + ax = axes[2] | |
| 254 | + ax.plot(d['trade_date'], d['vix_close'], color=fs.INK, lw=0.8) | |
| 255 | + ax.axhline(25, color=fs.SERIES[1], lw=0.8, ls='--') | |
| 256 | + ax.text(d['trade_date'].iloc[len(d) // 2], 27, 'VIX = 25 (stress threshold)', fontsize=7.5, color=fs.SERIES[1], va='bottom', ha='center', | |
| 257 | + bbox=dict(facecolor='white', edgecolor='none', pad=0.6)) | |
| 258 | + ax.set_ylabel('VIX close') | |
| 259 | + ax.set_xlabel('Trading day') | |
| 260 | + fs.panel_label(ax, '(c) VIX') | |
| 261 | + fig.legend(handles=[Patch(facecolor=fs.OI['grey'], alpha=0.25, label='Days with VIX close above 25')], loc='outside lower center') | |
| 262 | + fs.save(fig, 'fig_rq3_series') | |
| 263 | + n = len(d); clip0 = float((d['implied_corr_raw'] <= 0).mean() * 100) if 'implied_corr_raw' in d else np.nan | |
| 264 | + share_hi = float(hi_vix.mean() * 100) | |
| 265 | + cap = (rf"Implied and realized correlation among S\&P~500 constituents, 2010--2025 ({num(n)} trading days, series rebuilt from the shipped parquets by stage~1b of " | |
| 266 | + rf"\texttt{{scripts/04}}, 29 constituents, days with fewer than 20 constituents dropped). Panel (a): the CBOE-style implied correlation from the ATM implied " | |
| 267 | + rf"variances of SPX and of the constituents (eq.~(\ref{{eq:implcorr}}), clipped to $[0,1]$) and the 22-day pairwise realized correlation; panel (b): the log ratio " | |
| 268 | + rf"that enters the stress regressions (mean {_txt(d['log_ratio'].mean(),3)}, standard deviation {_txt(d['log_ratio'].std(),3)}); panel (c): the VIX close from the " | |
| 269 | + rf"complete daily series of the FinData OS lake. Shaded: days on which the VIX closes above 25 ({_txt(share_hi,1)}\% of days), the level component of the " | |
| 270 | + rf"post-origin stress definition. Source: \texttt{{data/processed/correlation\_divergence\_v2.parquet}} (regenerable; diagnostics in Table~\ref{{tab:rq3_diag}}).") | |
| 271 | + figure('fig_rq3_series', cap, 'fig:rq3_series') | |
| 272 | + | |
| 273 | + | |
| 274 | +# ---------------------------------------------------------------------------- | |
| 275 | +# core: IRFs with bands | |
| 276 | +# ---------------------------------------------------------------------------- | |
| 277 | +def fig_irf_bands(): | |
| 278 | + df = pd.read_csv(R / "irf_bands.csv") | |
| 279 | + df = df[df['spec'] == 'levels'] | |
| 280 | + fig, axes = fs.new(width=6.5, height=2.9, nrows=1, ncols=2, sharex=True) | |
| 281 | + panels = [('rv_to_ivshock', '(a) RV response to a one-SD IV shock'), ('iv_to_rvshock', '(b) IV response to a one-SD RV shock')] | |
| 282 | + for ax, (key, title) in zip(axes, panels): | |
| 283 | + for i, (ordering, label) in enumerate((('rvfirst', 'RV ordered first'), ('ivfirst', 'IV ordered first'))): | |
| 284 | + sub = df[df['ordering'] == ordering].sort_values('horizon') | |
| 285 | + ax.fill_between(sub['horizon'], sub[f'{key}_p05'], sub[f'{key}_p95'], color=fs.SERIES[i], alpha=0.18, lw=0, hatch=None if i == 0 else '///', edgecolor=fs.SERIES[i]) | |
| 286 | + ax.plot(sub['horizon'], sub[f'{key}_mean'], color=fs.SERIES[i], ls=fs.LINESTYLES[i], marker=fs.MARKERS[i], markevery=4, label=label) | |
| 287 | + ax.axhline(0, color=fs.INK, lw=0.7) | |
| 288 | + fs.panel_label(ax, title) | |
| 289 | + ax.set_xlabel('Horizon (trading days)') | |
| 290 | + ax.margins(x=0.02) | |
| 291 | + axes[0].set_ylabel('Response (standard-deviation units)') | |
| 292 | + axes[0].legend(loc='upper right') | |
| 293 | + fs.save(fig, 'fig_irf_bands') | |
| 294 | + n = int(df['n_tickers'].max()) | |
| 295 | + rv1 = df[(df['ordering'] == 'rvfirst')].sort_values('horizon') | |
| 296 | + peak = rv1.loc[rv1['rv_to_ivshock_mean'].idxmax()] | |
| 297 | + cap = (rf"Orthogonalized impulse responses of the bivariate VAR(5) in ATM implied volatility and realized variance (standardized levels), averaged across the " | |
| 298 | + rf"{n} tickers with at least 1{{,}}000 complete observations. Lines: cross-ticker mean response to a one-standard-deviation shock under the two recursive " | |
| 299 | + rf"orderings (solid, circles: RV ordered first, so that RV shocks move IV contemporaneously but not the reverse; dashed, squares: IV first); shaded bands " | |
| 300 | + rf"(plain / hatched): cross-ticker average of the per-ticker residual-bootstrap 5th and 95th percentiles (500 draws, seed 42). " | |
| 301 | + rf"Under the RV-first ordering the RV response to an IV shock peaks at {_txt(peak['rv_to_ivshock_mean'],3)} SD at horizon {int(peak['horizon'])} " | |
| 302 | + rf"(band [{_txt(peak['rv_to_ivshock_p05'],3)}, {_txt(peak['rv_to_ivshock_p95'],3)}]). Variance decompositions and stationarity diagnostics are in " | |
| 303 | + rf"Table~\ref{{tab:fevd_v2}}. Source: \texttt{{results/irf\_bands.csv}} (package L1-E8).") | |
| 304 | + figure('fig_irf_bands', cap, 'fig:irf_bands') | |
| 305 | + | |
| 306 | + | |
| 307 | +# ---------------------------------------------------------------------------- | |
| 308 | +# appendix: specification curve | |
| 309 | +# ---------------------------------------------------------------------------- | |
| 310 | +def fig_spec_curve(): | |
| 311 | + d = pd.read_csv(R / "rq1_spec_curve.csv").sort_values('t_stat').reset_index(drop=True) | |
| 312 | + d['x'] = np.arange(len(d)) | |
| 313 | + se_types = list(dict.fromkeys(d['se_type'])) | |
| 314 | + dims = [('se_type', 'Standard errors'), ('cutoff', 'Winsorization'), ('sample', 'Sample'), ('horizon_days', 'Horizon (days)'), ('controls', 'Controls'), ('return_basis', 'Return basis')] | |
| 315 | + levels = [(k, lab, [str(v) for v in dict.fromkeys(d[k].astype(str))]) for k, lab, in dims] | |
| 316 | + n_rows = sum(len(v) for _, _, v in levels) | |
| 317 | + fig, (ax, axd) = fs.new(width=6.5, height=5.2, nrows=2, ncols=1, sharex=True, gridspec_kw={'height_ratios': [2.2, n_rows / 6]}) | |
| 318 | + for i, se in enumerate(se_types): | |
| 319 | + s = d[d['se_type'] == se] | |
| 320 | + ax.scatter(s['x'], s['t_stat'], s=9, marker=fs.MARKERS[i % 6], color=fs.SERIES[i % 6], edgecolor=fs.INK, linewidth=0.2, label=se.replace('_', ' ')) | |
| 321 | + for yv, lab in ((1.96, '$|t| = 1.96$'), (-1.96, None), (3.02, '$|t| = 3.02$ (Bonferroni threshold of version 1.2)'), (-3.02, None)): | |
| 322 | + ax.axhline(yv, color=fs.INK2, lw=0.7, ls=':' if abs(yv) > 3 else '--') | |
| 323 | + if lab: | |
| 324 | + ax.text(len(d) * 0.01, yv + 0.1, lab, fontsize=7, color=fs.INK2, va='bottom') | |
| 325 | + ax.axhline(0, color=fs.INK, lw=0.7) | |
| 326 | + ax.set_ylabel('$t$-statistic of Skew$_{25\\delta}$') | |
| 327 | + ax.legend(loc='upper left', ncol=2, fontsize=7) | |
| 328 | + fs.panel_label(ax, f'(a) Sorted $t$-statistics over {len(d)} specifications') | |
| 329 | + yt, ytl = [], [] | |
| 330 | + row = 0 | |
| 331 | + for k, lab, vals in levels: | |
| 332 | + for v in vals: | |
| 333 | + s = d[d[k].astype(str) == v] | |
| 334 | + axd.scatter(s['x'], np.full(len(s), row), s=3, marker='|', color=fs.INK, linewidth=0.6) | |
| 335 | + yt.append(row); ytl.append(f"{lab}: {v.replace('_', ' ')}" if v == vals[0] else f"{v.replace('_', ' ')}") | |
| 336 | + row += 1 | |
| 337 | + axd.axhline(row - 0.5, color=fs.RULE, lw=0.5) | |
| 338 | + axd.set_yticks(yt); axd.set_yticklabels(ytl, fontsize=6.5) | |
| 339 | + axd.set_ylim(-0.6, row - 0.4); axd.invert_yaxis() | |
| 340 | + axd.set_xlabel('Specifications ordered by $t$-statistic') | |
| 341 | + axd.grid(False) | |
| 342 | + fs.panel_label(axd, '(b) Specification indicators') | |
| 343 | + fs.save(fig, 'fig_spec_curve') | |
| 344 | + cap = (rf"Specification curve of the 25-delta skew coefficient in the return regressions ({len(d)} specifications: winsorization cutoff $\times$ standard-error " | |
| 345 | + rf"estimator $\times$ sample $\times$ horizon $\times$ controls $\times$ return basis; Table~\ref{{tab:spec_curve}}). Panel (a): $t$-statistics sorted in increasing " | |
| 346 | + rf"order, marker and colour by standard-error estimator; dashed lines at $\pm1.96$, dotted at $\pm3.02$. Panel (b): the specification behind each point. " | |
| 347 | + rf"The median $t$ is {_txt(d['t_stat'].median())}; $|t|>1.96$ in {_txt(d['abs_t_gt_196'].mean()*100,1)}\% and $|t|>3.02$ in {_txt(d['abs_t_gt_302'].mean()*100,1)}\% " | |
| 348 | + rf"of specifications; among the Driscoll--Kraay specifications the share with $|t|>3.02$ is " | |
| 349 | + rf"{_txt(d[d['se_type'].str.startswith('DK')]['abs_t_gt_302'].mean()*100,1)}\%. Source: \texttt{{results/rq1\_spec\_curve.csv}} (package L1-C5).") | |
| 350 | + figure('fig_spec_curve', cap, 'fig:spec_curve') | |
| 351 | + | |
| 352 | + | |
| 353 | +# ---------------------------------------------------------------------------- | |
| 354 | +# appendix: price magnet | |
| 355 | +# ---------------------------------------------------------------------------- | |
| 356 | +def fig_magnet(): | |
| 357 | + d = pd.read_csv(R / "rq4_price_magnet_v2.csv") | |
| 358 | + fig, axes = fs.new(width=6.5, height=2.9, nrows=1, ncols=2) | |
| 359 | + for ax, (kind, title, strip) in zip(axes, [('distance', '(a) By distance of the open to the strike', 'Distance '), ('oi_quintile', '(b) By open-interest concentration quintile', 'OI quintile ')]): | |
| 360 | + s = d[d['kind'] == kind].reset_index(drop=True) | |
| 361 | + x = np.arange(len(s)) | |
| 362 | + ax.errorbar(x - 0.12, s['null_share_mean'] * 100, yerr=[(s['null_share_mean'] - s['null_p025']) * 100, (s['null_p975'] - s['null_share_mean']) * 100], | |
| 363 | + fmt=fs.MARKERS[1], color=fs.SERIES[1], capsize=3, lw=1, label='Permutation null (mean, 95% interval)', markeredgecolor=fs.INK, markeredgewidth=0.3) | |
| 364 | + ax.errorbar(x + 0.12, s['observed_share'] * 100, yerr=1.96 * s['se_observed_date_clustered'] * 100, | |
| 365 | + fmt=fs.MARKERS[0], color=fs.SERIES[0], capsize=3, lw=1, label='Observed share ($\\pm 1.96$ date-clustered SE)', markeredgecolor=fs.INK, markeredgewidth=0.3) | |
| 366 | + ax.plot(x, s['mirror_placebo_share'] * 100, ls='', marker='x', color=fs.INK2, label='Mirror-strike placebo') | |
| 367 | + ax.axhline(50, color=fs.INK, lw=0.7, ls=':') | |
| 368 | + ax.set_xticks(x); ax.set_xticklabels([g.replace(strip, '').replace('Q1_Low', 'Q1 (low)').replace('Q5_High', 'Q5 (high)') for g in s['group']], rotation=0, fontsize=7.5) | |
| 369 | + ax.set_ylabel('Days closing nearer the max-OI strike (%)') | |
| 370 | + fs.panel_label(ax, title) | |
| 371 | + ax.margins(y=0.2) | |
| 372 | + for xx, (n, p) in enumerate(zip(s['n_obs'], s['p_value_perm'])): | |
| 373 | + ax.annotate(f"n = {int(n):,}\np = {p:.2f}", xy=(xx, ax.get_ylim()[0]), xytext=(0, 2), textcoords='offset points', ha='center', va='bottom', fontsize=6.3, color=fs.INK2) | |
| 374 | + axes[0].legend(loc='center right', fontsize=6.8) | |
| 375 | + fs.save(fig, 'fig_magnet') | |
| 376 | + a = d[d['kind'] == 'all'].iloc[0] | |
| 377 | + cap = (rf"Price-magnet test on the shipped panel ({num(a['n_obs'])} ticker-days of eight stocks, 2010--2025, open within 0.1\%--10\% of the maximum-open-interest " | |
| 378 | + rf"strike among options expiring within 30 days). Circles: observed share of days on which the close is relatively closer to that strike than the open, with " | |
| 379 | + rf"$\pm1.96$ date-clustered standard errors; squares: the permutation null that re-pairs each (day, distance, strike side) with the open-to-close move of a random " | |
| 380 | + rf"other day of the same ticker ({int(a['n_perm'])} draws, seed {int(a['seed'])}), mean and 95\% interval; crosses: the mirror-strike placebo. The dotted line is the " | |
| 381 | + rf"50\% benchmark of version~1.2, which is not the correct null because a driftless move reaches a strike at distance $d$ only if it has the right sign and does not " | |
| 382 | + rf"overshoot. Over all ticker-days the observed share is {_txt(a['observed_share']*100,2)}\% against a null of {_txt(a['null_share_mean']*100,2)}\% " | |
| 383 | + rf"[{_txt(a['null_p025']*100,1)}, {_txt(a['null_p975']*100,1)}] (permutation $p$ = {a['p_value_perm']:.2f}; minimum detectable difference " | |
| 384 | + rf"{_txt(a['mde_80pct_power_pp'],1)} pp at 80\% power). Source: \texttt{{results/rq4\_price\_magnet\_v2.csv}} (Table~\ref{{tab:magnet_v2}}).") | |
| 385 | + figure('fig_magnet', cap, 'fig:magnet') | |
| 386 | + | |
| 387 | + | |
| 388 | +# ---------------------------------------------------------------------------- | |
| 389 | +# appendix: RQ5 clipping dependence | |
| 390 | +# ---------------------------------------------------------------------------- | |
| 391 | +def fig_rq5_clipping(): | |
| 392 | + c = pd.read_csv(R / "rq5_clipping_variant.csv") | |
| 393 | + c = c[c['period'] != '2016-2025'].copy() | |
| 394 | + c['year'] = c['period'].astype(int) | |
| 395 | + rules = list(dict.fromkeys(c['rule'])) | |
| 396 | + rule_lab = {'features clipped at training quantiles': 'HAR-RV, regressors clipped at training 1%/99% quantiles', | |
| 397 | + 'features unclipped': 'HAR-RV, regressors unclipped', 'log-RV model': 'log-RV HAR specification'} | |
| 398 | + fig, ax = fs.new(width=6.5, height=3.0) | |
| 399 | + years = sorted(c['year'].unique()) | |
| 400 | + for i, rule in enumerate(rules): | |
| 401 | + s = c[(c['rule'] == rule) & (c['model'] == 'OLS: HAR-RV')].sort_values('year') | |
| 402 | + ax.plot(s['year'], s['r2_oos'], color=fs.SERIES[i], marker=fs.MARKERS[i], ls=fs.LINESTYLES[i], label=rule_lab.get(rule, rule)) | |
| 403 | + v = c[(c['rule'] == rules[0]) & (c['model'] == 'OLS: HAR-RV')].sort_values('year') | |
| 404 | + ax.plot(v['year'], v['r2_oos_vix'], color=fs.INK, marker='x', ls='-', lw=1.0, label='VIX benchmark, $(VIX_t/100)^2/252$') | |
| 405 | + ax.axhline(0, color=fs.INK, lw=0.7) | |
| 406 | + # axis truncated below the lowest HAR value so that the HAR lines stay readable; off-scale VIX values are marked | |
| 407 | + floor = float(c[c['model'] == 'OLS: HAR-RV']['r2_oos'].min()) - 0.35 | |
| 408 | + ax.set_ylim(floor, max(float(c['r2_oos'].max()), float(v['r2_oos_vix'].max())) + 0.25) | |
| 409 | + for _, rr in v[v['r2_oos_vix'] < floor].iterrows(): | |
| 410 | + ax.annotate(f"VIX {int(rr['year'])}: {rr['r2_oos_vix']:.2f} (off scale)", xy=(rr['year'], floor), xytext=(0, 14), textcoords='offset points', | |
| 411 | + ha='center', va='bottom', fontsize=7, color=fs.INK2, arrowprops=dict(arrowstyle='-|>', color=fs.INK2, lw=0.7)) | |
| 412 | + ax.set_xticks(years) | |
| 413 | + ax.set_xlabel('Test year (models re-estimated on 2010 to the previous year); vertical axis truncated below the lowest HAR-RV value') | |
| 414 | + ax.set_ylabel('$R^2_{OOS}$, next-day SPX realized variance') | |
| 415 | + fig.legend(loc='outside lower center', ncol=2, fontsize=7.5) | |
| 416 | + fs.save(fig, 'fig_rq5_clipping') | |
| 417 | + full = pd.read_csv(R / "rq5_clipping_variant.csv"); full = full[full['period'] == '2016-2025'] | |
| 418 | + g = lambda rule, m='OLS: HAR-RV': full[(full['rule'] == rule) & (full['model'] == m)].iloc[0] | |
| 419 | + clip, unc = g(rules[0]), g('features unclipped') | |
| 420 | + y20c = c[(c['rule'] == rules[0]) & (c['model'] == 'OLS: HAR-RV') & (c['year'] == 2020)]['r2_oos'].iloc[0] | |
| 421 | + y20u = c[(c['rule'] == 'features unclipped') & (c['model'] == 'OLS: HAR-RV') & (c['year'] == 2020)]['r2_oos'].iloc[0] | |
| 422 | + cap = (rf"Dependence of the RQ5 benchmark comparison on the treatment of the test-year regressors (SPX, next-day realized variance, rolling annual re-estimation " | |
| 423 | + rf"2016--2025, FinData OS lake data, \texttt{{scripts/32}}). Lines: out-of-sample $R^2$ of the HAR-RV model by test year (vertical axis truncated below the lowest HAR-RV value; off-scale values of the VIX benchmark are marked with their value) when the test-year regressors are clipped at " | |
| 424 | + rf"the training window's 1\%/99\% quantiles (the scheme of the RQ2 design), left unclipped, or when the model is estimated in logs; crosses: the VIX benchmark. " | |
| 425 | + rf"Clipping caps $RV_t$ in 2020 far below its realized level: the HAR-RV $R^2_{{OOS}}$ for 2020 is {_txt(y20c,3)} clipped against {_txt(y20u,3)} unclipped, and over " | |
| 426 | + rf"2016--2025 it is {_txt(clip['r2_oos'],3)} clipped against {_txt(unc['r2_oos'],3)} unclipped, versus {_txt(clip['r2_oos_vix'],3)} for the VIX benchmark " | |
| 427 | + rf"(Diebold--Mariano $t$ of the unclipped HAR-RV against the VIX: {_txt(unc['dm_vs_vix_t'])}). Source: \texttt{{results/rq5\_clipping\_variant.csv}} (Table~\ref{{tab:rq5_v2}}, Panel C).") | |
| 428 | + figure('fig_rq5_clipping', cap, 'fig:rq5_clipping') | |
| 429 | + | |
| 430 | + | |
| 431 | +def main(): | |
| 432 | + print("=" * 70); print("MANUSCRIPT FIGURES (Stage E)"); print("=" * 70) | |
| 433 | + config.ensure_output_dirs() | |
| 434 | + for fn in [fig_coverage, fig_rq2_loss_ratios, fig_rq2_subperiod, fig_rq3_headline, fig_rq3_series, fig_irf_bands, | |
| 435 | + fig_spec_curve, fig_magnet, fig_rq5_clipping]: | |
| 436 | + print(fn.__name__) | |
| 437 | + fn() | |
| 438 | + print("\nFIGURES COMPLETE.") | |
| 439 | + | |
| 440 | + | |
| 441 | +if __name__ == "__main__": | |
| 442 | + main() | |
| 443 | ||