Dickey–Fuller, augmented Dickey–Fuller and Engle–Granger from scratch

"""M13A · Dickey–Fuller, augmented Dickey–Fuller and Engle–Granger from scratch. Read notes §13A.2–13A.5 first.

ADF regression with p lagged differences (notes eq. 13A.7):
    Delta y_t = [c] + [d t] + gamma y_{t-1} + sum_{j=1}^p delta_j Delta y_{t-j} + e_t
The ADF statistic is the OLS t-ratio of gamma. The effective sample is
t = p+2, ..., T in 1-based time (T - p - 1 equations), which is the convention
of statsmodels.tsa.stattools.adfuller(maxlag=p, autolag=None). The tests compare
your statistic with statsmodels to 1e-9.

Run (tests import this file; set COURSE_IMPL=solutions to run the reference):
    cd exercises && .venv/bin/python -m pytest modules/M13A-cointegration-stat-arb/tests/test_m13a_adf.py
"""

from __future__ import annotations

import numpy as np

from toolkit.stats import ols as _ols

TRENDS = ("n", "c", "ct")


def ols(y: np.ndarray, X: np.ndarray):
    """OLS coefficients, classical standard errors and residuals (toolkit.stats.ols). (Given.)

    Returns (beta, se, resid) for the design X as given (no constant is added).
    """
    fit = _ols(y, X)
    return fit.params, fit.bse, fit.resid


def adf_design(y: np.ndarray, lags: int = 0, trend: str = "c"):
    """Build the ADF regression (eq. 13A.7).

    Returns (dy, X) where dy has length T - lags - 1 and the FIRST column of X is
    y_{t-1}; then the `lags` lagged differences Delta y_{t-1}, ..., Delta y_{t-lags};
    then a constant if trend in {"c", "ct"}; then a linear time index 1, 2, ... if trend == "ct".
    Raise ValueError for a trend not in TRENDS.
    Hint: with dy_all = np.diff(y), use rows i = lags, ..., len(dy_all)-1; the lagged level for
    row i is y[i] and the j-th lagged difference is dy_all[i - j].
    """
    raise NotImplementedError


def adf_stat(y: np.ndarray, lags: int = 0, trend: str = "c") -> float:
    """ADF t-statistic on gamma (eq. 13A.7). lags = 0 gives the Dickey–Fuller test."""
    raise NotImplementedError


def df_gamma(y: np.ndarray, trend: str = "c") -> float:
    """OLS gamma_hat = rho_hat - 1 from the DF regression with no lagged differences."""
    raise NotImplementedError


def simulate_df_quantiles(T: int, nsim: int, trend: str = "c", seed: int = 0,
                          probs=(0.01, 0.05, 0.10)) -> np.ndarray:
    """Monte Carlo quantiles of the DF t-statistic under H0 (Prop. 13A.1).

    Draw nsim Gaussian random walks of length T started at 0 (cumulative sums of
    np.random.default_rng(seed).standard_normal((nsim, T))), compute adf_stat(., 0, trend)
    for each, and return np.quantile(stats, probs).
    """
    raise NotImplementedError


def spurious_regression(T: int, nsim: int, seed: int = 0):
    """Granger–Newbold experiment (notes §13A.3).

    nsim times: draw two independent Gaussian random walks y, x of length T, regress y on
    [1, x] by OLS. Return (reject_rate, median_r2): the share of |t| > 1.96 on the slope and
    the median R^2 (R^2 = 1 - SSR / sum (y - ybar)^2).
    """
    raise NotImplementedError


def engle_granger(y: np.ndarray, x: np.ndarray, lags: int = 1):
    """Engle–Granger two-step (notes §13A.5).

    Step 1: OLS of y on a constant and x; residuals u_t.
    Step 2: ADF t-statistic on u with `lags` lagged differences and NO deterministic terms
    (trend "n"): the residuals already have mean zero.
    Returns (alpha_hat, beta_hat, stat). The tests check stat against statsmodels' coint.
    """
    raise NotImplementedError