mathematical foundations

The mathematics of Drift

Full derivations for every formula. Each section explains the problem, derives the solution from first principles, and notes how it maps to source code.

Layer 1

Data — caching and provider abstraction

Why Parquet?

Parquet is a columnar format. For time-series data where you frequently slice by date range or select specific columns, columnar storage means only the requested columns are read off disk. A row-based format reads every column even when you want one. For a 500-ticker, 5-year OHLCV dataset this difference is 10–50×.

Snappy compression decompresses at ~500 MB/s — 3× faster than gzip at only ~10% size penalty. On repeated pipeline runs this matters more than cold storage size.

Cache key scheme

Each file is addressed by a content hash of (dataset, ticker, start, end, interval). This prevents stale reads across date ranges. The deliberate choice NOT to support partial-range reuse avoids a class of subtle bugs where a cached wide range silently returns stale data for a sub-range request.

Key invariant: the DataLoader is the only file in the codebase that knows OpenBB exists. All layers above get DataFrames. Swapping the data vendor is a one-file change.
Layer 2a

Factor engine — cross-sectional signal construction

The rank transform — derivation from first principles

We want a function mapping ranks r ∈ {1..N} linearly to scores s ∈ {−1,+1}. Write the linear form s = a·r + b. Two boundary conditions pin the two unknowns:

s(1) = −1 → a·1 + b = −1
s(N) = +1 → a·N + b = +1

Subtract the first equation from the second — b cancels:

a(N−1) = 2 → a = 2/(N−1)

Substitute back to find b, then combine over a common denominator:

sᵢ = (2rᵢ − N − 1) / (N − 1)
Verify: r=1 → (2−N−1)/(N−1) = −1 ✓ · r=N → (N−1)/(N−1) = +1 ✓
The rank transform is outlier-proof by construction. An outlier at +∞ gets the same score (+1.00) as a value at rank N. The z-score of that outlier would be +∞.

The seven factors

Fama-French 5 factors: Market beta (systematic risk), Size (−log mktcap), Value (B/M ratio), Profitability (ROE), Investment (−asset_growth). Plus Momentum (12-1 month return, skip most recent month per Jegadeesh & Titman 1993) and Quality (composite z-score of ROE, low leverage, low earnings volatility per Asness et al. 2019).

Layer 2b

HMM regime detection — the forward algorithm

The computational problem

For K=3 states and T=252 trading days, brute-force enumeration of all state sequences requires K^T = 3^252 ≈ 10^120 operations. The forward algorithm reduces this to O(K²T) via dynamic programming.

The forward variable

Define α_t(k) = P(o₁,...,o_t, s_t = k | λ) — the probability of seeing observations up to t AND being in state k at t.

Initialisation: α₁(k) = π_k · B_k(o₁)
Recursion: α_{t+1}(j) = B_j(o_{t+1}) · Σ_k [α_t(k) · A_{kj}]
Termination: P(O|λ) = Σ_k α_T(k)
B_j is outside the sum — emission depends only on the destination state, not the path. This is the Markov property.
Complexity: O(K²T) vs O(K^T). For K=3, T=252: 2,268 ops vs 10^120. Dynamic programming makes the intractable trivial.

Why StandardScaler is mandatory

log_return, realised_vol, vol_of_vol differ in scale by ~10×. Without scaling, Baum-Welch EM collapses all observations into one state. That state's covariance is well-estimated; the minority states have near-zero observations, their covariance matrices become rank-deficient, and the Cholesky decomposition inside the Gaussian emission fails.

Layer 2c

WST features — Cauchy wavelet scattering

Why the modulus is essential

A complex wavelet ψ_j is a bandpass filter. The raw convolution x★ψ_j oscillates at the centre frequency with zero mean — averaging gives nothing. The modulus |x★ψ_j| strips the oscillation and gives the instantaneous amplitude envelope — slowly varying, non-negative, convolvable at the next scale.

Without modulus: S₂ = x★ψ_j₁★ψ_j₂ is just another bandpass filter (composition of linear operators = linear operator). With modulus: S₂ captures how energy at scale j₁ modulates energy at scale j₂ — a genuinely non-linear feature.

Cauchy vs Morlet

Morlet: ψ̂(ω) ∝ exp(−(ω−ω₀)²/2σ²) (Gaussian decay, fast)
Cauchy: ψ̂(ω) ∝ ω^α · exp(−ω/ω_j) (algebraic decay, slower)
α=4 in Drift · higher α → closer to Morlet · lower α → heavier tails, more sensitive to large events

Financial returns are fat-tailed. A Morlet wavelet damps a 10σ return the same as a 3σ return at the same frequency. A Cauchy wavelet retains sensitivity to the 10σ event — exactly the information that matters for systematic risk.

The 37-dimensional coefficient vector

S₀ = ⟨|x|⟩ → 1 coefficient
S₁(j) = ⟨|x★ψ_j|⟩ → J=8 coefficients
S₂(j₁,j₂) = ⟨||x★ψ_j₁|★ψ_j₂|⟩ → C(8,2)=28 coefficients (j₁ < j₂)
Total: 1 + 8 + 28 = 37 · Each coefficient is translation-invariant: same volatility pattern at any point in time → same feature
Layer 3

Alpha engine — IC, ICIR, and signal combination

Why Spearman, not Pearson

A single earnings shock (+30% on one stock, ±2% everywhere else) produces a near-zero Pearson IC even if the factor correctly ranked all other stocks. Spearman rank correlation is Pearson correlation of the rank vectors — a single outlier gets rank N, not a 15σ influence on the correlation.

IC_t = Spearman(s_t, r_{t+h})
s_t are already rank-transformed → rank(s_t) preserves the same ordering

ICIR as signal Sharpe ratio

ICIR = IC̄ / σ_IC
t-stat = ICIR · √T
ICIR 0.3, T=252: t=4.8 (highly significant) · ICIR 0.05: t=0.8 (noise)

IC-weighted combination

Zᵢ = Σ_k w_k · s_{i,k} / Σ_k w_k
w_k = max(ICIR_k, 0) · negative ICIR → zero weight
Regime conditioning: momentum weight → 0 in bear regimes (Barroso & Santa-Clara 2015: momentum crashes occur specifically during market reversals following bear regimes).
Layer 4a

Black-Litterman — Bayesian posterior portfolios

Step 1 — Reverse optimisation

CAPM says market-cap weights w_mkt are mean-variance optimal. Running the formula backwards gives the implied equilibrium expected returns:

μ_eq = λ · Σ · w_mkt
λ = market risk aversion ≈ 2.5

Step 2 — Views from IC-weighted signals

P = I_N (absolute views on each stock)
q = Z · scale(μ_eq) (views in return units)
Ω = τ · diag(PΣPᵀ) (He-Litterman uncertainty)

Step 3 — Posterior (conjugate Gaussian update)

M = [(τΣ)⁻¹ + PᵀΩ⁻¹P]⁻¹ (posterior covariance)
μ_BL = M · [(τΣ)⁻¹μ_eq + PᵀΩ⁻¹q] (posterior mean)
Precision-weighted average: high-confidence views pull μ_BL toward q; low-confidence leaves it near μ_eq
Layer 4b

Hierarchical Risk Parity — clustering-based allocation

Step 1 — Distance matrix

d_{ij} = √((1 − ρ_{ij}) / 2)
ρ=1→d=0 · ρ=0→d=0.707 · ρ=−1→d=1 · satisfies triangle inequality (unlike d=1−ρ)

Step 2 — Ward linkage and quasi-diagonalisation

Ward linkage minimises within-cluster variance at each merge. The dendrogram leaf order (quasi-diagonalisation) reorders the covariance matrix so similar assets are adjacent.

Step 3 — Recursive bisection

Var(cluster) = w_c^T Σ w_c (equal weight within cluster)
α = Var(right) / (Var(left) + Var(right))
w_left ← α · w_parent · w_right ← (1−α) · w_parent
HRP is the default optimiser — it never fails numerically. BL requires a well-conditioned covariance; HRP works even when T < N.
Layer 4c

CVaR — the Rockafellar-Uryasev linear programme

The non-smooth problem

CVaR_α(w) = min_z [ z + 1/((1−α)T) · Σ_t max(−w^T r_t − z, 0) ]

The max() makes this non-smooth. Standard gradient methods fail. Rockafellar-Uryasev (2000) showed it reformulates as a linear programme:

min_{w,z,u} z + 1/((1−α)T) · Σ_t u_t
s.t. u_t ≥ −w^T r_t − z ∀ t
u_t ≥ 0, Σwᵢ = 1, 0 ≤ wᵢ ≤ max_weight
At optimality: z = VaR_α · objective = CVaR_α · u_t captures loss in excess of VaR
Layer 5

Risk model — Ledoit-Wolf shrinkage and BARRA decomposition

Why sample covariance fails

An N×N covariance matrix has N(N+1)/2 unique parameters. When T approaches N, the condition number (max/min eigenvalue) explodes. Inverting such a matrix amplifies noise by the condition number — 10³ or more for a 100-asset portfolio with 2 years of daily data.

Ledoit-Wolf shrinkage

Σ̂ = α · F + (1−α) · S
F = (tr(S)/N)·I (scaled identity target) · α computed analytically from S, no cross-validation

BARRA decomposition

rᵢ = Σ_k Xᵢ_k · f_k + uᵢ
Var(r_p) = w^T (X F X^T + Δ) w
= w^T X F X^T w + w^T Δ w
X = factor exposures (N×K) · F = factor covariance (K×K) · Δ = diagonal specific variance
Layer 6

Backtest — PSR, DSR, and walk-forward

Sharpe ratio variance under non-normality

Var(SR̂) ≈ (1/T) · [1 − γ₃ · SR̂ + (γ₄−1)/4 · SR̂²]
γ₃ = skewness · γ₄ = excess kurtosis · negative skew inflates variance → SR̂ less reliable

Probabilistic Sharpe Ratio

PSR(SR*) = Φ( (SR̂ − SR*) · √(T−1) / √(1 − γ₃·SR̂ + (γ₄−1)/4·SR̂²) )
Φ = standard normal CDF · negative skew → larger denominator → lower PSR → correctly penalises tail-risk harvesting

Deflated Sharpe Ratio

When N strategies are tested and the best is selected, E[max SR] under the null rises with N. DSR adjusts the benchmark SR* to account for selection bias:

SR*_DSR = SR*₀ + σ_{SR̂} · E[max(N)]
E[max(N)] ≈ (1−ρ̄)·Φ⁻¹(1−1/N) + ρ̄·Φ⁻¹(1−1/Ne) · for N=50 independent trials: E[max]≈2.0
PSR(88.6%) means there is an 88.6% probability the true Sharpe exceeds zero — not just that the observed Sharpe is positive. DSR further deducts for how many times you looked before selecting this strategy.