spb/anomaly-atlas Public License

1# =============================================================================2#  Project   : anomaly-atlas3#  File      : experiments/micro/expE_calendar_scan/benchmark.py4#  Purpose   : Calendar scan with pre-counted budget + permuted-calendar null5#  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 E — calendar scan (protocol + 8-test budget pre-specified in14hypothesis.md). Level 0 throughout. Also measures the H20 intraday artifact15profile (taxonomy input, not a hypothesis test)."""1617from __future__ import annotations1819import json20import sys21from datetime import UTC, date, datetime22from pathlib import Path2324import numpy as np2526REPO_ROOT = Path(__file__).resolve().parents[3]27sys.path.insert(0, str(REPO_ROOT / "benchmarks"))28sys.path.insert(0, str(REPO_ROOT / "src"))2930from hardware_manifest import collect_manifest  # noqa: E4023132from anomaly_atlas.data.hf_client import HFMarketDataClient  # noqa: E40233from anomaly_atlas.data.universe import LIQUID_ETF, LIQUID_STOCK, TRAIN  # noqa: E40234from anomaly_atlas.stats.multiple_testing import benjamini_hochberg  # noqa: E40235from anomaly_atlas.validation.artifacts import edge_spread  # noqa: E4023637N_PERM, PERM_SEED = 2000, 4238UMICRO_WINDOW = ("2014-01-01", "2016-01-01")  # H20 profile, cached from expC/D39HALF_HOURS = [f"{h:02d}:{m:02d}" for h, m in40              [(9, 30), (10, 0), (10, 30), (11, 0), (11, 30), (12, 0), (12, 30),41               (13, 0), (13, 30), (14, 0), (14, 30), (15, 0), (15, 30)]]424344def weekday(d: str) -> int:45    return date(int(d[:4]), int(d[5:7]), int(d[8:10])).weekday()464748def class_masks(dates: list[str]) -> dict[str, np.ndarray]:49    """The 8 pre-declared calendar classes as boolean masks over days."""50    n = len(dates)51    wd = np.array([weekday(d) for d in dates])52    masks = {name: wd == i for i, name in enumerate(["mon", "tue", "wed", "thu", "fri"])}5354    month = np.array([d[:7] for d in dates])55    tom = np.zeros(n, dtype=bool)56    for i in range(n):57        if i + 1 < n and month[i + 1] != month[i]:58            tom[i] = True                      # last trading day of month (-1)59        if i > 0 and month[i - 1] != month[i]:60            for j in range(i, min(i + 3, n)):  # first 3 trading days (+1..+3)61                if month[j] == month[i]:62                    tom[j] = True63    masks["turn_of_month"] = tom6465    # holidays: a non-weekend gap in the trading calendar66    pre = np.zeros(n, dtype=bool)67    post = np.zeros(n, dtype=bool)68    for i in range(n - 1):69        d0 = date(int(dates[i][:4]), int(dates[i][5:7]), int(dates[i][8:10]))70        d1 = date(int(dates[i + 1][:4]), int(dates[i + 1][5:7]), int(dates[i + 1][8:10]))71        gap_weekdays = np.busday_count(d0.isoformat(), d1.isoformat()) - 172        if gap_weekdays >= 1:73            pre[i] = True74            post[i + 1] = True75    masks["pre_holiday"] = pre76    masks["post_holiday"] = post77    return masks787980def stats_for(returns: np.ndarray, masks: dict[str, np.ndarray]) -> dict[str, float]:81    mu = returns.mean()82    return {name: float(returns[m].mean() - mu) if m.sum() >= 20 else float("nan")83            for name, m in masks.items()}848586def permutation_test(returns: np.ndarray, years: np.ndarray,87                     masks: dict[str, np.ndarray]) -> tuple[dict, dict]:88    """Within-year permutation: marginal p per class + family-wise max-stat p."""89    rng = np.random.default_rng(PERM_SEED)90    observed = stats_for(returns, masks)91    names = list(masks)92    exceed = dict.fromkeys(names, 0)93    fw_exceed = dict.fromkeys(names, 0)94    perm_dist: dict[str, list[float]] = {k: [] for k in names}95    year_idx = [np.where(years == y)[0] for y in np.unique(years)]96    for _ in range(N_PERM):97        perm = returns.copy()98        for idx in year_idx:99            perm[idx] = perm[idx][rng.permutation(len(idx))]100        s = stats_for(perm, masks)101        max_abs = max(abs(v) for v in s.values() if np.isfinite(v))102        for name in names:103            perm_dist[name].append(s[name])104            if np.isfinite(s[name]) and abs(s[name]) >= abs(observed[name]):105                exceed[name] += 1106            if np.isfinite(observed[name]) and max_abs >= abs(observed[name]):107                fw_exceed[name] += 1108    marg = {k: (exceed[k] + 1) / (N_PERM + 1) for k in names}109    fw = {k: (fw_exceed[k] + 1) / (N_PERM + 1) for k in names}110    band = {k: [round(float(np.percentile(perm_dist[k], q)) * 1e4, 3) for q in (2.5, 97.5)]111            for k in names}112    return {"observed_bp": {k: round(v * 1e4, 3) for k, v in observed.items()},113            "perm_band95_bp": band,114            "p_marginal": {k: round(v, 4) for k, v in marg.items()},115            "p_familywise": {k: round(v, 4) for k, v in fw.items()}}, observed116117118def h20_profile(client: HFMarketDataClient) -> dict:119    """Intraday half-hour profile of |return|, EDGE spread, staleness."""120    tickers = [("stock", t) for t in LIQUID_STOCK] + [("etf", t) for t in LIQUID_ETF]121    buckets = {hh: {"absret": [], "bars": [], "minutes": 0, "present": 0} for hh in HALF_HOURS}122    s, e = UMICRO_WINDOW123    for asset, ticker in tickers:124        bars = client.get_bars(asset, ticker, "1min", "adj_split", s, e)125        by_day: dict[str, list[dict]] = {}126        for b in bars:127            t = b["datetime"][11:16]128            if "09:30" <= t < "16:00":129                by_day.setdefault(b["datetime"][:10], []).append(b)130        n_days = len(by_day)131        for hh_i, hh in enumerate(HALF_HOURS):132            hi = HALF_HOURS[hh_i + 1] if hh_i + 1 < len(HALF_HOURS) else "16:00"133            sel = [b for bs in by_day.values() for b in bs if hh <= b["datetime"][11:16] < hi]134            closes = np.array([b["close"] for b in sel])135            if len(closes) > 100:136                r = np.abs(np.diff(np.log(closes)))137                buckets[hh]["absret"].append(float(np.median(r)))138                buckets[hh]["bars"].append(sel)139            buckets[hh]["minutes"] += 30 * n_days140            buckets[hh]["present"] += len(sel)141    profile = []142    for hh in HALF_HOURS:143        b = buckets[hh]144        spreads = []145        for sel in b["bars"]:146            o = np.array([x["open"] for x in sel])147            h = np.array([x["high"] for x in sel])148            lo = np.array([x["low"] for x in sel])149            c = np.array([x["close"] for x in sel])150            sp = edge_spread(o, h, lo, c)151            if np.isfinite(sp):152                spreads.append(sp)153        profile.append({154            "bucket": hh,155            "median_abs_1min_ret_bp": round(float(np.median(b["absret"])) * 1e4, 3)156            if b["absret"] else None,157            "median_edge_spread_bp": round(float(np.median(spreads)) * 1e4, 3)158            if spreads else None,159            "staleness": round(1 - b["present"] / b["minutes"], 4) if b["minutes"] else None,160        })161    return {"window": UMICRO_WINDOW, "universe": LIQUID_STOCK + LIQUID_ETF,162            "profile": profile}163164165def main() -> None:166    run_utc = datetime.now(UTC)167    client = HFMarketDataClient()168169    bars = client.get_bars("etf", "SPY", "1day", "adj_splitdiv", TRAIN[0], TRAIN[1])170    dates = [b["datetime"][:10] for b in bars]171    closes = np.array([b["close"] for b in bars])172    returns = np.diff(np.log(closes))173    dates = dates[1:]  # return dates174    years = np.array([d[:4] for d in dates])175    masks = class_masks(dates)176177    table, observed = permutation_test(returns, years, masks)178    fdr = benjamini_hochberg(np.array([table["p_marginal"][k] for k in masks]), alpha=0.05)179    table["fdr_marginal"] = {k: bool(r) for k, r in zip(masks, fdr, strict=True)}180181    # sub-period stability, descriptive182    subs = {}183    for name, (lo, hi) in {"2000-2007": ("2000", "2007"), "2008-2015": ("2008", "2015")}.items():184        sel = (years >= lo) & (years <= hi)185        subs[name] = {k: round(v * 1e4, 3)186                      for k, v in stats_for(returns[sel],187                                            {k: m[sel] for k, m in masks.items()}).items()}188189    results = {190        "experiment": "expE_calendar_scan",191        "run_utc": run_utc.isoformat(),192        "author": "Simon-Pierre Boucher",193        "contact": "contact@spboucher.ai",194        "data_source": "hfmarketdata.io",195        "confidence_level": 0,196        "protocol": {"instrument": "SPY 1day adj_splitdiv", "train": TRAIN,197                     "budget_tests": 8, "n_perm": N_PERM, "perm_seed": PERM_SEED,198                     "class_counts": {k: int(m.sum()) for k, m in masks.items()},199                     "n_days": int(len(returns))},200        "calendar_tests": table,201        "subperiods_bp": subs,202        "h20_intraday_profile": h20_profile(client),203        "client_stats": vars(client.stats) | {"refreshes": list(client.stats.refreshes)},204        "manifest": collect_manifest(),205    }206    out_dir = REPO_ROOT / "results" / "expE_calendar_scan" / run_utc.strftime("%Y%m%dT%H%M%SZ")207    out_dir.mkdir(parents=True)208    (out_dir / "results.json").write_text(json.dumps(results, indent=2) + "\n")209    print(f"wrote {out_dir.relative_to(REPO_ROOT)}/results.json")210    print(json.dumps({k: results["calendar_tests"][k] for k in211                      ("observed_bp", "p_marginal", "p_familywise", "fdr_marginal")}, indent=1))212    print("H20 profile:", json.dumps(results["h20_intraday_profile"]["profile"], indent=1)[:800])213214215if __name__ == "__main__":216    main()217