PythagorasGoal/scripts/research/xsec/xslib.py

"""xslib — SHARED CROSS-SECTIONAL research harness over the certified Hyperliquid alt panel.

Built for the "cerca altre strategie" wave (2026-06-20, follow-up to the 104-hypothesis BTC/ETH
sweep that exhausted the single-asset directional space). The frontier the prior synthesis pointed
to: CROSS-SECTIONAL / multi-asset mechanisms on the 51 certified Hyperliquid alts (1d, 2024-2026),
where the ~1.3 BTC/ETH-directional ceiling does NOT bind, and DISTINCT from XS01 (plain x-sec momentum).

Why a new harness: the panel is N assets × ~900 days. A strategy = a per-asset SCORE computed
causally (data <= close[i]); the harness ranks it cross-sectionally each rebalance, goes long the
top-k / short the bottom-k (market-neutral) or long-only top-k, vol-targets, charges fee on turnover,
and — crucially — the weight decided at bar i is applied to the return of bar i+1, so look-ahead is
structurally impossible (same convention as src.portfolio xs_book / sleeves._xsec_returns).

A candidate only matters if it (a) is robust (positive FULL + hold-out 2025+ + jackknife), AND
(b) is DISTINCT from XS01 (low correlation), AND (c) ADDS to the live TP01+XS01+VRP01 portfolio.
CAVEAT baked in: the panel is ~2.5 years — every result is SUGGESTIVE, not robust like 6y BTC/ETH.

Quick start (agent script):
    import sys; sys.path.insert(0, "/opt/docker/PythagorasGoal/scripts/research/xsec")
    import xslib as xs, numpy as np
    p = xs.load_panel("all")                       # or "majors", a list, or an int N (top-N liquidity)
    score = xs.past_return(p.close, 30)            # momentum: higher = long
    rep = xs.study_xs("MOM30", lambda P: xs.past_return(P.close, 30), H=10, k=5)
    print(xs.fmt(rep)); print("JSON:", xs.as_json(rep))
"""
from __future__ import annotations

import glob
import json
import sys
import warnings
from dataclasses import dataclass
from functools import lru_cache
from pathlib import Path

import numpy as np
import pandas as pd

# panel research has many all-NaN edge windows (rolling beta/vol on first rows) -> benign
warnings.filterwarnings("ignore", category=RuntimeWarning)

_ROOT = Path(__file__).resolve().parents[3]
if str(_ROOT) not in sys.path:
    sys.path.insert(0, str(_ROOT))
sys.path.insert(0, str(_ROOT / "scripts" / "research" / "alt"))
import altlib as al  # noqa: E402  (reuse _sh, _dd_ret, _to_daily, HOLDOUT, metric helpers)

RAW = _ROOT / "data" / "raw"
HOLDOUT = pd.Timestamp("2025-01-01", tz="UTC")
FEE = 0.001  # round-trip; charged /2 per side on turnover

MAJORS = ["BTC", "ETH", "SOL", "BNB", "XRP", "DOGE", "AVAX", "LINK", "LTC", "ADA",
          "ARB", "OP", "SUI", "APT", "INJ", "TIA", "SEI", "NEAR", "AAVE"]


# ===========================================================================
# PANEL
# ===========================================================================
@dataclass
class Panel:
    syms: list
    index: pd.DatetimeIndex
    close: np.ndarray
    open: np.ndarray
    high: np.ndarray
    low: np.ndarray
    vol: np.ndarray
    ret: np.ndarray   # daily simple returns, ret[0]=0


@lru_cache(maxsize=16)
def load_panel(universe="all", min_rows: int = 700) -> Panel:
    """Common-date OHLCV panel of the certified HL alts (1d). `universe`:
       'all'  -> every alt with >= min_rows of history (drops short ones e.g. ALGO/SAND),
       'majors' -> the 19 XS01 majors, a list of symbols, or an int N (top-N by median $-volume)."""
    close, vol, high, low, opn = {}, {}, {}, {}, {}
    for f in sorted(glob.glob(str(RAW / "hl_*_1d.parquet"))):
        sym = Path(f).stem.replace("hl_", "").replace("_1d", "").upper()
        d = pd.read_parquet(f)
        if len(d) < min_rows:
            continue
        idx = pd.to_datetime(d["timestamp"], unit="ms", utc=True)
        close[sym] = pd.Series(d["close"].values.astype(float), index=idx)
        vol[sym] = pd.Series(d["volume"].values.astype(float), index=idx)
        high[sym] = pd.Series(d["high"].values.astype(float), index=idx)
        low[sym] = pd.Series(d["low"].values.astype(float), index=idx)
        opn[sym] = pd.Series(d["open"].values.astype(float), index=idx)
    C = pd.concat(close, axis=1, join="inner").sort_index().dropna()
    syms = list(C.columns)
    if universe == "majors":
        syms = [s for s in MAJORS if s in syms]
    elif isinstance(universe, (list, tuple)):
        syms = [s for s in universe if s in syms]
    elif isinstance(universe, int):
        dollar = {s: float(np.nanmedian(C[s].values * pd.concat(vol, axis=1)[s].reindex(C.index).values))
                  for s in syms}
        syms = sorted(syms, key=lambda s: -dollar[s])[:universe]
    C = C[syms]
    idx = C.index

    def stack(dd):
        return pd.concat(dd, axis=1).reindex(index=idx)[syms].values.astype(float)
    cl = C.values
    ret = np.zeros_like(cl)
    ret[1:] = cl[1:] / cl[:-1] - 1.0
    return Panel(syms, idx, cl, stack(opn), stack(high), stack(low), stack(vol), ret)


# ===========================================================================
# CAUSAL CROSS-SECTIONAL HELPERS  (value at row i uses data <= i)
# ===========================================================================
def past_return(close, L):
    out = np.full_like(close, np.nan)
    out[L:] = close[L:] / close[:-L] - 1.0
    return out


def roll_std(mat, win):
    return pd.DataFrame(mat).rolling(win, min_periods=max(2, win // 2)).std().values


def roll_mean(mat, win):
    return pd.DataFrame(mat).rolling(win, min_periods=max(2, win // 2)).mean().values


def roll_skew(mat, win):
    return pd.DataFrame(mat).rolling(win, min_periods=max(3, win // 2)).skew().values


def ewm_mean(mat, span):
    return pd.DataFrame(mat).ewm(span=span, adjust=False).mean().values


def xs_zscore(mat):
    """Cross-sectional z-score per row (across assets). NaN-safe."""
    m = np.nanmean(mat, axis=1, keepdims=True)
    s = np.nanstd(mat, axis=1, keepdims=True)
    return (mat - m) / np.where(s > 0, s, np.nan)


def xs_rank(mat):
    """Cross-sectional rank in [0,1] per row (0=lowest)."""
    out = np.full_like(mat, np.nan, dtype=float)
    for i in range(mat.shape[0]):
        row = mat[i]
        ok = np.isfinite(row)
        if ok.sum() >= 2:
            r = pd.Series(row[ok]).rank().values
            out[i, ok] = (r - 1) / (ok.sum() - 1)
    return out


def market_ret(ret):
    """Equal-weight market return per day (n,)."""
    return np.nanmean(ret, axis=1)


def roll_beta(ret, win):
    """Rolling beta of each asset to the equal-weight market (n,A), causal."""
    mkt = market_ret(ret)
    ms = pd.Series(mkt)
    var = ms.rolling(win, min_periods=max(5, win // 2)).var()
    out = np.full_like(ret, np.nan)
    for a in range(ret.shape[1]):
        cov = pd.Series(ret[:, a]).rolling(win, min_periods=max(5, win // 2)).cov(ms)
        out[:, a] = (cov / var.replace(0, np.nan)).values
    return out


def residual_return(ret, win):
    """Idiosyncratic daily return = ret - beta*market (beta rolling, causal)."""
    beta = roll_beta(ret, win)
    mkt = market_ret(ret)[:, None]
    return ret - beta * mkt


def volume_z(vol, win):
    m = roll_mean(vol, win)
    s = roll_std(vol, win)
    return (vol - m) / np.where(s > 0, s, np.nan)


# ===========================================================================
# BACKTEST — generic cross-sectional book from a per-asset SCORE matrix.
# score[i] (data <= i) -> rank assets -> long top-k / short bottom-k; W[i] earns dret[i+1].
# ===========================================================================
def xs_backtest(panel: Panel, score, H=10, k=5, long_short=True, target_vol=0.20,
                fee=FEE, vt_cap=3.0):
    px = panel.close
    n, A = px.shape
    dret = panel.ret
    score = np.asarray(score, float)
    if score.shape != (n, A):
        raise ValueError(f"score shape {score.shape} != panel {(n, A)}")
    W = np.zeros((n, A))
    w = np.zeros(A)
    for i in range(n):
        if i % H == 0:
            row = score[i]
            fin = np.isfinite(row)
            if fin.sum() >= 2 * k:
                ranked = np.where(fin, row, -np.inf)
                order = np.argsort(ranked)
                order = order[np.isfinite(ranked[order])]
                lo, hi = order[:k], order[-k:]
                w = np.zeros(A)
                if long_short:
                    w[hi] = 0.5 / k
                    w[lo] = -0.5 / k
                else:
                    w[hi] = 1.0 / k
        W[i] = w
    gross = np.zeros(n)
    gross[1:] = np.sum(W[:-1] * dret[1:], axis=1)
    turn = np.zeros(n)
    turn[0] = np.abs(W[0]).sum()
    turn[1:] = np.abs(np.diff(W, axis=0)).sum(axis=1)
    net = gross - turn * (fee / 2.0)
    s = pd.Series(net, index=panel.index)
    rv = s.rolling(30, min_periods=15).std().shift(1) * np.sqrt(365.25)
    scale = np.clip(np.nan_to_num(target_vol / rv.replace(0, np.nan).values, nan=0.0), 0, vt_cap)
    return pd.Series(s.values * scale, index=panel.index)


# ===========================================================================
# BASELINES (live stack) + MARGINAL scoring
# ===========================================================================
@lru_cache(maxsize=1)
def baselines():
    """Daily returns of the LIVE stack: TP01, XS01, and the combined active portfolio."""
    from src.portfolio.portfolio import StrategyPortfolio, to_daily
    from src.portfolio.sleeves import _tp01_returns, _xsec_returns, active_sleeves
    tp = to_daily(_tp01_returns())
    xs01 = to_daily(_xsec_returns())
    active = StrategyPortfolio(active_sleeves()).combined_daily()
    return dict(tp01=tp, xs01=xs01, active=active)


def _corr(a, b):
    J = pd.concat({"a": a, "b": b}, axis=1, join="inner").dropna()
    return round(float(J["a"].corr(J["b"])), 3) if len(J) > 5 else None


def marginal_vs(cand, base, weights=(0.2, 0.35)):
    """Does `cand` improve `base`? blend uplift (full & hold-out), + OOS jackknife robustness."""
    J = pd.concat({"B": base, "C": cand}, axis=1, join="inner").dropna()
    if len(J) < 30:
        return dict(verdict="N/A", reason="overlap < 30d")
    JH = J[J.index >= HOLDOUT]
    has_h = len(JH) > 20
    out = dict(corr=_corr(J["B"], J["C"]), base_full=round(al._sh(J["B"]), 3),
               base_hold=round(al._sh(JH["B"]), 3) if has_h else None,
               cand_full=round(al._sh(J["C"]), 3), cand_hold=round(al._sh(JH["C"]), 3) if has_h else None,
               blends={})
    for w in weights:
        bf, bh = (1 - w) * J["B"] + w * J["C"], (1 - w) * JH["B"] + w * JH["C"]
        out["blends"][f"w{int(w * 100)}"] = dict(
            uplift_full=round(al._sh(bf) - al._sh(J["B"]), 3),
            uplift_hold=round(al._sh(bh) - al._sh(JH["B"]), 3) if has_h else None,
            dd=round(al._dd_ret(bf), 4))
    # OOS jackknife at w=0.2
    robust = False
    cu = jk = None
    if has_h:
        def _u(sub):
            return al._sh(0.8 * sub["B"] + 0.2 * sub["C"]) - al._sh(sub["B"])
        months = sorted(set(zip(JH.index.year, JH.index.month)))
        cu = round(_u(JH), 3)
        jk = round(min(_u(JH[~((JH.index.year == y) & (JH.index.month == m))]) for y, m in months), 3) \
            if len(months) > 1 else cu
        robust = bool(cu > 0.02 and jk > 0.0)
    out["holdout_uplift_w20"] = cu
    out["jackknife_min_uplift"] = jk
    out["robust_oos"] = robust
    up = out["blends"][f"w{int(weights[0] * 100)}"]["uplift_hold"]
    cc = out["corr"] if out["corr"] is not None else 0.0
    if cc is not None and cc > 0.85 and (up is None or abs(up) < 0.05):
        out["verdict"] = "REDUNDANT"
    elif up is not None and up >= 0.05 and robust:
        out["verdict"] = "ADDS"
    elif up is not None and up <= -0.05:
        out["verdict"] = "DILUTES"
    else:
        out["verdict"] = "NEUTRAL"
    return out


# ===========================================================================
# DRIVER
# ===========================================================================
def study_xs(name, score_fn, universe="all", H=10, k=5, long_short=True,
             target_vol=0.20, min_rows=700) -> dict:
    """Backtest one cross-sectional hypothesis and score it honestly:
       FULL + hold-out 2025+ + yearly, correlation to TP01 & XS01 (distinctness),
       and marginal contribution to the LIVE active portfolio. `score_fn(panel) -> (n,A)`
       per-asset score (higher = long), computed CAUSALLY (data <= close[i])."""
    p = load_panel(universe, min_rows=min_rows)
    score = score_fn(p)
    daily = al._to_daily(xs_backtest(p, score, H=H, k=k, long_short=long_short, target_vol=target_vol))
    net = daily.values
    idx = daily.index
    full = al._metrics_from_net(net, idx)
    hmask = idx >= HOLDOUT
    hold = al._metrics_from_net(net[hmask], idx[hmask]) if hmask.sum() > 20 else dict(sharpe=0.0, n=int(hmask.sum()))
    bl = baselines()
    marg = marginal_vs(daily, bl["active"])
    earns_slot = (full["sharpe"] > 0 and hold.get("sharpe", 0) > 0
                  and marg.get("verdict") == "ADDS"
                  and (_corr(daily, bl["xs01"]) or 0) < 0.6)  # distinct from existing x-sec
    return dict(
        name=name, universe=str(universe), H=H, k=k, long_short=long_short,
        n_assets=len(p.syms), n_days=int(len(idx)),
        full=full, holdout=hold, yearly=al._yearly(net, idx),
        corr_tp01=_corr(daily, bl["tp01"]), corr_xs01=_corr(daily, bl["xs01"]),
        corr_active=_corr(daily, bl["active"]),
        marginal=marg, earns_slot=earns_slot,
        caveat="panel ~2.5y (2024-26): suggestive, not robust",
    )


def _clean(o):
    if isinstance(o, dict):
        return {k: _clean(v) for k, v in o.items()}
    if isinstance(o, (list, tuple)):
        return [_clean(x) for x in o]
    if isinstance(o, (np.floating,)):
        return round(float(o), 4)
    if isinstance(o, (np.integer,)):
        return int(o)
    if isinstance(o, (np.bool_,)):
        return bool(o)
    return o


def as_json(rep):
    return json.dumps(_clean(rep), default=str)


def fmt(rep):
    m = rep["marginal"]
    yr = " ".join(f"{y}:{d['ret'] * 100:+.0f}%" for y, d in rep["yearly"].items())
    return (f"=== {rep['name']} [{rep['universe']} H{rep['H']} k{rep['k']} "
            f"{'LS' if rep['long_short'] else 'LO'}] EARNS_SLOT={rep['earns_slot']}\n"
            f"  FULL Sh {rep['full']['sharpe']:+.2f} DD {rep['full']['maxdd'] * 100:.0f}% "
            f"ret {rep['full']['ret'] * 100:+.0f}% | HOLD Sh {rep['holdout'].get('sharpe', 0):+.2f} "
            f"| corr TP01 {rep['corr_tp01']} XS01 {rep['corr_xs01']}\n"
            f"  marginal vs active: {m.get('verdict')} (corr {m.get('corr')}, "
            f"holdUplift_w20 {m.get('holdout_uplift_w20')}, jackknife {m.get('jackknife_min_uplift')}, "
            f"robust_oos {m.get('robust_oos')}) | {yr}")


if __name__ == "__main__":
    print("--- SMOKE TEST xslib ---")
    # 1) x-sec momentum (should resemble XS01 ballpark) ; 2) short-term reversal ; 3) low-vol
    print(fmt(study_xs("MOM30-90", lambda P: xs_zscore(past_return(P.close, 30)) + xs_zscore(past_return(P.close, 90)), H=10, k=5)))
    print(fmt(study_xs("REV5", lambda P: -past_return(P.close, 5), H=5, k=5)))
    print(fmt(study_xs("LOWVOL", lambda P: -roll_std(P.ret, 30), H=10, k=5)))
    print("\nJSON sample:", as_json(study_xs("MOM30", lambda P: past_return(P.close, 30)))[:240])