🤖 Phase 5 — Quant Strategies, ML & Backtesting

import pandas as pd
import statsmodels.api as sm
 
# Fama-French factor regression: decompose a stock's returns into factor exposures
# (Factor data available from Kenneth French's data library)
stock_returns = pd.Series(...)  # your stock's excess returns
factors = pd.DataFrame({
    'MKT': ..., 'SMB': ..., 'HML': ..., 'UMD': ...
})
 
X = sm.add_constant(factors)
model = sm.OLS(stock_returns, X).fit()
print(model.params)
# alpha (const): return NOT explained by factor exposures -- this is what active managers claim to add
# betas: how much exposure the stock has to each factor
print(f"Alpha (annualized): {model.params['const'] * 252:.2%}")
print(f"R-squared: {model.rsquared:.2%}")  # how much of the variance factors explain

from statsmodels.tsa.stattools import coint
import numpy as np
 
def find_cointegrated_pairs(price_df, significance=0.05):
    """Test all pairs for cointegration."""
    n = price_df.shape[1]
    pairs = []
    for i in range(n):
        for j in range(i+1, n):
            asset1, asset2 = price_df.columns[i], price_df.columns[j]
            score, pvalue, _ = coint(price_df[asset1], price_df[asset2])
            if pvalue < significance:
                pairs.append((asset1, asset2, pvalue))
    return sorted(pairs, key=lambda x: x[2])
 
# Pairs trading signal: z-score of the spread
def pairs_trading_signal(price1, price2, lookback=60, entry_z=2.0, exit_z=0.5):
    """
    Generate long/short signals based on the spread's z-score.
    Spread = price1 - hedge_ratio * price2 (hedge ratio from rolling OLS)
    """
    spread = price1 - price2  # simplified; in practice use rolling OLS hedge ratio
    rolling_mean = spread.rolling(lookback).mean()
    rolling_std = spread.rolling(lookback).std()
    z_score = (spread - rolling_mean) / rolling_std
 
    signal = pd.Series(0, index=spread.index)
    signal[z_score > entry_z] = -1   # spread too wide: short asset1, long asset2
    signal[z_score < -entry_z] = 1   # spread too narrow: long asset1, short asset2
    signal[z_score.abs() < exit_z] = 0  # close position
    return signal.replace(to_replace=0, method='ffill').fillna(0)  # hold position until exit

from sklearn.model_selection import TimeSeriesSplit
import numpy as np
 
# WRONG: random train/test split shuffles time -- catastrophic look-ahead bias
# from sklearn.model_selection import train_test_split
# X_train, X_test, y_train, y_test = train_test_split(X, y, shuffle=True)  # NEVER DO THIS
 
# CORRECT: time series split -- train always precedes test
tscv = TimeSeriesSplit(n_splits=5)
for train_idx, test_idx in tscv.split(X):
    X_train, X_test = X.iloc[train_idx], X.iloc[test_idx]
    y_train, y_test = y.iloc[train_idx], y.iloc[test_idx]
    # train and evaluate model here
 
# BETTER for overlapping labels: purged k-fold cross-validation (Lopez de Prado)
# Removes training samples whose label periods OVERLAP with the test period
def purged_train_test_split(X, y, label_horizon, test_start, test_end):
    """Remove training samples whose label window overlaps with test period."""
    test_mask = (X.index >= test_start) & (X.index <= test_end)
    purge_start = test_start - pd.Timedelta(days=label_horizon)
    train_mask = (X.index < purge_start) | (X.index > test_end)
    return X[train_mask], X[test_mask], y[train_mask], y[test_mask]

# Problem: raw prices are non-stationary (bad for ML), but full differencing (returns)
# destroys ALL memory/predictive structure.
# Fractional differentiation: achieve stationarity while preserving maximum memory
 
def frac_diff(series, d, threshold=1e-4):
    """Fractionally differentiate a series (Lopez de Prado method)."""
    weights = [1.0]
    for k in range(1, len(series)):
        w = -weights[-1] / k * (d - k + 1)
        if abs(w) < threshold:
            break
        weights.append(w)
    weights = np.array(weights[::-1])
 
    result = pd.Series(index=series.index, dtype=float)
    width = len(weights)
    for i in range(width, len(series)):
        result.iloc[i] = np.dot(weights, series.iloc[i-width:i])
    return result.dropna()

class EventDrivenBacktester:
    """
    Simplified event-driven backtester.
    Processes one bar at a time: signal -> order -> fill -> portfolio update
    """
    def __init__(self, initial_capital=100_000, commission=0.001, slippage=0.0005):
        self.cash = initial_capital
        self.positions = {}  # symbol -> shares
        self.commission = commission
        self.slippage = slippage
        self.equity_curve = []
        self.trades = []
 
    def execute_order(self, symbol, target_shares, price, timestamp):
        """Execute a trade with realistic costs."""
        current_shares = self.positions.get(symbol, 0)
        shares_to_trade = target_shares - current_shares
        if shares_to_trade == 0:
            return
 
        # Apply slippage: buys execute at higher price, sells at lower
        execution_price = price * (1 + self.slippage * np.sign(shares_to_trade))
        trade_value = abs(shares_to_trade) * execution_price
        commission_cost = trade_value * self.commission
 
        self.cash -= shares_to_trade * execution_price + commission_cost
        self.positions[symbol] = target_shares
 
        self.trades.append({
            'timestamp': timestamp, 'symbol': symbol,
            'shares': shares_to_trade, 'price': execution_price,
            'commission': commission_cost
        })
 
    def mark_to_market(self, prices: dict, timestamp):
        """Compute total equity given current prices."""
        position_value = sum(
            shares * prices.get(sym, 0)
            for sym, shares in self.positions.items()
        )
        equity = self.cash + position_value
        self.equity_curve.append({'timestamp': timestamp, 'equity': equity})
        return equity
 
    def run(self, data: pd.DataFrame, strategy_fn):
        """
        data: DataFrame with MultiIndex (timestamp, symbol) and OHLCV columns
        strategy_fn(history, current_positions) -> dict of {symbol: target_shares}
        """
        for timestamp in data.index.get_level_values(0).unique():
            bar = data.loc[timestamp]
            history = data.loc[:timestamp]  # only data UP TO now -- prevents look-ahead
 
            target_positions = strategy_fn(history, self.positions)
 
            for symbol, target in target_positions.items():
                if symbol in bar.index:
                    self.execute_order(symbol, target, bar.loc[symbol, 'close'], timestamp)
 
            prices = {sym: bar.loc[sym, 'close'] for sym in bar.index}
            self.mark_to_market(prices, timestamp)
 
        return pd.DataFrame(self.equity_curve).set_index('timestamp')

def strategy_performance_report(equity_curve: pd.Series, risk_free_rate=0.02):
    """Compute standard performance metrics from an equity curve."""
    returns = equity_curve.pct_change().dropna()
    n_years = len(returns) / 252
 
    total_return = equity_curve.iloc[-1] / equity_curve.iloc[0] - 1
    cagr = (1 + total_return) ** (1/n_years) - 1
    annual_vol = returns.std() * np.sqrt(252)
    sharpe = (returns.mean() * 252 - risk_free_rate) / annual_vol
 
    # Sortino: like Sharpe but only penalizes downside volatility
    downside_returns = returns[returns < 0]
    downside_vol = downside_returns.std() * np.sqrt(252)
    sortino = (returns.mean() * 252 - risk_free_rate) / downside_vol
 
    # Maximum drawdown
    cumulative = (1 + returns).cumprod()
    running_max = cumulative.cummax()
    drawdown = (cumulative - running_max) / running_max
    max_drawdown = drawdown.min()
 
    # Calmar ratio: CAGR / |max drawdown|
    calmar = cagr / abs(max_drawdown)
 
    return {
        'CAGR': f"{cagr:.2%}",
        'Annual Volatility': f"{annual_vol:.2%}",
        'Sharpe Ratio': f"{sharpe:.2f}",
        'Sortino Ratio': f"{sortino:.2f}",
        'Max Drawdown': f"{max_drawdown:.2%}",
        'Calmar Ratio': f"{calmar:.2f}",
        'Win Rate': f"{(returns > 0).mean():.2%}"
    }

def walk_forward_backtest(data, strategy_class, param_grid, train_window=252*2, test_window=63):
    """
    Walk-forward optimization:
    1. Optimize parameters on a training window
    2. Test (out-of-sample) on the following window
    3. Roll forward, repeat
    This simulates how a strategy would ACTUALLY be deployed and re-tuned over time.
    """
    results = []
    start = 0
    while start + train_window + test_window <= len(data):
        train_data = data.iloc[start : start + train_window]
        test_data = data.iloc[start + train_window : start + train_window + test_window]
 
        # Find best params on training data ONLY
        best_params, best_sharpe = None, -np.inf
        for params in param_grid:
            strategy = strategy_class(**params)
            train_equity = strategy.backtest(train_data)
            sharpe = strategy_performance_report(train_equity)['Sharpe Ratio']
            if float(sharpe) > best_sharpe:
                best_sharpe, best_params = float(sharpe), params
 
        # Apply best params to UNSEEN test data
        strategy = strategy_class(**best_params)
        test_equity = strategy.backtest(test_data)
        results.append({
            'period_start': test_data.index[0],
            'params': best_params,
            'oos_performance': strategy_performance_report(test_equity)
        })
 
        start += test_window  # roll forward
 
    return results