""" WP 6 -- Intraday periodicity and trading-session analysis (revision, 2026-10-17). The original 01_analysis.py labelled FirstRate time stamps as UTC, whereas they are in U.S. Eastern time (the ES series has its daily Globex maintenance break at 17:00 and the VIX series starts at 03:15 / 09:30 Eastern). This script rebuilds the intraday analysis with an explicit America/New_York -> UTC conversion from one-minute FirstRate bars held locally (sibling repository WP 12, same vendor; path in WP6_MINUTE_DIR), aggregated to the paper's 5-minute grid over the common sample 2016-03-11 to 2026-03-11: * hourly profile of average absolute 5-minute returns, in UTC, for ES, NQ, VIX, BTC; * session realized variances (Asian 00:00-08:00, European 08:00-14:30, U.S. 14:30-21:00 UTC) and their shares of the daily total; * session-ordered predictive regressions: does Bitcoin's Asian/European-session volatility help predict the U.S.-session volatility of ES (and vice versa), and does weekend Bitcoin volatility predict Monday equity volatility? Outputs: results/rev_intraday_*.csv, paper/tables/tab_sessions*.tex, figures/fig5_hourly_volatility.png. Run time: about 2 minutes (loading 1-minute parquet files). """ import os import time import warnings import numpy as np import pandas as pd import matplotlib matplotlib.use('Agg') import matplotlib.pyplot as plt import statsmodels.api as sm warnings.filterwarnings('ignore') T0 = time.time() ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) RESULTS = os.path.join(ROOT, 'results') FIGURES = os.path.join(ROOT, 'figures') TABLES = os.path.join(ROOT, 'paper', 'tables') PAPER_FIG = os.path.join(ROOT, 'paper', 'figures') MINUTE_DIR = os.environ.get('WP6_MINUTE_DIR', os.path.join(ROOT, '..', 'realized-volatility-horse-race', 'data', 'raw')) FILES = {'ES': 'futures_ES_1min.parquet', 'NQ': 'futures_NQ_1min.parquet', 'VIX': 'index_VIX_1min.parquet', 'BTC': 'crypto_BTC_1min.parquet'} START, END = '2016-03-11', '2026-03-11' COLORS = {'ES': '#1f77b4', 'NQ': '#aec7e8', 'VIX': '#8c564b', 'BTC': '#9467bd'} SESSIONS = {'Asian': (0.0, 8.0), 'European': (8.0, 14.5), 'US': (14.5, 21.0)} def log(m): print(f"[{time.time() - T0:6.1f}s] {m}", flush=True) def load_5min_utc(sym): df = pd.read_parquet(os.path.join(MINUTE_DIR, FILES[sym]), columns=['datetime', 'close']) s = df.set_index('datetime')['close'].sort_index() s = s.loc[START:END] # FirstRate stamps are naive U.S. Eastern (bar-start convention); convert to UTC idx = s.index.tz_localize('America/New_York', ambiguous='NaT', nonexistent='NaT') s = s[~idx.isna()] s.index = idx[~idx.isna()].tz_convert('UTC') s5 = s.resample('5min').last().dropna() r = np.log(s5).diff().dropna() # drop returns spanning a gap longer than one hour (session breaks, weekends) gaps = s5.index.to_series().diff().dt.total_seconds().reindex(r.index) r = r[gaps <= 3600] return r rets = {} for sym in FILES: rets[sym] = load_5min_utc(sym) log(f"{sym}: {len(rets[sym]):,} five-minute returns, {rets[sym].index.min().date()} -> {rets[sym].index.max().date()}") # ─── 1. Hourly profile (UTC) ────────────────────────────────────────────────── hourly = pd.DataFrame({sym: r.abs().groupby(r.index.hour).mean() * 100 for sym, r in rets.items()}) hourly.index.name = 'hour_utc' hourly.to_csv(os.path.join(RESULTS, 'rev_intraday_hourly_absret_utc.csv'), float_format='%.5f') hourly_et = pd.DataFrame({sym: r.abs().groupby(r.index.tz_convert('America/New_York').hour).mean() * 100 for sym, r in rets.items()}) hourly_et.index.name = 'hour_eastern' hourly_et.to_csv(os.path.join(RESULTS, 'rev_intraday_hourly_absret_eastern.csv'), float_format='%.5f') fig, axes = plt.subplots(2, 2, figsize=(13, 7.5)) for ax, sym in zip(axes.flatten(), ['VIX', 'ES', 'NQ', 'BTC']): h = hourly[sym] ax.bar(h.index, h.values, color=COLORS[sym], alpha=0.85, edgecolor='white', linewidth=0.5) ax.axvspan(-0.5, 7.5, alpha=0.07, color='blue', label='Asian 00:00–08:00 UTC') ax.axvspan(7.5, 14.0, alpha=0.07, color='green', label='European 08:00–14:30 UTC') ax.axvspan(14.0, 20.5, alpha=0.07, color='orange', label='U.S. 14:30–21:00 UTC') ax.set_title(f'{sym} — mean |5-min return| by hour (UTC)') ax.set_xlabel('Hour of day (UTC)'); ax.set_ylabel('Mean |return| (%)') ax.set_xticks(range(0, 24, 2)); ax.set_xlim(-0.6, 23.6) if sym == 'VIX': ax.legend(fontsize=7.5, loc='upper right') plt.suptitle('Intraday periodicity of 5-minute absolute returns, 2016-03 to 2026-03 (time stamps converted from U.S. Eastern)', fontsize=11) plt.tight_layout() for d in (FIGURES, PAPER_FIG): fig.savefig(os.path.join(d, 'fig5_hourly_volatility.png'), dpi=150, bbox_inches='tight') plt.close(fig) log("hourly figure written") # ─── 2. Session realized variances per UTC day ──────────────────────────────── def session_rv(r): hour = r.index.hour + r.index.minute / 60 day = r.index.tz_convert('UTC').normalize().tz_localize(None) out = {} for name, (a, b) in SESSIONS.items(): m = (hour >= a) & (hour < b) rv = (r[m] ** 2).groupby(day[m]).sum() n = r[m].groupby(day[m]).size() rv = rv[n >= 12] # at least one hour of bars out[name] = rv out['Day'] = (r ** 2).groupby(day).sum() return pd.DataFrame(out) sess = {sym: session_rv(r) for sym, r in rets.items()} rows = [] for sym, df in sess.items(): wd = df[df.index.dayofweek < 5] for name in SESSIONS: s = wd[name].dropna() rows.append({'Asset': sym, 'Session': name, 'Days': len(s), 'Ann_vol_pct': 100 * np.sqrt(252 * s).mean(), 'Mean_absret_pct': np.nan, 'Share_of_day_pct': 100 * (wd[name] / wd[['Asian', 'European', 'US']].sum(axis=1)).mean()}) sess_tbl = pd.DataFrame(rows) # mean absolute return per session for sym, r in rets.items(): hour = r.index.hour + r.index.minute / 60 for name, (a, b) in SESSIONS.items(): m = (hour >= a) & (hour < b) & (r.index.dayofweek < 5) sess_tbl.loc[(sess_tbl.Asset == sym) & (sess_tbl.Session == name), 'Mean_absret_pct'] = 100 * r[m].abs().mean() sess_tbl.to_csv(os.path.join(RESULTS, 'rev_intraday_sessions.csv'), index=False, float_format='%.4f') log("session table done") # weekend share for BTC: average over ISO weeks of (Saturday + Sunday variance) / (weekly variance), # which is robust to the extreme 2017 days; day-of-week profile as mean annualized daily RV btc_day = sess['BTC']['Day'] wk_var = btc_day.groupby(btc_day.index.to_period('W')).sum() wk_end = btc_day[btc_day.index.dayofweek >= 5].groupby(btc_day[btc_day.index.dayofweek >= 5].index.to_period('W')).sum() ratio = (wk_end / wk_var).dropna() ratio = ratio[btc_day.groupby(btc_day.index.to_period('W')).size() == 7] wk_share = ratio.mean() btc_rv_day = 100 * np.sqrt(252 * btc_day) btc_dow = btc_rv_day.groupby(btc_rv_day.index.dayofweek).mean() btc_dow.to_csv(os.path.join(RESULTS, 'rev_intraday_btc_dayofweek.csv'), header=['ann_vol_pct']) log(f"BTC weekend share of weekly variance: {100 * wk_share:.1f}%") # ─── 3. Session-ordered predictive regressions ─────────────────────────────── def predictive(dep, pred_same_day=None, controls_lag=4, label=''): """OLS of log dep RV on its own lags (+ optional same-day predictor and its lags); HAC Wald on predictor block.""" df = pd.DataFrame({'y': np.log(dep.where(dep > 0))}) X = [df['y'].shift(l).rename(f'y_l{l}') for l in range(1, controls_lag + 1)] pred_cols = [] if pred_same_day is not None: x = np.log(pred_same_day.where(pred_same_day > 0)).reindex(df.index) for l in range(0, controls_lag): X.append(x.shift(l).rename(f'x_l{l}')) pred_cols.append(f'x_l{l}') X = sm.add_constant(pd.concat(X, axis=1)) data = pd.concat([df['y'], X], axis=1).dropna() m = sm.OLS(data['y'], data.drop(columns='y')).fit(cov_type='HAC', cov_kwds={'maxlags': 10}) R = np.zeros((len(pred_cols), m.params.shape[0])) for k, c in enumerate(pred_cols): R[k, list(m.params.index).index(c)] = 1 w = m.wald_test(R, use_f=False, scalar=True) base = sm.OLS(data['y'], data.drop(columns=['y'] + pred_cols)).fit() return {'Test': label, 'T': int(m.nobs), 'coef_x_l0': float(m.params['x_l0']), 'se_x_l0': float(m.bse['x_l0']), 'Wald_HAC': float(w.statistic), 'p': float(w.pvalue), 'R2_with': m.rsquared, 'R2_without': base.rsquared} wd = {sym: df[df.index.dayofweek < 5] for sym, df in sess.items()} tests = [] # (i) Bitcoin earlier sessions -> ES / NQ / VIX U.S. session (same day) for tgt in ['ES', 'NQ', 'VIX']: for s in ['Asian', 'European']: tests.append(predictive(wd[tgt]['US'].dropna(), wd['BTC'][s], label=f'BTC {s} (t) $\\rightarrow$ {tgt} U.S. (t)')) # (ii) ES / VIX U.S. session (t-1) -> Bitcoin Asian session (t): shift predictor by one day for src in ['ES', 'NQ', 'VIX']: prev = wd[src]['US'].dropna().copy() prev.index = prev.index + pd.offsets.BDay(1) tests.append(predictive(wd['BTC']['Asian'].dropna(), prev, label=f'{src} U.S. (t$-$1) $\\rightarrow$ BTC Asian (t)')) # (iii) ES European (t) -> BTC U.S. (t) and BTC U.S. (t) -> ES U.S.(t) skipped: contemporaneous overlap # (iv) ES Asian (Globex overnight, t) -> BTC U.S. (t) for symmetry tests.append(predictive(wd['BTC']['US'].dropna(), wd['ES']['Asian'], label='ES Asian (t) $\\rightarrow$ BTC U.S. (t)')) tests.append(predictive(wd['ES']['US'].dropna(), wd['ES']['Asian'], label='ES Asian (t) $\\rightarrow$ ES U.S. (t) [own overnight]')) pred_df = pd.DataFrame(tests) pred_df.to_csv(os.path.join(RESULTS, 'rev_intraday_session_predictive.csv'), index=False, float_format='%.4f') log("session predictive regressions done") # (v) weekend Bitcoin variance -> Monday ES U.S.-session variance, controlling for Friday ES and Friday BTC es_us = wd['ES']['US'].dropna() mon = es_us[es_us.index.dayofweek == 0] rows = [] for d in mon.index: fri = d - pd.Timedelta(days=3) sat, sun = d - pd.Timedelta(days=2), d - pd.Timedelta(days=1) bwk = btc_day.get(sat, np.nan) + btc_day.get(sun, np.nan) rows.append({'date': d, 'es_mon': mon[d], 'es_fri': es_us.get(fri, np.nan), 'btc_wk': bwk, 'btc_fri': btc_day.get(fri, np.nan), 'btc_mon_asia': wd['BTC']['Asian'].get(d, np.nan)}) wk = pd.DataFrame(rows).set_index('date') wk = np.log(wk.where(wk > 0)).dropna() m1 = sm.OLS(wk['es_mon'], sm.add_constant(wk[['es_fri', 'btc_fri']])).fit(cov_type='HAC', cov_kwds={'maxlags': 4}) m2 = sm.OLS(wk['es_mon'], sm.add_constant(wk[['es_fri', 'btc_fri', 'btc_wk']])).fit(cov_type='HAC', cov_kwds={'maxlags': 4}) m3 = sm.OLS(wk['es_mon'], sm.add_constant(wk[['es_fri', 'btc_fri', 'btc_wk', 'btc_mon_asia']])).fit(cov_type='HAC', cov_kwds={'maxlags': 4}) wk_res = pd.DataFrame({ 'model': ['Friday controls only', '+ weekend BTC', '+ weekend BTC + Monday Asian BTC'], 'N': [int(m1.nobs), int(m2.nobs), int(m3.nobs)], 'coef_btc_wk': [np.nan, m2.params['btc_wk'], m3.params['btc_wk']], 'se_btc_wk': [np.nan, m2.bse['btc_wk'], m3.bse['btc_wk']], 'p_btc_wk': [np.nan, m2.pvalues['btc_wk'], m3.pvalues['btc_wk']], 'coef_btc_mon_asia': [np.nan, np.nan, m3.params['btc_mon_asia']], 'p_btc_mon_asia': [np.nan, np.nan, m3.pvalues['btc_mon_asia']], 'coef_es_fri': [m1.params['es_fri'], m2.params['es_fri'], m3.params['es_fri']], 'R2': [m1.rsquared, m2.rsquared, m3.rsquared]}) wk_res.to_csv(os.path.join(RESULTS, 'rev_intraday_weekend_btc_monday_es.csv'), index=False, float_format='%.4f') log("weekend regression done") # ─── 4. LaTeX tables ────────────────────────────────────────────────────────── def fmt(x, d=2, sign=False): if pd.isna(x): return '--' s = f"{x:+.{d}f}" if sign else f"{x:.{d}f}" return s.replace('-', '$-$').replace('+', '$+$') if sign else s def stars(p): return r'$^{***}$' if p < 0.01 else r'$^{**}$' if p < 0.05 else r'$^{*}$' if p < 0.10 else '' L = [r"\begin{tabular}{@{}llrrrr@{}}", r"\toprule", r"Asset & Session (UTC) & Days & Mean $|r|$ (\%) & Annualized vol.\ (\%) & Share of day (\%) \\", r"\midrule"] for sym in ['ES', 'NQ', 'VIX', 'BTC']: sub = sess_tbl[sess_tbl.Asset == sym] for k, (_, r) in enumerate(sub.iterrows()): lab = {'Asian': 'Asian 00:00--08:00', 'European': 'European 08:00--14:30', 'US': 'U.S. 14:30--21:00'}[r['Session']] L.append(f"{sym if k == 0 else ''} & {lab} & {int(r['Days'])} & {fmt(r['Mean_absret_pct'], 3)} & " f"{fmt(r['Ann_vol_pct'], 1)} & {fmt(r['Share_of_day_pct'], 1)} \\\\") if sym != 'BTC': L.append(r"\addlinespace") L += [r"\bottomrule", r"\end{tabular}"] with open(os.path.join(TABLES, 'tab_sessions.tex'), 'w') as f: f.write("\n".join(L)) L = [r"\begin{tabular}{@{}lrrrrr@{}}", r"\toprule", r"Predictive relation & $T$ & Coef.\ (same day) & Wald (HAC) & $p$ & $\Delta R^2$ \\", r"\midrule"] for _, r in pred_df.iterrows(): L.append(f"{r['Test']} & {int(r['T'])} & {fmt(r['coef_x_l0'], 3, True)} ({fmt(r['se_x_l0'], 3)}) & " f"{fmt(r['Wald_HAC'], 2)}{stars(r['p'])} & {fmt(r['p'], 3)} & {fmt(r['R2_with'] - r['R2_without'], 3)} \\\\") L += [r"\bottomrule", r"\end{tabular}"] with open(os.path.join(TABLES, 'tab_sessions_predictive.tex'), 'w') as f: f.write("\n".join(L)) L = [r"\begin{tabular}{@{}lrrrrr@{}}", r"\toprule", r"Specification & $N$ & $\ln \mathrm{RV}^{\mathrm{ES}}_{\mathrm{Fri}}$ & $\ln \mathrm{RV}^{\mathrm{BTC}}_{\mathrm{weekend}}$ & $\ln \mathrm{RV}^{\mathrm{BTC}}_{\mathrm{Mon,\,Asian}}$ & $R^2$ \\", r"\midrule"] for _, r in wk_res.iterrows(): bw = '--' if pd.isna(r['coef_btc_wk']) else f"{fmt(r['coef_btc_wk'], 3, True)}{stars(r['p_btc_wk'])} ({fmt(r['se_btc_wk'], 3)})" ba = '--' if pd.isna(r['coef_btc_mon_asia']) else f"{fmt(r['coef_btc_mon_asia'], 3, True)}{stars(r['p_btc_mon_asia'])}" L.append(f"{r['model']} & {int(r['N'])} & {fmt(r['coef_es_fri'], 3, True)} & {bw} & {ba} & {fmt(r['R2'], 3)} \\\\") L += [r"\bottomrule", r"\end{tabular}"] with open(os.path.join(TABLES, 'tab_weekend_btc.tex'), 'w') as f: f.write("\n".join(L)) key = {'btc_weekend_share_pct': 100 * wk_share} for sym in FILES: h = hourly[sym] key[f'{sym}_peak_hour_utc'] = int(h.idxmax()); key[f'{sym}_peak_absret'] = h.max() key[f'{sym}_min_hour_utc'] = int(h.idxmin()); key[f'{sym}_min_absret'] = h.min() key[f'{sym}_peak_to_trough'] = h.max() / h.min() for d, v in btc_dow.items(): key[f'BTC_dow{d}_annvol'] = v pd.Series(key).to_csv(os.path.join(RESULTS, 'rev_intraday_key_numbers.csv'), header=['value']) print(pd.Series(key).to_string()) print(sess_tbl.to_string()) print(pred_df.to_string()) print(wk_res.to_string()) log("done")