"""High-frequency event study around macro events. Implements the event-study methodology of MacKinlay (1997) with Patell (1976) and Boehmer-Musumeci-Poulsen (1991) test statistics. """ from __future__ import annotations import pandas as pd import numpy as np from scipy import stats from dataclasses import dataclass # --------------------------------------------------------------------------- # Data classes # --------------------------------------------------------------------------- @dataclass class EventStudyResult: """Container for event-study output.""" cars: pd.DataFrame # [event_date, window, CAR] test_stats: pd.DataFrame # [window, mean_CAR, t_stat, p_value] event_type: str n_events: int # --------------------------------------------------------------------------- # Core event study # --------------------------------------------------------------------------- def run_event_study( prices: pd.DataFrame, event_dates: list, event_type: str = "FOMC", window: tuple[int, int] = (-5, 5), estimation_window: int = 60, ) -> EventStudyResult: """Event study methodology following MacKinlay (1997). 1. For each event date, find the row in *prices* closest to that date. 2. Estimation window: ``[t - estimation_window, t + window[0] - 1]`` for the normal-return (constant-mean) model. 3. Event window: ``[t + window[0], t + window[1]]``. 4. Abnormal return = actual return - expected return. 5. CAR = cumulative sum of abnormal returns in the event window. 6. Test statistics: cross-sectional *t*-test of CARs. Parameters ---------- prices : pd.DataFrame Each column is a market's price series with a DatetimeIndex. event_dates : list Dates (or date-like objects) of the events. event_type : str Label for the event category (e.g. ``"FOMC"``). window : tuple[int, int] (pre, post) trading-day offsets around each event. estimation_window : int Number of trading days used to estimate the normal-return model. Returns ------- EventStudyResult """ # Compute simple returns returns = prices.pct_change().dropna(how="all") # Ensure DatetimeIndex if not isinstance(returns.index, pd.DatetimeIndex): returns.index = pd.to_datetime(returns.index) event_dates_dt = pd.to_datetime(event_dates) all_cars: list[dict] = [] all_ars: list[pd.DataFrame] = [] for ev_date in event_dates_dt: # Find closest available trading day idx = returns.index.searchsorted(ev_date) if idx >= len(returns.index): idx = len(returns.index) - 1 t0 = idx # position index est_start = t0 - estimation_window + window[0] - 1 est_end = t0 + window[0] - 1 ev_start = t0 + window[0] ev_end = t0 + window[1] # Boundary checks if est_start < 0 or ev_end >= len(returns): continue # Estimation-window returns → constant-mean model est_returns = returns.iloc[est_start:est_end + 1] expected = est_returns.mean() # per column est_std = est_returns.std() # Event-window abnormal returns ev_returns = returns.iloc[ev_start:ev_end + 1] ar = ev_returns.subtract(expected, axis=1) car = ar.cumsum() # Store CAR at end of window (last row) final_car = car.iloc[-1] for col in final_car.index: all_cars.append({ "event_date": ev_date, "market": col, "window": f"[{window[0]},{window[1]}]", "CAR": final_car[col], }) all_ars.append(ar) cars_df = pd.DataFrame(all_cars) n_events = int(cars_df["event_date"].nunique()) if len(cars_df) > 0 else 0 # Cross-sectional test statistics per market test_rows: list[dict] = [] if len(cars_df) > 0: for market, grp in cars_df.groupby("market"): car_vals = grp["CAR"].dropna() if len(car_vals) < 2: continue t_stat, p_val = stats.ttest_1samp(car_vals, 0) test_rows.append({ "market": market, "window": f"[{window[0]},{window[1]}]", "mean_CAR": float(car_vals.mean()), "t_stat": float(t_stat), "p_value": float(p_val), "n_events": len(car_vals), }) test_df = pd.DataFrame(test_rows) return EventStudyResult( cars=cars_df, test_stats=test_df, event_type=event_type, n_events=n_events, ) # --------------------------------------------------------------------------- # Patell (1976) test # --------------------------------------------------------------------------- def compute_patell_test( abnormal_returns: pd.DataFrame, std_errors: pd.Series, ) -> float: """Patell (1976) standardised-residual test. Parameters ---------- abnormal_returns : pd.DataFrame Abnormal returns (events x time). std_errors : pd.Series Estimation-window standard errors for each event. Returns ------- float Patell *z*-statistic. """ # Standardise each event's AR by its estimation-window std error sar = abnormal_returns.div(std_errors, axis=0) # Cumulative standardised AR csar = sar.sum(axis=1) # sum across event window per event n = len(csar) if n == 0: return np.nan z = csar.mean() / (csar.std() / np.sqrt(n)) return float(z) # --------------------------------------------------------------------------- # BMP (1991) test # --------------------------------------------------------------------------- def compute_bmp_test( standardized_ars: pd.DataFrame, ) -> tuple[float, float]: """Boehmer, Musumeci, Poulsen (1991) test for event-induced variance. Uses cross-sectional variance of standardised abnormal returns to produce a *t*-statistic robust to event-induced variance changes. Parameters ---------- standardized_ars : pd.DataFrame Standardised abnormal returns (events x time-steps in window). Returns ------- tuple[float, float] (BMP *t*-statistic, two-sided *p*-value). """ # Average SAR across the event window for each event avg_sar = standardized_ars.mean(axis=1) n = len(avg_sar) if n < 2: return (np.nan, np.nan) mean_sar = avg_sar.mean() se = avg_sar.std(ddof=1) / np.sqrt(n) t_stat = mean_sar / se if se > 0 else np.nan p_value = float(2 * stats.t.sf(abs(t_stat), df=n - 1)) if not np.isnan(t_stat) else np.nan return (float(t_stat), p_value)