Earnings surprises and post-earnings-announcement drift (PEAD)

"""M06A · Earnings surprises and post-earnings-announcement drift (PEAD).

Exercise stubs. Replace each ``raise NotImplementedError`` with your code, then run the tests in
../tests. See notes.md §6A.2 (Definition 6A.1, Propositions 6A.2–6A.3).

Two parts:
1. SUE from a quarterly EPS series with the seasonal-random-walk-with-drift expectation model.
2. An event study of drift on a synthetic panel in which the drift is injected, so you can check
   that your estimator recovers what was put in.
"""

from __future__ import annotations

import numpy as np
import statsmodels.api as sm  # noqa: F401  (for drift_regression)


# ---------------------------------------------------------------------------------------------
# Provided helper (not a stub): synthetic event panel with injected drift.
# ---------------------------------------------------------------------------------------------
def make_event_panel(n_events: int = 20000, horizon: int = 60, day0: float = 0.015,
                     drift: float = 0.012, sigma: float = 0.02, seed: int = 0):
    """Synthetic abnormal returns around earnings announcements.

    SUE_i ~ N(0, 1) independently.  Abnormal returns, for event-time tau = 0..horizon:
        AR_{i,0}   = day0 * SUE_i + eps_{i,0}
        AR_{i,tau} = (drift / horizon) * SUE_i + eps_{i,tau},   tau = 1..horizon
    with eps iid N(0, sigma^2). So the expected post-announcement CAR(1, horizon) of an event with
    surprise s is drift * s.
    Returns (sue, ar) with sue of shape (n_events,) and ar of shape (n_events, horizon + 1);
    column j of ar is event day tau = j.
    """
    rng = np.random.default_rng(seed)
    sue = rng.standard_normal(n_events)
    ar = rng.normal(0.0, sigma, size=(n_events, horizon + 1))
    ar[:, 0] += day0 * sue
    ar[:, 1:] += (drift / horizon) * sue[:, None]
    return sue, ar


# ---------------------------------------------------------------------------------------------
# Exercises
# ---------------------------------------------------------------------------------------------
def sue(eps, window: int = 8, min_obs: int = 4) -> np.ndarray:
    """Standardised unexpected earnings, seasonal random walk with drift (Definition 6A.1).

    eps: quarterly EPS, oldest first. Let d_q = eps_q - eps_{q-4} (seasonal difference).
    For quarter q, take the previous `window` seasonal differences d_{q-window}, ..., d_{q-1}
    (ignore NaNs); if fewer than `min_obs` are available, or their std is 0, SUE_q = NaN.
    Otherwise delta_q = their mean, s_q = their sample std (ddof=1), SUE_q = (d_q - delta_q) / s_q.
    Returns an array the same length as eps (NaN where undefined).
    """
    raise NotImplementedError


def assign_groups(x, n_groups: int = 10) -> np.ndarray:
    """Rank-based group labels 0..n_groups-1 (0 = lowest x). Use ranks 0..n-1 (ties broken by
    order of appearance, i.e. a stable sort) and label = floor(rank * n_groups / n)."""
    raise NotImplementedError


def car(ar, start: int, end: int) -> np.ndarray:
    """Cumulative abnormal return per event over event days start..end inclusive (columns of ar)."""
    raise NotImplementedError


def group_car(ar, sue_values, start: int, end: int, n_groups: int = 10) -> np.ndarray:
    """Mean CAR(start, end) for each SUE group (length n_groups, group 0 = most negative SUE)."""
    raise NotImplementedError


def drift_spread(ar, sue_values, start: int = 1, end: int = 60, n_groups: int = 10):
    """Top-minus-bottom group difference in mean CAR(start, end), with its standard error and t.

    se = sqrt(var_top / n_top + var_bottom / n_bottom) with sample variances (ddof=1).
    Returns (spread, se, t).
    """
    raise NotImplementedError


def drift_regression(ar, sue_values, start: int = 1, end: int = 60):
    """OLS of CAR(start, end) on a constant and SUE. Returns (slope, t_slope) with the usual
    homoskedastic standard error. Use statsmodels: sm.OLS(y, sm.add_constant(x)).fit()."""
    raise NotImplementedError


def expected_extreme_spread(n_groups: int = 10) -> float:
    """E[Z | top group] - E[Z | bottom group] for Z ~ N(0,1) and equal-probability groups.
    Derive it (Proposition 6A.3, Problem 6A.3): E[Z | Z > c] = phi(c) / P(Z > c). Use scipy.stats.norm."""
    raise NotImplementedError