We backtested 7 popular trading rules on US stocks, Europe and forex, after costs

THUZTRA GLOBAL RULES STUDY v1 - FROZEN SPEC (written 2026-11-08, before any result was computed)

Purpose: test widely taught retail trading rules on US stocks, European stocks and three major
FX pairs, after costs, with textbook parameters only (no tuning).

DATA (frozen in raw/, hashes in raw/SHA256SUMS)
- US: raw/F-F_Research_Data_Factors_daily.csv (Kenneth R. French Data Library). Columns are percent.
  Daily total market return R_t = (Mkt-RF + RF)/100. Cash return RF_t = RF/100.
  This is the CRSP value-weighted US market (NYSE, AMEX, NASDAQ) including dividends. Not the S&P 500.
- Europe: raw/Europe_3_Factors_Daily.csv (same library). R_t = (Mkt-RF + RF)/100 in USD, RF_t = RF/100
  (US one-month T-bill). Rows with -99.99 are dropped. Europe includes the UK.
- FX: raw/DEXUSEU.csv (EUR/USD, USD per EUR), raw/DEXUSUK.csv (GBP/USD, USD per GBP),
  raw/DEXJPUS.csv (USD/JPY, JPY per USD). Federal Reserve H.10 noon buying rates via FRED.
  Blank values (holidays) are dropped. Pair return r_t = P_t / P_{t-1} - 1.

PRICE SERIES FOR SIGNALS
- Stocks: P_t = cumulative product of (1 + R_t) over the whole file (total-return index).
- FX: P_t = the quoted rate.
Indicators use the whole file history (so warm-up comes before the evaluation window).

TIMING (no look-ahead)
- A signal s_t uses data up to and including close t. The position held during day t+1 is pos_{t+1} = s_t.
- The day-t return uses pos_t.

RULES (stocks: long = 1, flat = 0. FX: long = +1, short = -1, flat = 0)
1. BH  Buy and hold: s = 1 always (FX: long the pair always).
2. SMA200  Price above 200-day simple moving average: s = 1 if P > SMA200 else 0 (FX: +1 / -1).
3. GC  50/200 moving-average crossover: s = 1 if SMA50 > SMA200 else 0 (FX: +1 / -1).
4. MACD  12/26/9: MACD = EMA12(P) - EMA26(P); signal = EMA9(MACD); s = 1 if MACD > signal else 0 (FX: +1/-1).
   EMA_n: alpha = 2/(n+1), recursive, seeded with the first value (pandas ewm(adjust=False)).
5. RSI  RSI(14) 30/70 (Wilder): d = P_t - P_{t-1}; gain = max(d,0); loss = max(-d,0);
   avg = ewm(alpha=1/14, adjust=False) of gain and loss starting at the first diff; RSI = 100 - 100/(1+avgGain/avgLoss)
   (avgLoss = 0 -> RSI = 100). State machine, starts at 0:
   stocks: go to 1 when RSI < 30; go to 0 when RSI > 70; else keep state.
   FX: go to +1 when RSI < 30; go to -1 when RSI > 70; else keep state.
6. BOLL  Bollinger(20, 2) reversion: mid = SMA20, sd = population std (ddof=0) of last 20 P.
   stocks: enter 1 when P < mid - 2 sd; exit to 0 when P >= mid; else keep.
   FX: as stocks for longs; also enter -1 when P > mid + 2 sd; exit -1 to 0 when P <= mid.
   (From flat, a long entry is tested before a short entry.)
7. DON  Donchian 20/10 breakout (close-based Turtle rule): hi20 = max of P over the previous 20 days (excluding today),
   lo10 = min over previous 10 days, lo20 / hi10 likewise.
   stocks: enter 1 when P > hi20; exit to 0 when P < lo10; else keep.
   FX: also enter -1 when P < lo20; exit -1 to 0 when P > hi10. From flat, long tested first.
8. MOM12  12-month time-series momentum: on the last trading day of each calendar month,
   s = 1 if P_t / P_{t-252} - 1 > 0 else 0 (FX: +1 / -1). Held until the next month-end decision.
   Before the first month-end decision with 252 days of history: 0.
Undefined indicator values (warm-up) -> s = 0.

COSTS (per side = per unit of |change in position|)
- US: 5 bps per side (index-ETF spread + slippage). ETF fee 0.09% a year, charged as 0.0009/252 per day times |pos|.
- Europe: 10 bps per side. ETF fee 0.20% a year likewise.
- FX: 2 bps per side (retail spread + slippage, major pairs). Financing mark-up 1.00% a year, charged 0.01/252 per day times |pos|.
  The interest-rate differential (carry) is NOT modelled. Taxes are NOT modelled.
- A move from -1 to +1 is 2 units of change.
Daily net return:
- stocks: ret_t = pos_t * R_t + (1 - pos_t) * RF_t - cost_side * |pos_t - pos_{t-1}| - fee_daily * pos_t
- FX:     ret_t = pos_t * r_t - cost_side * |pos_t - pos_{t-1}| - fin_daily * |pos_t|
  (FX results are unlevered excess returns: cash collateral interest is not added.)
At the start of the evaluation window the previous position is taken as 0 (so the first entry pays its cost).
Gross = the same with all costs and fees set to 0.

EVALUATION WINDOWS (fixed in advance)
- US headline: 1993-02-01 to 2026-09-30 (first US-listed index ETF launched 1993-01-29).
- US long history (robustness, costs are optimistic before 1993): 1927-01-03 to 2026-09-30.
- Europe: 1991-07-01 to 2026-09-30.
- FX: 2000-01-03 to 2026-10-30 (all three pairs).
- Split check: the evaluation window is cut at its calendar midpoint; Sharpe is reported for each half.

METRICS (on daily net returns inside the window; equity starts at 1 on the last date before the window)
- years = (last date - date before first window date in the file).days / 365.25
- CAGR = equity_end ** (1/years) - 1
- k = N_days / years (bars per year)
- Volatility = std(ret, ddof=1) * sqrt(k)
- Sharpe: stocks use x = ret - RF; FX uses x = ret. Sharpe = mean(x) / std(x, ddof=1) * sqrt(k)
- Max drawdown = min over t of (equity_t / running max equity - 1), running max includes the starting 1.
- Stocks: % of full calendar years in the window where the rule's compounded return > buy and hold's (strictly).
- FX: % of full calendar years with a positive compounded return.
- Switches per year = count of days with pos_t != pos_{t-1} (including the first entry) / years.
- Exposure = mean of |pos_t|.

"""Thuztra global rules study v1. Implements SPEC.txt exactly.

The page /reference/trading-rules-backtested (backend/src/trading_rules_study.ts)
renders out/page_data.json, copied to backend/src/trading_rules_study.json.

The raw files are not in this repo: the French library holds their copyright, so
only derived statistics are published. To reproduce, download into raw/ and check
the hashes in SHA256SUMS (taken 2026-11-08; the French library revises its files):
  https://mba.tuck.dartmouth.edu/pages/faculty/ken.french/ftp/F-F_Research_Data_Factors_daily_CSV.zip
  https://mba.tuck.dartmouth.edu/pages/faculty/ken.french/ftp/Europe_3_Factors_Daily_CSV.zip
  https://fred.stlouisfed.org/graph/fredgraph.csv?id=DEXUSEU   (saved as DEXUSEU.csv)
  https://fred.stlouisfed.org/graph/fredgraph.csv?id=DEXUSUK   (saved as DEXUSUK.csv)
  https://fred.stlouisfed.org/graph/fredgraph.csv?id=DEXJPUS   (saved as DEXJPUS.csv)
Unzip both zips in raw/, then run: python3 study.py  (needs pandas and numpy).
check.py is the independent second calculation of two rules (GC, RSI) on US and EUR/USD,
written from SPEC.txt alone; run.txt is study.py's printed output.
"""
import csv
import datetime as dt
import json
import math
import os
import sys

import numpy as np
import pandas as pd

HERE = os.path.dirname(os.path.abspath(__file__))
RAW = os.path.join(HERE, "raw")
OUT = os.path.join(HERE, "out")

RULES = ["BH", "SMA200", "GC", "MACD", "RSI", "BOLL", "DON", "MOM12"]


def load_french(name):
    rows = []
    with open(os.path.join(RAW, name)) as f:
        for line in f:
            parts = [p.strip() for p in line.split(",")]
            if len(parts) == 5 and parts[0].isdigit() and len(parts[0]) == 8:
                mkt_rf, rf = float(parts[1]), float(parts[4])
                if mkt_rf <= -99.98 or rf <= -99.98:
                    continue
                rows.append((pd.Timestamp(parts[0]), (mkt_rf + rf) / 100.0, rf / 100.0))
    df = pd.DataFrame(rows, columns=["date", "R", "RF"]).set_index("date").sort_index()
    df["P"] = (1.0 + df["R"]).cumprod()
    return df


def load_fred(name):
    df = pd.read_csv(os.path.join(RAW, name))
    col = df.columns[1]
    df = df[pd.to_numeric(df[col], errors="coerce").notna()].copy()
    df["date"] = pd.to_datetime(df["observation_date"])
    df["P"] = df[col].astype(float)
    df = df.set_index("date")[["P"]].sort_index()
    df["R"] = df["P"].pct_change()
    return df


def ema(s, n):
    return s.ewm(alpha=2.0 / (n + 1), adjust=False).mean()


def rsi(p, n=14):
    d = p.diff()
    gain = d.clip(lower=0)
    loss = (-d).clip(lower=0)
    ag = gain.iloc[1:].ewm(alpha=1.0 / n, adjust=False).mean()
    al = loss.iloc[1:].ewm(alpha=1.0 / n, adjust=False).mean()
    out = pd.Series(np.nan, index=p.index)
    with np.errstate(divide="ignore", invalid="ignore"):
        r = np.where(al.values == 0, 100.0, 100.0 - 100.0 / (1.0 + ag.values / al.values))
    out.iloc[1:] = r
    return out


def signals(p, fx):
    """Return dict rule -> signal series s_t (decided at close t)."""
    neg = -1.0 if fx else 0.0
    idx = p.index
    sig = {}
    sig["BH"] = pd.Series(1.0, index=idx)
    sma200 = p.rolling(200).mean()
    sma50 = p.rolling(50).mean()
    sig["SMA200"] = pd.Series(np.where(sma200.isna(), 0.0, np.where(p > sma200, 1.0, neg)), index=idx)
    sig["GC"] = pd.Series(np.where(sma200.isna(), 0.0, np.where(sma50 > sma200, 1.0, neg)), index=idx)
    macd = ema(p, 12) - ema(p, 26)
    sigl = ema(macd, 9)
    sig["MACD"] = pd.Series(np.where(macd > sigl, 1.0, neg), index=idx)

    rs = rsi(p).values
    st, cur = np.zeros(len(p)), 0.0
    for i, v in enumerate(rs):
        if not np.isnan(v):
            if v < 30:
                cur = 1.0
            elif v > 70:
                cur = -1.0 if fx else 0.0
        st[i] = cur
    sig["RSI"] = pd.Series(st, index=idx)

    mid = p.rolling(20).mean().values
    sd = p.rolling(20).std(ddof=0).values
    pv = p.values
    st, cur = np.zeros(len(p)), 0.0
    for i in range(len(p)):
        if not np.isnan(mid[i]):
            lo, hi = mid[i] - 2 * sd[i], mid[i] + 2 * sd[i]
            if cur == 0.0:
                if pv[i] < lo:
                    cur = 1.0
                elif fx and pv[i] > hi:
                    cur = -1.0
            elif cur == 1.0:
                if pv[i] >= mid[i]:
                    cur = 0.0
            elif cur == -1.0:
                if pv[i] <= mid[i]:
                    cur = 0.0
        st[i] = cur
    sig["BOLL"] = pd.Series(st, index=idx)

    hi20 = p.shift(1).rolling(20).max().values
    lo20 = p.shift(1).rolling(20).min().values
    hi10 = p.shift(1).rolling(10).max().values
    lo10 = p.shift(1).rolling(10).min().values
    st, cur = np.zeros(len(p)), 0.0
    for i in range(len(p)):
        if not np.isnan(hi20[i]):
            if cur == 0.0:
                if pv[i] > hi20[i]:
                    cur = 1.0
                elif fx and pv[i] < lo20[i]:
                    cur = -1.0
            elif cur == 1.0:
                if pv[i] < lo10[i]:
                    cur = 0.0
            elif cur == -1.0:
                if pv[i] > hi10[i]:
                    cur = 0.0
        st[i] = cur
    sig["DON"] = pd.Series(st, index=idx)

    ret12 = p / p.shift(252) - 1.0
    month = idx.to_period("M")
    is_me = np.r_[month[1:] != month[:-1], True]
    st, cur = np.zeros(len(p)), 0.0
    r12 = ret12.values
    for i in range(len(p)):
        if is_me[i] and not np.isnan(r12[i]):
            cur = 1.0 if r12[i] > 0 else neg
        st[i] = cur
    sig["MOM12"] = pd.Series(st, index=idx)
    return sig


def run_rule(df, s, start, end, fx, side, fee_yr, gross=False):
    pos_all = s.shift(1).fillna(0.0)
    m = (df.index >= start) & (df.index <= end)
    sub = df.loc[m]
    pos = pos_all.loc[m].values
    prev = np.r_[0.0, pos[:-1]]
    c = 0.0 if gross else side
    f = 0.0 if gross else fee_yr / 252.0
    if fx:
        ret = pos * sub["R"].values - c * np.abs(pos - prev) - f * np.abs(pos)
        x = ret
    else:
        ret = pos * sub["R"].values + (1 - pos) * sub["RF"].values - c * np.abs(pos - prev) - f * pos
        x = ret - sub["RF"].values
    d0 = df.index[df.index < sub.index[0]][-1]
    return pd.Series(ret, index=sub.index), x, pos, prev, d0


def metrics(ret, x, pos, prev, d0):
    eq = (1 + ret).cumprod()
    years = (ret.index[-1] - d0).days / 365.25
    k = len(ret) / years
    peak = np.maximum.accumulate(np.r_[1.0, eq.values])[1:]
    return {
        "cagr": float(eq.iloc[-1] ** (1 / years) - 1),
        "vol": float(np.std(ret.values, ddof=1) * math.sqrt(k)),
        "sharpe": float(np.mean(x) / np.std(x, ddof=1) * math.sqrt(k)),
        "maxdd": float(np.min(eq.values / peak - 1)),
        "switches_per_year": float(np.sum(pos != prev) / years),
        "exposure": float(np.mean(np.abs(pos))),
        "final_equity": float(eq.iloc[-1]),
        "years": years,
        "days": int(len(ret)),
    }


def yearly(ret):
    y = (1 + ret).groupby(ret.index.year).prod() - 1
    first, last = ret.index[0], ret.index[-1]
    full = [yr for yr in y.index if not (yr == first.year and first > pd.Timestamp(yr, 1, 1) + pd.Timedelta(days=7))
            and not (yr == last.year and last < pd.Timestamp(yr, 12, 24))]
    return y.loc[full]


def study(name, df, start, end, fx, side, fee_yr):
    start, end = pd.Timestamp(start), pd.Timestamp(end)
    end = min(end, df.index[-1])
    sig = signals(df["P"], fx)
    mid = start + (end - start) / 2
    out = {}
    bh_y = None
    for rule in RULES:
        ret, x, pos, prev, d0 = run_rule(df, sig[rule], start, end, fx, side, fee_yr)
        g = run_rule(df, sig[rule], start, end, fx, side, fee_yr, gross=True)
        mt = metrics(ret, x, pos, prev, d0)
        mt["gross_cagr"] = metrics(*g)["cagr"]
        halves = []
        for a, b in ((start, mid), (mid + pd.Timedelta(days=1), end)):
            r2, x2, p2, pv2, d2 = run_rule(df, sig[rule], a, b, fx, side, fee_yr)
            halves.append(metrics(r2, x2, p2, pv2, d2)["sharpe"])
        mt["sharpe_first_half"], mt["sharpe_second_half"] = halves
        y = yearly(ret)
        if rule == "BH":
            bh_y = y
        if fx:
            mt["pct_years_positive"] = float((y > 0).mean())
        else:
            mt["pct_years_beat_bh"] = float((y > bh_y).mean()) if rule != "BH" else None
        mt["full_years"] = int(len(y))
        mt["yearly"] = {str(k): float(v) for k, v in y.items()}
        out[rule] = mt
    return {
        "market": name,
        "start": str(ret.index[0].date()),
        "end": str(ret.index[-1].date()),
        "split_date": str(mid.date()),
        "fx": fx,
        "cost_side_bps": side * 1e4,
        "annual_fee_pct": fee_yr * 100,
        "rules": out,
    }


def main():
    us = load_french("F-F_Research_Data_Factors_daily.csv")
    eu = load_french("Europe_3_Factors_Daily.csv")
    fx = {k: load_fred(f) for k, f in (("EURUSD", "DEXUSEU.csv"), ("GBPUSD", "DEXUSUK.csv"), ("USDJPY", "DEXJPUS.csv"))}
    res = {
        "generated_by": "study.py",
        "markets": {
            "US": study("US stock market", us, "1993-02-01", "2026-09-30", False, 0.0005, 0.0009),
            "US_LONG": study("US stock market, long history", us, "1927-01-03", "2026-09-30", False, 0.0005, 0.0009),
            "EUROPE": study("European stock market", eu, "1991-07-01", "2026-09-30", False, 0.0010, 0.0020),
        },
    }
    for k, d in fx.items():
        res["markets"][k] = study(k, d, "2000-01-03", "2026-10-30", True, 0.0002, 0.01)
    os.makedirs(OUT, exist_ok=True)
    with open(os.path.join(OUT, "results.json"), "w") as f:
        json.dump(res, f, indent=1, sort_keys=True)
    page = {"markets": {mk: {**{k: v for k, v in m.items() if k != "rules"},
                             "rules": {r: {k: v for k, v in d.items() if k != "yearly"} for r, d in m["rules"].items()}}
                        for mk, m in res["markets"].items()}}
    with open(os.path.join(OUT, "page_data.json"), "w") as f:
        json.dump(page, f, indent=1, sort_keys=True)
        f.write("\n")
    with open(os.path.join(OUT, "results.csv"), "w", newline="") as f:
        w = csv.writer(f)
        w.writerow(["market", "rule", "start", "end", "cagr_net", "cagr_gross", "sharpe", "volatility", "max_drawdown",
                    "pct_years_beat_buy_and_hold", "pct_years_positive", "switches_per_year", "exposure",
                    "sharpe_first_half", "sharpe_second_half", "cost_per_side_bps", "annual_fee_or_financing_pct"])
        for mk, m in res["markets"].items():
            for rule, r in m["rules"].items():
                w.writerow([mk, rule, m["start"], m["end"], f"{r['cagr']:.4f}", f"{r['gross_cagr']:.4f}",
                            f"{r['sharpe']:.3f}", f"{r['vol']:.4f}", f"{r['maxdd']:.4f}",
                            "" if r.get("pct_years_beat_bh") is None else f"{r['pct_years_beat_bh']:.3f}",
                            "" if r.get("pct_years_positive") is None else f"{r['pct_years_positive']:.3f}",
                            f"{r['switches_per_year']:.2f}", f"{r['exposure']:.3f}",
                            f"{r['sharpe_first_half']:.3f}", f"{r['sharpe_second_half']:.3f}",
                            f"{m['cost_side_bps']:.0f}", f"{m['annual_fee_pct']:.2f}"])
    for mk, m in res["markets"].items():
        print(f"== {mk} {m['start']}..{m['end']} split {m['split_date']}")
        for rule, r in m["rules"].items():
            extra = r.get("pct_years_beat_bh") if not m["fx"] else r.get("pct_years_positive")
            print(f"  {rule:7s} cagr {r['cagr']*100:6.2f}% gross {r['gross_cagr']*100:6.2f}% sharpe {r['sharpe']:5.2f} "
                  f"mdd {r['maxdd']*100:6.1f}% yrs {('' if extra is None else f'{extra*100:4.0f}%'):>5s} "
                  f"sw/yr {r['switches_per_year']:5.1f} exp {r['exposure']:.2f} halves {r['sharpe_first_half']:.2f}/{r['sharpe_second_half']:.2f}")


if __name__ == "__main__":
    sys.exit(main())

53598ec5af305d321ee93db492a32076fc882f976e67d27eba724d0266d20ad2  Europe_3_Factors_Daily_CSV.zip
96c23b88e651cbfff3740a53f8225e59ee39893371e96385a1d8125ec2cd9c89  F-F_Research_Data_Factors_daily_CSV.zip
60ba955456192c4a1cbf6a1f4dee25fb040b2120bdcea0050a4e55b04cf8e997  DEXJPUS.csv
3d95f828a707b568522fc90b3efb08a82169777536b3f9a131ad2df065f3a791  DEXUSEU.csv
75c6228bdb54f33f37dec3059b688f4e43dd6c5e88ce5905e69ed2a3b02567c5  DEXUSUK.csv
9d76d27f486206667a6bcdf18dde487b11505ba637e8159f23810ecdf3e265b7  Europe_3_Factors_Daily.csv
3f7910502dd6f188748210e5145d95642de1d8cff56c2cf0fa85985299610ff0  F-F_Research_Data_Factors_daily.csv