""" Factor census: every stock factor the lab has, measured one way on the current US panel (PIT S&P 500 after the Tiingo merge and the pre-2011 membership fix). Earlier per-factor numbers (factor_scan_us.csv, qqq_scan_factors.csv, ranks_us.pkl) predate that panel. PYTHONPATH=. uv run --with numba python -m research.lab.factor_census Factors: every research/lab/factors.py factor for "us", plus rec_63 (the recovery_mom input); the EDGAR fundamentals (research/lab/fundamentals.py; size = -log market cap); the alt-data panels (insider Form 4, Lazy Prices, earnings-announcement return, FINRA short-volume ratio lagged 2 sessions, SEC fails-to-deliver lagged 21 calendar days); and three house composites (mean member percentile ranks): house_overnight_recovery overnight_mom + rec_126 (overnight_recovery_top10) house_recovery_mom rec_63 + mom_12_1 (recovery_mom_top10) house_low_risk ivol + low_beta (lowrisk_top20) Every score is "higher = better" as the library defines it. A negative spread means the factor works the other way round. Measures (members with a score and a forward return, >= 50 names a day): ic daily Spearman rank IC against the open(t+1) -> close(t+21) return ls top-quintile minus bottom-quintile 21-day return, annualised x 12 top_x top quintile minus the member mean, annualised t Newey-West, 21 lags, on the daily overlapping series turn share of the top quintile replaced after 21 days top20 long-only top 20, every 21 days, buffer 2, next-open fills, 10 bp a side, phase-averaged (us_sweep.run): excess Sharpe, CAGR, MaxDD, and alpha against the equal-weight member basket (Newey-West t), from the factor's coverage start corr correlation of the ls series with house_overnight_recovery's and mom_12_1's; span_t is the t of the intercept when the factor's ls is regressed on house_overnight_recovery's ls Periods: 2006-10..2016 and 2017..2026-09 (halves), 2014..2026-09 (after the survivor-biased years). A factor whose data start later uses what it has; a half with less than 2 years is left out. Classes, fixed before the run: A works alone: ls t >= 2 over the span, ls > 0 in every usable half, and top-20 alpha vs EW > 0 in every usable half B long-short only: ls t >= 2 and > 0 in every half, but the long-only top 20 has no alpha in some half C unstable or weak: ls t >= 2 with a half <= 0, or 1 <= t < 2 with every half > 0 R works reversed: ls t <= -2 N nothing: the rest About 70 factors are tested at once, so roughly 3 will pass t >= 2 by chance. The t >= 3 bar (Harvey, Liu & Zhu 2016) is reported too. Descriptive census, no strategy gate. The top-20 runs go to trials_us.jsonl. The script refuses to overwrite census.csv. --part construction (reported): the house inputs and the census's strongest long-only rows at the live overnight_recovery_top10 construction (top 10, every 10 days, buffer 2), vs the census construction. Writes construction.csv. Result (2026-10-12), 86 factors: - Classes: 1 A, 2 B, 11 C, 3 R, 69 N. Six spreads reach |t| >= 2, about what chance gives at ~80 tests. No spread reaches |t| >= 3. - A low_svr21 (2018-26 only): IC t 3.0, spread +4.2%/yr (t 2.5), top-20 alpha vs EW +5.5% (t 2.8), corr 0.27 with the house composite, span t 2.7. - B accruals: spread +2.8% (t 2.35), positive in both halves. Top-20 alpha +4.2% (t 1.75), but 2010-16 is -0.8%. - B rev_5: spread +4.3% (t 2.8), both halves. 77% of the top quintile turns over each month, and the net top-20 alpha is -0.7%. - R max_ret -6.9% (t -2.0): lottery names won. overnight_mom_1 and overnight_20 -5.2% / -4.5% (t -2.3 / -2.1): last month's overnight gains reverse. - Momentum and the house composites have no quintile spread (mom_12_1 -4.5%, t -1.05; 2006-16 -10%, 2017-26 +1.5%; overnight_mom -1.5%, t -0.56). Their edge is only in the extreme top names. Top 10 / every 10 / buffer 2 vs EW: house_overnight_recovery +6.4% (t 1.40), halves -0.2 / +13.6 overnight_mom alone +4.4% (t 0.97), halves +4.4 / +4.7 rec_126 alone +2.0% mom_12_1 -0.4% At top 20 / every 21: +2.9 / +1.3 / -0.3 / -2.2. - Low-risk factors have negative raw spreads (low_vol -6.4%, low_beta -6.1%) because high beta won the bull market. On a beta-adjusted basis house_low_risk is +0.9% (t 0.57). - Long-only alpha with both halves positive (fundamentals start in 2010): roe +4.0% (t 2.75), net_issuance +3.1% (2.0), sp +3.0% (1.4; top 10 +7.1%, t 2.5, found after looking). Their quintile spreads are flat. - CN volume factors invert in US large caps: pv_corr top-20 alpha -3.8% (t -3.1), turn_vol -3.5% (-2.9), vol_cv -3.4% (-2.8). - Boundary: S&P 500 only, with pre-2014 survivor bias. The fundamental and short-side data start later, so their 2006-16 half is shorter or missing. """ from __future__ import annotations import argparse import os import pickle import time import warnings import numpy as np import pandas as pd from research.lab import fundamentals as fu from research.lab.alpha import _nw_se, alpha_stats from research.lab.factor_scan import RESULTS, load from research.lab.factors import FACTORS from research.lab.levered_book import stats, tbill from research.lab.us_sweep import run warnings.filterwarnings("ignore") OUT = os.path.join(RESULTS, "factor_census") END = "2026-09-24" SPAN = ("2006-10-02", END) HALVES = {"h1": ("2006-10-02", "2016-12-30"), "h2": ("2017-01-03", END)} POST = ("2014-01-02", END) H = 21 FAMILY = { **dict.fromkeys(["mom_12_1", "mom_6_1", "mom_12_7", "mom_3", "resid_mom", "resid_mom_6", "sharpe_mom", "high_52w", "fip", "mom_fip", "rec_126", "rec_63", "ma200_dist", "ts_trend", "sector_mom", "trend_hzz", "amp_cut_mom", "stmom_hiturn", "seasonality"], "momentum/trend"), **dict.fromkeys(["overnight_mom", "overnight_mom_6", "overnight_mom_1", "overnight_20", "intraday_rev"], "overnight/intraday"), **dict.fromkeys(["rev_5", "rev_21", "rev_resid_21", "rev_sector_21"], "reversal"), **dict.fromkeys(["ivol", "ivol_20", "low_vol", "low_beta", "max_ret", "skew", "ideal_amp"], "low risk/lottery"), **dict.fromkeys(["amihud", "abn_volume", "low_turnover", "abn_turnover", "turn_vol", "pv_corr", "upvol_ratio", "abn_vol_ratio", "vol_cv"], "volume/liquidity"), **dict.fromkeys(["ep", "cfp", "sp", "bm", "fcfp", "shareholder_yield", "size"], "value/size"), **dict.fromkeys(["gp_a", "op_a", "cfo_a", "roe", "roa", "gross_margin", "fscore", "leverage", "accruals"], "quality/profitability"), **dict.fromkeys(["asset_growth", "net_issuance"], "investment/issuance"), **dict.fromkeys(["rev_growth", "gm_change", "roa_change", "sue", "sue_rev", "ear", "ear_z", "ear_count"], "growth/earnings"), } def nw_mean_t(x: pd.Series, lags: int = H) -> tuple[float, float, int]: v = x.dropna().to_numpy() n = len(v) if n < 250: return float("nan"), float("nan"), n e = v - v.mean() s = e @ e / n for k in range(1, lags + 1): s += 2 * (1 - k / (lags + 1)) * (e[k:] @ e[:-k]) / n return float(v.mean()), float(v.mean() / np.sqrt(s / n)), n def cross_section(score: pd.DataFrame, fwd: pd.DataFrame, member: pd.DataFrame): valid = member & fwd.notna() & score.notna() n = valid.sum(axis=1) S, F = score.where(valid), fwd.where(valid) q = S.rank(axis=1, pct=True) fr = F.rank(axis=1, pct=True) qa, fa = q.sub(q.mean(axis=1), axis=0), fr.sub(fr.mean(axis=1), axis=0) ic = (qa * fa).sum(axis=1) / np.sqrt((qa ** 2).sum(axis=1) * (fa ** 2).sum(axis=1)) top, bot = q > 0.8, q <= 0.2 ok = (n >= 50) & (top.sum(axis=1) >= 10) & (bot.sum(axis=1) >= 10) ls = (F.where(top).mean(axis=1) - F.where(bot).mean(axis=1)).where(ok) tx = (F.where(top).mean(axis=1) - F.mean(axis=1)).where(ok) stay = (top & top.shift(-H, fill_value=False)).sum(axis=1) / top.sum(axis=1).replace(0, np.nan) return ic.where(ok), ls, tx, (1 - stay).where(ok), n def span_t(y: pd.Series, x: pd.Series) -> float: d = pd.concat([y, x], axis=1).dropna() if len(d) < 250: return float("nan") X = np.c_[np.ones(len(d)), d.iloc[:, 1].to_numpy()] b, *_ = np.linalg.lstsq(X, d.iloc[:, 0].to_numpy(), rcond=None) e = d.iloc[:, 0].to_numpy() - X @ b return float(b[0] / _nw_se(X, e, H)[0]) def pct_rank(df: pd.DataFrame, member: pd.DataFrame) -> pd.DataFrame: return df.where(member).rank(axis=1, pct=True) def scores(p) -> dict[str, pd.DataFrame]: cols, dates = p.close.columns, p.dates S = {} for n, f in FACTORS.items(): if "us" in f.markets: S[n] = f(p).reindex(index=dates, columns=cols) c = p.close.ffill() S["rec_63"] = c / c.rolling(63, min_periods=63).min() - 1 S["house_overnight_recovery"] = (pct_rank(S["overnight_mom"], p.member) + pct_rank(S["rec_126"], p.member)) / 2 S["house_recovery_mom"] = (pct_rank(S["rec_63"], p.member) + pct_rank(S["mom_12_1"], p.member)) / 2 S["house_low_risk"] = (pct_rank(S["ivol"], p.member) + pct_rank(S["low_beta"], p.member)) / 2 F = fu.load(p) for k, v in F.items(): v = v.reindex(index=dates, columns=cols).astype(float) if k == "mcap": S["size"] = -np.log(v) else: S[k] = v for path, keep in (("insider_factors.pkl", None), ("lazy_prices_factors.pkl", None), ("earnings_factors.pkl", ["ear", "ear_z", "ear_count"])): with open(os.path.join("data", "research", path), "rb") as fh: d = pickle.load(fh) for k, v in d.items(): if keep is None or k in keep: S[k] = v.reindex(index=dates, columns=cols).astype(float) from research.lab.short_volume import factors as svr_factors, load_svr for k, v in svr_factors(load_svr(p.close)).items(): if k in ("svr21", "svr_z"): S[f"low_{k}"] = -v.shift(2).reindex(index=dates, columns=cols) from research.lab.fails import PUBLISH_LAG_DAYS, load_ftd ftd = load_ftd(p.close, p.volume)["ftd_int"] ftd.index = ftd.index + pd.Timedelta(days=PUBLISH_LAG_DAYS) S["low_ftd_int"] = -(ftd.reindex(index=dates.union(ftd.index), columns=cols).ffill(limit=45).reindex(dates)) for k in list(S): if k.startswith(("insider_", "lazy_")): FAMILY.setdefault(k, "alt: insider" if k.startswith("insider_") else "alt: text (Lazy Prices)") elif k.startswith("low_svr") or k == "low_ftd_int": FAMILY.setdefault(k, "alt: short side") elif k.startswith("house_"): FAMILY.setdefault(k, "house composite") return S def classify(r: dict) -> str: halves = [r[f"ls_{h}"] for h in HALVES if np.isfinite(r.get(f"ls_{h}", np.nan))] alphas = [r[f"alpha_{h}"] for h in HALVES if np.isfinite(r.get(f"alpha_{h}", np.nan))] t = r["ls_t"] if not np.isfinite(t): return "N" if t <= -2: return "R" if t >= 2 and halves and all(h > 0 for h in halves): return "A" if alphas and all(a > 0 for a in alphas) else "B" if t >= 2 or (t >= 1 and halves and all(h > 0 for h in halves)): return "C" return "N" CONSTRUCTION = ["house_overnight_recovery", "overnight_mom", "rec_126", "house_recovery_mom", "rec_63", "mom_12_1", "house_low_risk", "low_svr21", "roe", "net_issuance", "accruals", "sp"] def construction(p, S) -> None: out_csv = os.path.join(OUT, "construction.csv") if os.path.exists(out_csv): raise SystemExit(f"{out_csv} exists; this census does not overwrite results") ret = p.close.pct_change(fill_method=None) ew = ret.where(p.member.shift(1, fill_value=False)).mean(axis=1) rf = tbill(p.dates).fillna(0.0) sp = p.slice("2006-06-01", END) rows = [] for name in CONSTRUCTION: sc = S[name].where(p.member) n = sc.notna().sum(axis=1) a = max(pd.Timestamp(SPAN[0]), n[n >= 100].index[0] + pd.Timedelta(days=31)) for label, kw in (("top10_e10_b2", dict(top_n=10, every=10, buffer=2.0)), ("top20_e21_b2", dict(top_n=20, every=H, buffer=2.0))): rr = run(sp, sc, [(a.date().isoformat(), END)], **kw)["_rets"].loc[a:END] s, al = stats(rr, rf), alpha_stats(rr, ew.reindex(rr.index), rf.reindex(rr.index)) r = {"factor": name, "construction": label, "start": a.date().isoformat(), "cagr": float(s["cagr"]), "sharpe": float(s["sharpe"]), "mdd": float(s["mdd"]), "alpha": float(al["alpha"]), "alpha_t": float(al["t"]), "beta": float(al["beta"])} for h, (ha, hb) in HALVES.items(): seg = rr.loc[ha:hb] r[f"alpha_{h}"] = (float(alpha_stats(seg, ew.reindex(seg.index), rf.reindex(seg.index))["alpha"]) if len(seg) >= 2 * 252 else float("nan")) rows.append(r) print(f"{name:26s} {label:13s} cagr {r['cagr']*100:5.1f}% sh {r['sharpe']:.2f} mdd {r['mdd']*100:6.1f}% " f"alpha {r['alpha']*100:+5.1f}% (t {r['alpha_t']:+.2f}) beta {r['beta']:.2f} halves " f"{r['alpha_h1']*100:+5.1f} / {r['alpha_h2']*100:+5.1f}", flush=True) pd.DataFrame(rows).to_csv(out_csv, index=False) from research.lab.trials import append_trials append_trials(os.path.join(RESULTS, "trials_us.jsonl"), [{"name": f"factor_census:{r['factor']}:{r['construction']}", **{k: r[k] for k in ("cagr", "sharpe", "mdd", "alpha", "alpha_t")}} for r in rows], {"study": "factor_census construction", "role": "reported; previously inspected history"}) def main() -> None: ap = argparse.ArgumentParser() ap.add_argument("--part", choices=["census", "construction"], default="census") part = ap.parse_args().part os.makedirs(OUT, exist_ok=True) if part == "construction": p, _ = load("us") construction(p, scores(p)) return out_csv = os.path.join(OUT, "census.csv") if os.path.exists(out_csv): raise SystemExit(f"{out_csv} exists; this census does not overwrite results") p, _ = load("us") t0 = time.time() S = scores(p) print(f"{len(S)} factors built [{time.time() - t0:.0f}s]", flush=True) o = p.open.where(p.open > 0) fwd = p.close.ffill().shift(-H) / o.shift(-1) - 1 # a delisted name keeps its last print member = p.member ret = p.close.pct_change(fill_method=None) ew = ret.where(member.shift(1, fill_value=False)).mean(axis=1) rf = tbill(p.dates).fillna(0.0) sp = p.slice("2006-06-01", END) LS, rows = {}, [] for name, sc in S.items(): t = time.time() ic, ls, tx, turn, n = cross_section(sc, fwd, member) ic, ls, tx, turn, n = (x.loc[SPAN[0]:SPAN[1]] for x in (ic, ls, tx, turn, n)) LS[name] = ls cov = n[n >= 100] start = cov.index[0] if len(cov) else None r = {"factor": name, "family": FAMILY.get(name, "other"), "start": start.date().isoformat() if start is not None else None, "names_2025": float(n.loc["2025"].mean()) if len(n.loc["2025"]) else float("nan")} r["ic"], r["ic_t"], _ = nw_mean_t(ic) m, r["ls_t"], _ = nw_mean_t(ls) r["ls"] = m * 12 m, r["top_x_t"], _ = nw_mean_t(tx) r["top_x"] = m * 12 r["turn"] = float(turn.mean()) for h, (a, b) in {**HALVES, "post2014": POST}.items(): m, tt, nn = nw_mean_t(ls.loc[a:b]) ok = nn >= 2 * 252 r[f"ls_{h}"], r[f"ls_t_{h}"] = (m * 12, tt) if ok else (float("nan"), float("nan")) if start is not None: a = max(pd.Timestamp(SPAN[0]), start + pd.Timedelta(days=31)) st = run(sp, sc.where(member), [(a.date().isoformat(), END)], top_n=20, buffer=2.0, every=H) rr = st["_rets"].loc[a:END] s = stats(rr, rf) al = alpha_stats(rr, ew.reindex(rr.index), rf.reindex(rr.index)) r.update({"top20_sharpe": float(s["sharpe"]), "top20_cagr": float(s["cagr"]), "top20_mdd": float(s["mdd"]), "alpha": float(al["alpha"]), "alpha_t": float(al["t"]), "top20_turn": float(st["turn"])}) for h, (ha, hb) in HALVES.items(): seg = rr.loc[ha:hb] if len(seg) >= 2 * 252: r[f"alpha_{h}"] = float(alpha_stats(seg, ew.reindex(seg.index), rf.reindex(seg.index))["alpha"]) rows.append(r) print(f"{name:26s} ic {r['ic']:+.3f} (t {r['ic_t']:+.1f}) ls {r['ls']*100:+6.1f}% t {r['ls_t']:+5.2f} " f"[{r.get('ls_h1', np.nan)*100:+5.1f} / {r.get('ls_h2', np.nan)*100:+5.1f}] " f"top20 sh {r.get('top20_sharpe', np.nan):.2f} alpha {r.get('alpha', np.nan)*100:+5.1f}% " f"(t {r.get('alpha_t', np.nan):+.1f}) [{time.time() - t:.0f}s]", flush=True) L = pd.DataFrame(LS) for r in rows: y = L[r["factor"]] r["corr_house"] = float(y.corr(L["house_overnight_recovery"])) r["corr_mom"] = float(y.corr(L["mom_12_1"])) r["span_t_house"] = span_t(y, L["house_overnight_recovery"]) r["class"] = classify(r) r["t3"] = bool(abs(r["ls_t"]) >= 3) if np.isfinite(r["ls_t"]) else False df = pd.DataFrame(rows) df.to_csv(out_csv, index=False) L.to_pickle(os.path.join(OUT, "ls_returns.pkl")) with pd.option_context("display.width", 250, "display.max_rows", 200): show = df.sort_values(["family", "ls_t"], ascending=[True, False]) cols = ["family", "factor", "class", "start", "ic", "ic_t", "ls", "ls_t", "ls_h1", "ls_h2", "ls_t_post2014", "top20_sharpe", "alpha", "alpha_t", "alpha_h1", "alpha_h2", "corr_house", "span_t_house", "turn"] print(show[cols].round(3).to_string(index=False)) from research.lab.trials import append_trials append_trials(os.path.join(RESULTS, "trials_us.jsonl"), [{"name": f"factor_census:{r['factor']}:top20", **{k: r[k] for k in ("top20_sharpe", "top20_cagr", "top20_mdd", "alpha", "alpha_t") if k in r}} for r in rows], {"study": "factor_census", "role": "descriptive census; previously inspected history"}) if __name__ == "__main__": main()