Save the following complete program as residual_momentum.py and run it in a Python environment with NumPy installed (python -m pip install numpy). Then run python residual_momentum.py. The program includes every raw return and factor row, converts percentages once, uses SVD least squares with an intercept, rejects invalid inputs, and prints coefficients, residual vectors, scores, ranks and timestamps. It needs no network access or proprietary factor library.
Synthetic data only. Three assets, twelve monthly rows, eight-month regressions and a three-month signal with December skipped. The function is intentionally limited to this dated fixture.python
"""Synthetic teaching example. Requires Python 3 and NumPy."""
import calendar
import csv
import io
import json
from datetime import datetime, timezone
import numpy as np
RAW = """month,Mkt_RF,SMB,HML,RF,A,B,C
2025-01,1.00,0.30,-0.20,0.20,2.58,0.24,1.49
2025-02,-2.00,0.50,0.40,0.20,-4.02,-0.10,-1.77
2025-03,3.00,-0.60,0.20,0.20,6.00,1.12,3.38
2025-04,0.50,0.20,-0.70,0.20,1.89,0.26,0.38
2025-05,-1.00,-0.40,0.60,0.20,-2.34,-0.12,-0.62
2025-06,2.00,0.70,0.10,0.20,4.75,0.81,1.99
2025-07,-0.50,-0.20,-0.40,0.20,-1.16,0.09,0.00
2025-08,1.50,0.10,0.80,0.20,2.90,0.93,1.47
2025-09,4.00,0.40,-0.20,0.20,7.42,3.22,4.60
2025-10,3.00,-0.30,0.50,0.20,5.43,2.41,3.47
2025-11,5.00,0.20,-0.60,0.20,8.96,4.36,5.70
2025-12,2.00,0.60,0.30,0.20,4.55,0.63,2.42
"""
UTC = timezone.utc
IDS = ('A', 'B', 'C')
def stamp(year, month, day=2, hour=0):
return datetime(year, month, day, hour, tzinfo=UTC)
def data():
rows = list(csv.DictReader(io.StringIO(RAW)))
months = [row['month'] for row in rows]
# All numeric source fields are percent per month; convert exactly once.
values = np.array([[float(row[k]) / 100 for k in
('Mkt_RF', 'SMB', 'HML', 'RF', *IDS)] for row in rows])
ends, released = [], []
for month in months:
y, m = map(int, month.split('-'))
ends.append(stamp(y, m, calendar.monthrange(y, m)[1]))
ny, nm = (y + 1, 1) if m == 12 else (y, m + 1)
released.append(stamp(ny, nm, 1, 1))
return months, values, ends, released
def ols(factors, excess):
x = np.column_stack((np.ones(len(factors)), factors))
if len(x) <= x.shape[1] or not np.isfinite(x).all() or not np.isfinite(excess).all():
raise ValueError('missing, non-finite or insufficient observations')
coefficients, _, rank, singular = np.linalg.lstsq(x, excess, rcond=None)
if rank != 4 or singular[0] / singular[-1] > 1e8:
raise ValueError('singular or ill-conditioned factor matrix')
return coefficients
def score(residuals):
sd = np.std(residuals, axis=0, ddof=1)
if not np.isfinite(sd).all() or np.any(sd <= 1e-12):
raise ValueError('zero or invalid residual volatility')
return residuals.sum(axis=0) / sd
def rank(values):
return sorted(IDS, key=lambda name: (-values[IDS.index(name)], name))
def calculate(values, factor_ends, asset_ends, factor_release, asset_release,
decision=stamp(2026, 1), fit_function=ols):
if values.shape != (12, 7) or len(factor_ends) != 12 or len(asset_ends) != 12:
raise ValueError('expected twelve complete monthly rows')
expected = data()[2]
if factor_ends != expected or asset_ends != expected:
raise ValueError('different or missing calendar endpoints')
if len(factor_release) != 12 or len(asset_release) != 12:
raise ValueError('missing availability timestamps')
if not np.isfinite(values).all() or np.any(values[:, 4:] < -1):
raise ValueError('invalid returns or missing values')
f, rf, returns = values[:, :3], values[:, 3], values[:, 4:]
y = returns - rf[:, None]
def known(indices, cutoff):
if any(max(factor_release[j], asset_release[j]) > cutoff for j in indices):
raise ValueError('data unavailable at declared cutoff')
known(range(4, 12), decision)
# At January decision: one May-December fit, then September-November scores.
full = fit_function(f[4:12], y[4:12])
published_style = y[8:11] - f[8:11] @ full[1:]
refitted_errors = published_style - full[0]
frozen, audit = [], []
for j in range(8, 11):
fit_at = stamp(2025, j + 1)
residual_at = stamp(2025, j + 2)
known(range(j - 8, j), fit_at)
known([j], residual_at)
if residual_at > decision:
raise ValueError('residual not yet available')
b = fit_function(f[j - 8:j], y[j - 8:j])
frozen.append(y[j] - b[0] - f[j] @ b[1:])
audit.append({'month': j + 1, 'fit_at': fit_at.isoformat(),
'residual_at': residual_at.isoformat(),
'coefficients': b.tolist()})
frozen = np.array(frozen)
raw = np.prod(1 + returns[8:11], axis=0) - 1
signals = {'raw': raw, 'published_style': score(published_style),
'refitted_ols_error': score(refitted_errors),
'frozen_ols_error': score(frozen)}
return {'scores': {k: v.tolist() for k, v in signals.items()},
'ranks': {k: rank(v) for k, v in signals.items()},
'published_style_returns': published_style.tolist(),
'frozen_errors': frozen.tolist(), 'final_coefficients': full.tolist(),
'audit': audit}
if __name__ == '__main__':
_, values, ends, released = data()
print(json.dumps(calculate(values, ends, ends, released, released), indent=2))