fithia2/apps/intraday_bt/scripts/diag_orb_obv_acceleration.py

"""
V34 Diagnostic: OBV Slope Acceleration

Tests whether SHORT-TERM OBV accumulation (5-day slope) provides signal orthogonal to
V24's existing 20-day OBV slope. Hypothesis: the most recent accumulation (last 5 days)
before a gap event captures fresher institutional positioning than the 20-day average.

Features:
  obv_slope_5    : 5-day OBV accumulation slope (same formula as obv_slope_20, shorter window)
  obv_slope_accel: obv_slope_5 - obv_slope_20 (positive = recent acceleration of accumulation)

The key test: does obv_slope_5 pass G5a (|ρ(slope_5, slope_20)| < 0.70)?
If highly correlated (ρ > 0.70), it's redundant with V24's signal. If orthogonal, it could
complement V24's score.

Context: V24 already uses obv_slope_20 at weight=0.05. V25-V33 all failed.
This is the final attempt using existing data (obv_slope_5 uses same formula, different window).
"""
from __future__ import annotations

import concurrent.futures
import datetime as dt
import os
import sys
from pathlib import Path

import pyarrow.parquet as pq
import yaml
from zoneinfo import ZoneInfo

sys.path.insert(0, os.path.abspath(os.path.join(os.path.dirname(__file__), "../../..")))

from libs.common.time_utils import trading_days_between
from libs.intraday.domain import ORBStrategyParams
from libs.intraday.features import compute_obv_slope_approx, enrich_daily_bars
from libs.intraday.orb_simulator import ORBSimulationState, run_orb_simulation_with_state
from libs.intraday.screener import orb_pre_screen_candidates

V24_CONFIG = "configs/intraday/strategies/orb_gainers_v24_quality_overlay.yaml"
UNIVERSE_FILE = "configs/symbols_midlarge_snapshot_exact.yaml"
INTRADAY_CACHE_DIR = "data/cache/intraday"
LOOKBACK_DAYS = 400

_ET = ZoneInfo("America/New_York")
_MKT_OPEN = dt.time(9, 30)
_MKT_CLOSE = dt.time(16, 0)


def _parse_ts(ts_raw: object) -> dt.datetime:
    s = str(ts_raw)
    if s.endswith("Z"):
        s = s[:-1] + "+00:00"
    return dt.datetime.fromisoformat(s).astimezone(_ET)


def _build_daily_bar(path: Path, date: str) -> dict | None:
    try:
        table = pq.read_table(str(path))
        rows = table.to_pydict()
    except Exception:
        return None
    opens, highs, lows, closes, vols = [], [], [], [], []
    for i, ts_raw in enumerate(rows.get("timestamp", [])):
        try:
            ts = _parse_ts(ts_raw)
        except Exception:
            continue
        if _MKT_OPEN <= ts.time() < _MKT_CLOSE:
            opens.append(float(rows["open"][i] or 0))
            highs.append(float(rows["high"][i] or 0))
            lows.append(float(rows["low"][i] or 0))
            closes.append(float(rows["close"][i] or 0))
            vols.append(float(rows["volume"][i] or 0))
    if not opens:
        return None
    return {
        "date": date, "open": opens[0], "high": max(highs),
        "low": min(lows), "close": closes[-1], "volume": sum(vols),
    }


def build_daily_bars(tickers: list[str], dates: list[str], workers: int = 8) -> dict[str, list[dict]]:
    root = Path(INTRADAY_CACHE_DIR)
    def _load(ticker: str) -> tuple[str, list[dict]]:
        d_path = root / ticker
        if not d_path.is_dir():
            return ticker, []
        bars: list[dict] = []
        for date in dates:
            p = d_path / f"{date}.parquet"
            if not p.exists():
                continue
            bar = _build_daily_bar(p, date)
            if bar and bar["close"] > 0:
                bars.append(bar)
        return ticker, bars
    result: dict[str, list[dict]] = {}
    with concurrent.futures.ThreadPoolExecutor(max_workers=workers) as ex:
        for ticker, bars in ex.map(_load, tickers):
            if bars:
                result[ticker] = bars
    return result


def load_intraday_bulk(candidates: dict[str, list[str]]) -> dict[str, dict[str, list[dict]]]:
    import pandas as pd
    result: dict[str, dict[str, list[dict]]] = {}
    for date, tickers in candidates.items():
        day_bars: dict[str, list[dict]] = {}
        for ticker in tickers:
            p = Path(INTRADAY_CACHE_DIR) / ticker / f"{date}.parquet"
            if not p.exists():
                continue
            try:
                df = pd.read_parquet(str(p))
                if not df.empty and len(df) >= 5:
                    day_bars[ticker] = df.to_dict("records")
            except Exception:
                pass
        if day_bars:
            result[date] = day_bars
    return result


def pearson(xs: list[float], ys: list[float]) -> float | None:
    n = len(xs)
    if n < 2:
        return None
    xm, ym = sum(xs) / n, sum(ys) / n
    num = sum((xs[i] - xm) * (ys[i] - ym) for i in range(n))
    dx = sum((x - xm) ** 2 for x in xs) ** 0.5
    dy = sum((y - ym) ** 2 for y in ys) ** 0.5
    if dx <= 0 or dy <= 0:
        return None
    return num / (dx * dy)


def tercile_stats(vals: list[float], rs: list[float]) -> dict:
    if len(vals) < 9:
        return {}
    pairs = sorted(zip(vals, rs), key=lambda p: p[0])
    n = len(pairs)
    t = n // 3
    def stats(sub):
        ys = [p[1] for p in sub]
        return {"n": len(ys), "wr": sum(1 for y in ys if y > 0) / len(ys), "avg_r": sum(ys) / len(ys)}
    return {"low": stats(pairs[:t]), "mid": stats(pairs[t:2*t]), "high": stats(pairs[2*t:])}


def report_feature(label: str, vals: list[float], rs: list[float], obv20: list[float]) -> tuple[bool, bool, bool]:
    n = len(vals)
    p = pearson(vals, rs)
    ts = tercile_stats(vals, rs)
    rho_obv20 = pearson(vals, obv20) if obv20 else None
    if not ts or p is None:
        print(f"  {label}: insufficient data")
        return False, False, False
    low, mid, high = ts["low"], ts["mid"], ts["high"]
    avg_r_gap = abs(high["avg_r"] - low["avg_r"])
    wr_gap = abs(high["wr"] - low["wr"])
    g1 = abs(p) >= 0.07 and n >= 120
    g2 = avg_r_gap >= 0.30
    g3 = wr_gap >= 0.05
    g5a = rho_obv20 is None or abs(rho_obv20) < 0.70
    print(f"\n  [{label}]  n={n}  Pearson={p:+.3f}")
    print(f"    G1: {'PASS' if g1 else 'FAIL'} (|{p:.3f}| {'≥' if abs(p)>=0.07 else '<'} 0.07, n={n})")
    print(f"    G2: {'PASS' if g2 else 'FAIL'} (avg_R gap = {avg_r_gap:.3f}R  [threshold 0.30R])")
    print(f"    G3: {'PASS' if g3 else 'FAIL'} (WR gap = {wr_gap*100:.1f}pp  [threshold 5pp])")
    print(f"    G5a: {'PASS' if g5a else 'FAIL'} (|ρ(feature, obv_slope_20)| = {abs(rho_obv20):.3f}  [threshold 0.70])" if rho_obv20 is not None else "    G5a: N/A")
    print(f"    Bottom tercile: WR {low['wr']*100:.1f}%  avg_R {low['avg_r']:+.3f}  n={low['n']}")
    print(f"    Middle tercile: WR {mid['wr']*100:.1f}%  avg_R {mid['avg_r']:+.3f}  n={mid['n']}")
    print(f"    Top tercile:    WR {high['wr']*100:.1f}%  avg_R {high['avg_r']:+.3f}  n={high['n']}")
    all_pass = g1 and g2 and g3 and g5a
    print(f"    → {'ALL GATES PASS ✓' if all_pass else 'FAIL'}")
    return g1, g2, g3


def main() -> None:
    print("=== V34 OBV Slope Acceleration Diagnostic ===\n")

    with open(V24_CONFIG) as f:
        raw = yaml.safe_load(f)
    params = ORBStrategyParams(**raw["orb_strategy"])

    today = dt.date(2026, 4, 21)
    all_td = trading_days_between(today - dt.timedelta(days=700), today)
    trading_days_list = [d.isoformat() for d in all_td[-LOOKBACK_DAYS:]]
    print(f"Window: {trading_days_list[0]} → {trading_days_list[-1]} ({LOOKBACK_DAYS} trading days)")

    extended_td = [d.isoformat() for d in all_td[-(LOOKBACK_DAYS + 35):]]
    first_cal = dt.date.fromisoformat(extended_td[0])
    last_cal = dt.date.fromisoformat(trading_days_list[-1])
    needed_dates: list[str] = []
    d = first_cal
    while d <= last_cal:
        needed_dates.append(d.isoformat())
        d += dt.timedelta(days=1)

    with open(UNIVERSE_FILE) as f:
        udata = yaml.safe_load(f)
    universe = udata.get("symbols", udata) if isinstance(udata, dict) else udata
    if "QQQ" not in universe:
        universe = list(universe) + ["QQQ"]
    print(f"Universe: {len(universe)} tickers")

    print(f"Building daily bars ({len(needed_dates)} calendar days)...")
    daily_bars = build_daily_bars(universe, needed_dates)
    print(f"Built: {len(daily_bars)} tickers")
    print("Computing enrichment...")
    enrichment = enrich_daily_bars(daily_bars, trading_days_list)

    candidates = orb_pre_screen_candidates(
        daily_bars, trading_days_list, enrichment,
        min_price=params.min_price, min_atr=params.min_atr_14,
        min_avg_dollar_vol=params.min_avg_dollar_volume, max_per_day=None,
    )
    print(f"Pre-screened: {sum(len(v) for v in candidates.values())} ticker-days")
    all_intraday = load_intraday_bulk(candidates)

    print("Running V24 simulation...")
    state = ORBSimulationState(equity=params.initial_capital)
    day_results, _ = run_orb_simulation_with_state(
        all_intraday, trading_days_list, params, enrichment, state=state,
    )
    all_trades = [t for dr in day_results for t in dr.trades]
    trades_with_r = [t for t in all_trades if getattr(t, "r_multiple_at_exit", None) is not None]
    print(f"Total trades: {len(all_trades)}, with r_multiple: {len(trades_with_r)}\n")

    sorted_daily: dict[str, list[dict]] = {
        t: sorted(bars, key=lambda b: b["date"]) for t, bars in daily_bars.items()
    }

    # Compute per-trade features
    slope5_vals, slope20_vals, accel_vals, r_mults = [], [], [], []
    missing5, missing20, missingaccel = 0, 0, 0

    for trade in trades_with_r:
        ticker = trade.ticker
        date = str(trade.date)[:10]
        r = float(trade.r_multiple_at_exit)

        bars_t = sorted_daily.get(ticker, [])
        prev_bars = [b for b in bars_t if b["date"][:10] < date]

        s5 = compute_obv_slope_approx(prev_bars, lookback=5) if len(prev_bars) >= 7 else None
        s20 = compute_obv_slope_approx(prev_bars, lookback=20) if len(prev_bars) >= 22 else None

        if s5 is None:
            missing5 += 1
        if s20 is None:
            missing20 += 1
        if s5 is None or s20 is None:
            missingaccel += 1
            continue

        slope5_vals.append(s5)
        slope20_vals.append(s20)
        accel_vals.append(s5 - s20)
        r_mults.append(r)

    n_valid = len(r_mults)
    print(f"Valid trades (both slope5 & slope20): {n_valid} / {len(trades_with_r)}")
    print(f"Missing slope5: {missing5}  slope20: {missing20}")

    print("\n" + "=" * 60)
    print("GATE RESULTS (G1: |P|≥0.07, n≥120; G2: avg_R≥0.30R; G3: WR≥5pp; G5a: ρ<0.70 vs obv_slope_20)")

    report_feature("obv_slope_5", slope5_vals, r_mults, slope20_vals)
    report_feature("obv_slope_accel (slope5 - slope20)", accel_vals, r_mults, slope20_vals)

    p_5_20 = pearson(slope5_vals, slope20_vals)
    print(f"\n  ρ(slope_5, slope_20) = {p_5_20:.3f}  (G5a threshold: 0.70)")

    print("\n=== Summary ===")
    if p_5_20 is not None and abs(p_5_20) >= 0.70:
        print("G5a: obv_slope_5 is REDUNDANT with obv_slope_20 (ρ ≥ 0.70). Cannot add independent signal.")
        print("Conclusion: V24's obv_slope_20 captures the full OBV axis. No improvement possible via window change.")
    else:
        print("G5a: obv_slope_5 has orthogonal component vs obv_slope_20. Check individual feature gates above.")


if __name__ == "__main__":
    main()