Warrior_EA/research/test_anomalies.py

"""The families this project has NOT tested, chosen by strength of published prior.

Everything measured here so far asked one question - can recent price or flow predict the
NEXT bar's DIRECTION - and answered no four different ways. That is one hypothesis family,
and it happens to be the one most thoroughly arbitraged. These are the others, each with a
real literature behind it rather than a hunch:

  1 TIME-SERIES MOMENTUM (Moskowitz/Ooi/Pedersen 2012). The sign of the past k months
    predicts the next month, across every asset class they tested, for 100+ years. It
    trades RARELY, so the spread that killed the flow signal is amortised over a huge
    holding period - at D1 with a 60-day hold, a 0.3bp spread is ~0.5% of a typical move
    instead of 10% of it. This is the single best untested idea available here.

  2 SEASONALITY / SESSION. Hour-of-day and day-of-week effects in FX are documented and
    persistent (fixing flows, session opens). Cheap to test, and it is a CONDITIONER: even
    a weak directional rule can become tradeable if it only fires when the drift is with it.

  3 OVERNIGHT vs INTRADAY. In equities most of the premium accrues overnight, not during
    the session. SP500 here is a CFD on exactly that.

  4 VOLATILITY-MANAGED DRIFT (Moreira/Muir 2017). The one that needs NO directional edge
    at all: for an asset with positive drift, scaling exposure by inverse recent variance
    raises the Sharpe ratio, because volatility is far more forecastable than returns. If
    nothing else here works, THIS is the honest path to a system - it monetises drift and
    predictable vol rather than a directional forecast that does not exist.

DISCIPLINE, same as test_flow.py
--------------------------------
Costs charged on every entry and exit. Non-overlapping holds so trades are independent.
Buy-and-hold is reported alongside every directional rule, because on a drifting asset a
rule that is merely long most of the time will look skilful and is not. Nulls preserve
what the rule is not claiming and destroy only what it is.
"""
import numpy as np, sys, os, datetime as dt
sys.stdout.reconfigure(encoding='utf-8', errors='replace')

BARS = 'c:/Users/admin/Documents/Workspaces/Market Data/bars/'
SYMS = ['EURUSD', 'USDJPY', 'XAUUSD', 'SP500']


def load(sym, tf):
    z = np.load(f"{BARS}{sym}_{tf}_ticks.npz", allow_pickle=True)
    a = z['bars']
    cols = [str(c) for c in z['columns']]
    I = {c: k for k, c in enumerate(cols)}
    return a, I


def sharpe(r, per_year):
    r = np.asarray(r, float)
    if len(r) < 3 or r.std(ddof=1) == 0:
        return 0.0, 0.0
    s = r.mean() / r.std(ddof=1) * np.sqrt(per_year)
    t = r.mean() / (r.std(ddof=1) / np.sqrt(len(r)))
    return s, t


# ---------------------------------------------------------------- 1 momentum
def momentum(sym, tf='D1', lookbacks=(5, 10, 20, 60, 120, 250), nperm=5000, seed=7):
    a, I = load(sym, tf)
    c = a[:, I['close']]
    spread = a[:, I['spread_mean']]
    n = len(c)
    if n < 400:
        print(f"  {sym} {tf}: only {n} bars, skipping")
        return
    lr = np.diff(np.log(c))
    per_year = {'D1': 252, 'H4': 252 * 6, 'H1': 252 * 24}.get(tf, 252)
    rng = np.random.default_rng(seed)
    #--- round-trip cost as a log return: cross the spread in and out
    cost = float(np.nanmedian(spread / c))
    print(f"\n--- {sym} {tf} time-series momentum  n={n:,}  round-trip cost {2*cost*1e4:.2f} bp ---")
    bh_r = lr
    s_bh, t_bh = sharpe(bh_r, per_year)
    print(f"    buy & hold: Sharpe {s_bh:+.2f}  t {t_bh:+.2f}  "
          f"ann.ret {bh_r.mean()*per_year*100:+.2f}%")
    print(f"    {'look':>5}{'hold':>6}{'trades':>8}{'ann.ret':>9}{'Sharpe':>8}{'t':>7}"
          f"{'p(perm)':>9}{'vs B&H':>8}")
    for L in lookbacks:
        H = L                                    # hold == lookback, the MOP convention
        if n < L + H + 50:
            continue
        entries = np.arange(L, n - H - 1, H)     # NON-OVERLAPPING
        if len(entries) < 30:
            continue
        past = np.log(c[entries]) - np.log(c[entries - L])
        fwd = np.log(c[entries + H]) - np.log(c[entries])
        sig = np.sign(past)
        r = sig * fwd - 2 * cost                 # one round trip per trade
        s, t = sharpe(r, per_year / H)
        # permutation null: shuffle the SIGNAL across trades. Keeps the return series and
        # the long/short mix, destroys only the pairing of signal to period - which is
        # exactly and only what momentum claims.
        obs = r.mean()
        pm = np.empty(nperm)
        for b in range(nperm):
            pm[b] = (rng.permutation(sig) * fwd - 2 * cost).mean()
        p = (1 + np.sum(pm >= obs)) / (nperm + 1)
        bh_match = fwd.mean() - 2 * cost         # always-long over the same periods
        print(f"    {L:>5}{H:>6}{len(entries):>8}{r.mean()*per_year/H*100:>+9.2f}"
              f"{s:>+8.2f}{t:>+7.2f}{p:>9.4f}"
              f"{(r.mean()-bh_match)*per_year/H*100:>+8.2f}")


# ------------------------------------------------------- 2 hour / day-of-week
def seasonality(sym, tf='H1', nperm=5000, seed=11):
    a, I = load(sym, tf)
    c = a[:, I['close']]
    t = a[:, I['time']]
    lr = np.diff(np.log(c))
    hrs = np.array([dt.datetime.fromtimestamp(x / 1000, dt.UTC).hour for x in t[1:]])
    dows = np.array([dt.datetime.fromtimestamp(x / 1000, dt.UTC).weekday() for x in t[1:]])
    rng = np.random.default_rng(seed)
    print(f"\n--- {sym} {tf} seasonality  n={len(lr):,} ---")
    for label, key, K in (('hour', hrs, 24), ('dow', dows, 7)):
        means = np.array([lr[key == k].mean() if (key == k).sum() > 30 else 0.0
                          for k in range(K)])
        cnts = np.array([(key == k).sum() for k in range(K)])
        obs = np.max(np.abs(means))              # max-statistic over the K buckets
        pm = np.empty(nperm)
        for b in range(nperm):
            sh = rng.permutation(lr)
            pm[b] = np.max(np.abs([sh[key == k].mean() if cnts[k] > 30 else 0.0
                                   for k in range(K)]))
        p = (1 + np.sum(pm >= obs)) / (nperm + 1)
        best = int(np.argmax(np.abs(means)))
        print(f"    {label:>5}: strongest bucket {best:>2}  mean {means[best]*1e4:+.2f} bp "
              f"(n={cnts[best]:,})   family-wise p={p:.4f}"
              f"{'  *' if p < 0.05 else ''}")


# ----------------------------------------------------------- 4 vol-managed drift
def vol_managed(sym, tf='D1', win=20):
    """Moreira-Muir. Needs no directional forecast: scale exposure by inverse recent
    variance. Only meaningful where drift is positive, so buy&hold is printed beside it."""
    a, I = load(sym, tf)
    c = a[:, I['close']]
    lr = np.diff(np.log(c))
    n = len(lr)
    if n < win + 50:
        return
    per_year = {'D1': 252, 'H4': 252 * 6}.get(tf, 252)
    #--- variance over the PREVIOUS win bars only, shifted so a bar never scales itself
    var = np.full(n, np.nan)
    csum = np.concatenate([[0.0], np.cumsum(lr ** 2)])
    for i in range(win, n):
        var[i] = (csum[i] - csum[i - win]) / win
    ok = ~np.isnan(var) & (var > 0)
    w = np.zeros(n)
    w[ok] = 1.0 / var[ok]
    w[ok] /= np.nanmean(w[ok])                   # unit average exposure => comparable scale
    w = np.clip(w, 0, 5)                         # no unbounded leverage on a quiet patch
    r_vm = w[:-1] * lr[1:]                       # weight known BEFORE the return it scales
    r_bh = lr[1:]
    s_vm, t_vm = sharpe(r_vm[ok[:-1]], per_year)
    s_bh, t_bh = sharpe(r_bh[ok[:-1]], per_year)
    print(f"    {sym:>7} {tf}: buy&hold Sharpe {s_bh:+.2f}   vol-managed {s_vm:+.2f}   "
          f"delta {s_vm-s_bh:+.2f}   turnover-free")


if __name__ == '__main__':
    which = sys.argv[1] if len(sys.argv) > 1 else 'all'
    if which in ('all', 'mom'):
        print("=" * 78 + "\nTIME-SERIES MOMENTUM (Moskowitz/Ooi/Pedersen)\n" + "=" * 78)
        for s in SYMS:
            for tf in ('D1', 'H4'):
                if os.path.exists(f"{BARS}{s}_{tf}_ticks.npz"):
                    momentum(s, tf)
    if which in ('all', 'seas'):
        print("\n" + "=" * 78 + "\nSEASONALITY\n" + "=" * 78)
        for s in SYMS:
            seasonality(s, 'H1')
    if which in ('all', 'vol'):
        print("\n" + "=" * 78 +
              "\nVOLATILITY-MANAGED DRIFT (Moreira/Muir) - needs no directional edge\n"
              + "=" * 78)
        for s in SYMS:
            for tf in ('D1', 'H4'):
                if os.path.exists(f"{BARS}{s}_{tf}_ticks.npz"):
                    vol_managed(s, tf)