spb/anomaly-atlas Public License

1# =============================================================================2#  Project   : anomaly-atlas3#  File      : experiments/micro/expB_artifact_baselines/benchmark.py4#  Purpose   : Measure the artifact nulls: bounce, staleness, LOCF lead-lag5#  Author    : Simon-Pierre Boucher6#  Contact   : contact@spboucher.ai7#  Data src  : hfmarketdata.io (sole data source)8#  Created   : 2026-08-129#  Modified  : 2026-08-1210#  Platform  : macOS / Apple Silicon (arm64)11#  License   : All rights reserved (research code)12# =============================================================================13"""Experiment B — artifact baselines on real data (pre-specified protocol in14hypothesis.md; detectors gated on synthetic ground truth first, §8.1).1516Every number produced here is a NULL LEVEL (Level 0 by construction): the17fake-signal magnitude that later experiments must exceed before claiming18anything. Universe, window, seeds are pre-specified; all data flows through19the cached hf_client.20"""2122from __future__ import annotations2324import json25import sys26from datetime import UTC, datetime27from pathlib import Path2829import numpy as np3031REPO_ROOT = Path(__file__).resolve().parents[3]32sys.path.insert(0, str(REPO_ROOT / "benchmarks"))33sys.path.insert(0, str(REPO_ROOT / "src"))3435from hardware_manifest import collect_manifest  # noqa: E4023637from anomaly_atlas.data.hf_client import HFMarketDataClient  # noqa: E40238from anomaly_atlas.stats.bootstrap import moving_block_bootstrap, percentile_ci  # noqa: E40239from anomaly_atlas.stats.reversion import ac1, variance_ratio  # noqa: E40240from anomaly_atlas.validation.artifacts import roll_spread  # noqa: E4024142LIQUID_STOCK = ["AAPL", "MSFT", "NVDA", "AMZN", "GOOGL", "META", "TSLA", "JPM", "XOM", "UNH"]43LIQUID_ETF = ["SPY", "QQQ"]44N_RANDOM, RANDOM_SEED = 30, 4245START, END = "2024-01-02", "2024-04-01"46ADJ = "adj_split"47BOOT_N, BOOT_SEED = 300, 4248MAX_LAG = 34950RTH_MINUTES = [f"{h:02d}:{m:02d}" for h in range(9, 16) for m in range(60)]51RTH_MINUTES = [t for t in RTH_MINUTES if "09:30" <= t < "16:00"]  # 390 slots52SLOT = {t: i for i, t in enumerate(RTH_MINUTES)}535455def rth_day_grids(bars: list[dict]) -> dict[str, np.ndarray]:56    """day -> 390-slot array of log close prices (NaN where no print)."""57    days: dict[str, np.ndarray] = {}58    for b in bars:59        dt = b["datetime"]60        t = dt[11:16]61        if not ("09:30" <= t < "16:00"):62            continue63        grid = days.setdefault(dt[:10], np.full(390, np.nan))64        grid[SLOT[t]] = np.log(b["close"])65    return days666768def trade_time_returns(days: dict[str, np.ndarray]) -> np.ndarray:69    """Within-day log returns between consecutive PRINTS (no grid, no LOCF)."""70    out = []71    for day in sorted(days):72        p = days[day]73        obs = p[np.isfinite(p)]74        if len(obs) >= 2:75            out.append(np.diff(obs))76    return np.concatenate(out) if out else np.array([])777879def locf_grid_returns(days: dict[str, np.ndarray], day_list: list[str]) -> np.ndarray:80    """Concatenated per-day LOCF grid returns, NaN before first print and at81    day boundaries — the join that MANUFACTURES the stale-price artifact."""82    out = []83    for day in day_list:84        p = days.get(day)85        if p is None:86            out.append(np.full(389, np.nan))87            continue88        filled = p.copy()89        for i in range(1, 390):90            if not np.isfinite(filled[i]):91                filled[i] = filled[i - 1]92        out.append(np.diff(filled))  # NaN propagates before first print93    return np.concatenate(out)949596def nan_xcorr(x: np.ndarray, y: np.ndarray, max_lag: int) -> dict[int, float]:97    """corr(x_{t-k}, y_t) over finite pairs only; k>0 = x leads y."""98    n = min(len(x), len(y))99    x, y = x[:n], y[:n]100    out: dict[int, float] = {}101    for k in range(-max_lag, max_lag + 1):102        a = x[: n - k] if k >= 0 else x[-k:]103        b = y[k:] if k >= 0 else y[: n + k]104        m = np.isfinite(a) & np.isfinite(b)105        if m.sum() < 100 or a[m].std() == 0 or b[m].std() == 0:106            out[k] = float("nan")107            continue108        out[k] = float(np.corrcoef(a[m], b[m])[0, 1])109    return out110111112def analyze_ticker(days: dict[str, np.ndarray], day_list: list[str]) -> dict | None:113    present = (114        np.concatenate([np.isfinite(days[d]) for d in day_list if d in days])115        if any(d in days for d in day_list)116        else np.array([])117    )118    n_days_covered = sum(d in days for d in day_list)119    if n_days_covered < 30:120        return None121    r = trade_time_returns(days)122    if len(r) < 2_000:123        return None124    block = max(50, len(r) // max(n_days_covered, 1))125    boot = moving_block_bootstrap(r, ac1, block=block, n_boot=BOOT_N, seed=BOOT_SEED)126    lo, hi = percentile_ci(boot)127    spread = roll_spread(r)128    return {129        "days_covered": n_days_covered,130        "staleness": round(1.0 - present.mean() * len(present) / (390 * n_days_covered), 4)131        if n_days_covered132        else None,133        "rth_fill_ratio": round(present.sum() / (390 * n_days_covered), 4),134        "n_trade_returns": int(len(r)),135        "ac1": round(ac1(r), 5),136        "ac1_ci95": [round(lo, 5), round(hi, 5)],137        "roll_rel_spread": round(spread, 6) if np.isfinite(spread) else None,138        "vr5": round(variance_ratio(r, 5), 4),139        "vr30": round(variance_ratio(r, 30), 4),140    }141142143def main() -> None:144    run_utc = datetime.now(UTC)145    client = HFMarketDataClient()146147    # deterministic random universe (seed pre-specified)148    all_stock = client.tickers("stock", timeframe="1min", adjustment=ADJ)149    rng = np.random.default_rng(RANDOM_SEED)150    random_universe = sorted(rng.choice(sorted(all_stock), N_RANDOM, replace=False))151    universe = (152        [("stock", t, "liquid") for t in LIQUID_STOCK]153        + [("etf", t, "liquid") for t in LIQUID_ETF]154        + [("stock", t, "random") for t in random_universe]155    )156157    # fetch + grid everything158    grids: dict[str, dict[str, np.ndarray]] = {}159    for asset, ticker, _ in universe:160        bars = client.get_bars(asset, ticker, "1min", ADJ, START, END)161        grids[ticker] = rth_day_grids(bars)162        print(f"{ticker}: {sum(len(v[np.isfinite(v)]) for v in grids[ticker].values())} RTH bars")163    day_list = sorted(grids["SPY"].keys())  # trading calendar := SPY days164165    # B1-B3: per-ticker artifact levels166    per_ticker: dict[str, dict] = {}167    for asset, ticker, bucket in universe:168        m = analyze_ticker(grids[ticker], day_list)169        if m is not None:170            m["bucket"] = bucket171            m["asset"] = asset172            per_ticker[ticker] = m173174    # B4: LOCF lead-lag vs SPY175    spy_r = locf_grid_returns(grids["SPY"], day_list)176    for ticker, m in per_ticker.items():177        if ticker == "SPY":178            continue179        r = locf_grid_returns(grids[ticker], day_list)180        xc = nan_xcorr(spy_r, r, MAX_LAG)181        m["xcorr_vs_spy"] = {str(k): round(v, 5) if np.isfinite(v) else None for k, v in xc.items()}182        m["spy_leads_+1"] = round(xc[1], 5) if np.isfinite(xc[1]) else None183184    # SPX (index) vs SPY — the non-synchronous-session case185    spx_bars = client.get_bars("index", "SPX", "1min", None, START, END)186    spx_grid = rth_day_grids(spx_bars)187    spx_r = locf_grid_returns(spx_grid, day_list)188    spx_xc = nan_xcorr(spx_r, spy_r, MAX_LAG)189190    # staleness -> artifact monotonicity (Spearman)191    pairs = [192        (m["staleness"], m["spy_leads_+1"])193        for m in per_ticker.values()194        if m.get("spy_leads_+1") is not None and m["staleness"] is not None195    ]196    xs = np.array([p[0] for p in pairs])197    ys = np.array([p[1] for p in pairs])198    rx = np.argsort(np.argsort(xs)).astype(float)199    ry = np.argsort(np.argsort(ys)).astype(float)200    spearman = float(np.corrcoef(rx, ry)[0, 1]) if len(pairs) > 5 else float("nan")201202    # aggregates by staleness tercile203    stale_vals = sorted(m["staleness"] for m in per_ticker.values())204    t1, t2 = np.percentile(stale_vals, [33.3, 66.7])205206    def tercile(s: float) -> str:207        return "fresh" if s <= t1 else "mid" if s <= t2 else "stale"208209    agg: dict[str, dict] = {}210    for name in ("fresh", "mid", "stale"):211        rows = [m for m in per_ticker.values() if tercile(m["staleness"]) == name]212        if not rows:213            continue214        agg[name] = {215            "n": len(rows),216            "median_staleness": round(float(np.median([m["staleness"] for m in rows])), 4),217            "median_ac1": round(float(np.median([m["ac1"] for m in rows])), 5),218            "median_roll_spread": round(219                float(np.median([m["roll_rel_spread"] for m in rows if m["roll_rel_spread"]])), 6220            ),221            "median_vr5": round(float(np.median([m["vr5"] for m in rows])), 4),222            "median_vr30": round(float(np.median([m["vr30"] for m in rows])), 4),223            "median_spy_leads_+1": round(224                float(225                    np.median(226                        [m["spy_leads_+1"] for m in rows if m.get("spy_leads_+1") is not None]227                    )228                ),229                5,230            ),231        }232233    results = {234        "experiment": "expB_artifact_baselines",235        "run_utc": run_utc.isoformat(),236        "author": "Simon-Pierre Boucher",237        "contact": "contact@spboucher.ai",238        "data_source": "hfmarketdata.io",239        "protocol": {240            "window": [START, END],241            "adjustment": ADJ,242            "rth": "09:30-16:00",243            "liquid": LIQUID_STOCK + LIQUID_ETF,244            "random_universe": list(random_universe),245            "random_seed": RANDOM_SEED,246            "boot": [BOOT_N, BOOT_SEED],247            "confidence_level": 0,248            "note": "artifact NULL levels — descriptive, in-sample by design",249        },250        "per_ticker": per_ticker,251        "terciles": {"cuts": [round(float(t1), 4), round(float(t2), 4)], "agg": agg},252        "staleness_vs_spy_lead_spearman": round(spearman, 4),253        "spx_vs_spy_xcorr": {254            str(k): round(v, 5) if np.isfinite(v) else None for k, v in spx_xc.items()255        },256        "client_stats": vars(client.stats) | {"refreshes": list(client.stats.refreshes)},257        "manifest": collect_manifest(),258    }259260    out_dir = REPO_ROOT / "results" / "expB_artifact_baselines" / run_utc.strftime("%Y%m%dT%H%M%SZ")261    out_dir.mkdir(parents=True)262    (out_dir / "results.json").write_text(json.dumps(results, indent=2) + "\n")263    print(f"\nwrote {out_dir.relative_to(REPO_ROOT)}/results.json")264    print("terciles:", json.dumps(agg, indent=1))265    print("spearman(staleness, SPY leads +1):", round(spearman, 4))266    print("SPX vs SPY xcorr:", results["spx_vs_spy_xcorr"])267268269if __name__ == "__main__":270    main()271