fithia2/libs/backtest/scenarios/coupling.py

"""Event-price coupling for synthetic scenario backtesting.

Controls the signal-to-noise ratio between event quality features and
subsequent price movements. This is the core mechanism for overfitting detection:

  signal_strength=0.0  → pure noise  → strategy should return ~0  (false positive check)
  signal_strength=0.35 → realistic SNR → strategy captures genuine alpha
  signal_strength=0.60 → strong signal → verify strategy responds to alpha

The coupling injects a drift into bars AFTER the execution date based on the
event's score and reaction features, while preserving OHLCV consistency.
"""
from __future__ import annotations

import datetime as dt
import math
from typing import Any

import numpy as np


def couple_events_to_prices(
    candidates: dict[dt.date, list[dict[str, Any]]],
    bars_by_symbol: dict[str, dict[dt.date, dict[str, Any]]],
    trading_dates: list[dt.date],
    signal_strength: float,
    signal_decay_days: int = 10,
    false_positive_rate: float = 0.15,
    rng: np.random.Generator | None = None,
) -> None:
    """Inject signal-driven drift into bars after each event's execution date.

    Modifies bars_by_symbol in-place. Also updates each candidate's
    entry_price / entry_price_est / event_close to match the actual
    bar close on its reaction_date (so stop prices are consistent).

    Args:
        candidates: candidates_by_exec_date dict from generate_events().
        bars_by_symbol: OHLCV bars to modify in-place.
        trading_dates: Ordered list of NYSE trading dates.
        signal_strength: 0.0 = pure noise, 0.35 = realistic, 0.6 = strong alpha.
        signal_decay_days: Days over which signal drift decays to zero.
        false_positive_rate: Fraction of qualifying events that produce negative
                             returns (traps / false positives).
        rng: NumPy random generator.
    """
    if rng is None:
        rng = np.random.default_rng()

    date_to_idx = {d: i for i, d in enumerate(trading_dates)}

    for exec_date, rows in candidates.items():
        exec_idx = date_to_idx.get(exec_date)
        if exec_idx is None:
            continue

        # reaction_date is the day before execution
        react_idx = exec_idx - 1
        if react_idx < 0:
            continue
        reaction_date = trading_dates[react_idx]

        for row in rows:
            symbol = str(row.get("symbol", ""))
            sym_bars = bars_by_symbol.get(symbol)
            if sym_bars is None:
                continue

            # Sync entry_price with actual bar close on reaction_date
            react_bar = sym_bars.get(reaction_date)
            if react_bar and react_bar.get("close", 0) > 0:
                actual_close = float(react_bar["close"])
                row["entry_price"] = round(actual_close, 4)
                row["entry_price_est"] = round(actual_close, 4)
                row["event_close"] = round(actual_close, 4)
                # Recompute ATR based on actual price
                atr_pct = float(row.get("atr_14", actual_close * 0.022)) / max(float(row.get("entry_price", actual_close)), 1e-4)
                row["atr_14"] = round(actual_close * atr_pct, 4)

            # Skip coupling if signal_strength == 0 (pure noise scenario)
            if signal_strength <= 1e-9:
                continue

            # Compute expected drift from event features
            score = float(row.get("score", 0.5))
            reaction_return = float(row.get("reaction_day_return", 0.0))
            volume_ratio = float(row.get("volume_ratio_20d", 1.5))

            expected_drift_5d = signal_strength * (
                0.030 * (score - 0.5)
                + 0.020 * reaction_return
                + 0.008 * max(0.0, volume_ratio - 1.5)
            )

            # False positive: flip signal direction
            if rng.random() < false_positive_rate:
                expected_drift_5d = -expected_drift_5d * 0.7

            if abs(expected_drift_5d) < 1e-6:
                continue

            # Distribute drift over signal_decay_days using a decay schedule
            total_drift = expected_drift_5d
            decay = _compute_decay_schedule(total_drift, signal_decay_days)

            # Inject drift into bars starting at exec_date + 1 (first day we hold)
            apply_start_idx = exec_idx + 1
            for k, daily_adj in enumerate(decay):
                bar_idx = apply_start_idx + k
                if bar_idx >= len(trading_dates):
                    break
                bar_date = trading_dates[bar_idx]
                bar = sym_bars.get(bar_date)
                if bar is None:
                    continue

                _apply_drift_to_bar(bar, daily_adj)


def _compute_decay_schedule(total_drift: float, decay_days: int) -> list[float]:
    """Distribute total_drift over decay_days using exponential decay.

    Returns a list of per-day drift adjustments that sum to total_drift.
    """
    if decay_days <= 0:
        return [total_drift]

    # Exponential decay weights
    weights = [math.exp(-0.5 * k / max(decay_days, 1)) for k in range(decay_days)]
    total_weight = sum(weights)
    return [total_drift * w / total_weight for w in weights]


def _apply_drift_to_bar(bar: dict[str, Any], daily_drift: float) -> None:
    """Multiply all OHLCV price fields by (1 + daily_drift), preserving consistency.

    Applies uniform multiplicative adjustment so that OHLC relationships are
    maintained exactly. Volume is unchanged.
    """
    factor = 1.0 + daily_drift
    factor = max(0.5, min(2.0, factor))  # guard against extreme values

    for field in ("open", "high", "low", "close"):
        val = bar.get(field)
        if val is not None and float(val) > 0:
            bar[field] = round(float(val) * factor, 4)

    # Ensure OHLCV consistency after adjustment
    o = bar.get("open", 0)
    h = bar.get("high", 0)
    lo = bar.get("low", 0)
    c = bar.get("close", 0)

    if o and h and lo and c:
        bar["high"] = round(max(float(h), float(o), float(c)), 4)
        bar["low"] = round(max(0.01, min(float(lo), float(o), float(c))), 4)