""" ml/robustness.py — Robustness testing: Monte Carlo, stress tests, multi-scenario validation. Capabilities: 1. Monte Carlo equity curve simulation (shuffled returns, bootstrapped trades) 2. Stress tests: higher costs, slippage spikes, execution delays, gap risk 3. Parameter sensitivity analysis (grid search over key params) 4. Walk-forward stability scoring """ from __future__ import annotations import logging from dataclasses import dataclass, field import numpy as np logger = logging.getLogger("trading_system.ml.robustness") @dataclass class MonteCarloResult: """Monte Carlo simulation results.""" n_simulations: int median_return: float p5_return: float # 5th percentile (worst case) p25_return: float p75_return: float p95_return: float # 95th percentile (best case) prob_positive: float # probability of positive total return median_sharpe: float median_max_dd: float p95_max_dd: float # worst-case max drawdown (95th pct) simulated_returns: np.ndarray # (n_sims,) final total returns @dataclass class StressTestResult: """Results from a single stress scenario.""" scenario_name: str original_return: float stressed_return: float original_sharpe: float stressed_sharpe: float original_max_dd: float stressed_max_dd: float degradation_pct: float # % change in return @dataclass class RobustnessReport: """Full robustness test results.""" monte_carlo: MonteCarloResult | None = None stress_tests: list[StressTestResult] = field(default_factory=list) parameter_sensitivity: dict = field(default_factory=dict) overall_score: float = 0.0 # 0-100 robustness score class MonteCarloSimulator: """Monte Carlo equity curve simulation via return bootstrapping.""" def __init__( self, n_simulations: int = 1000, block_size: int = 5, # block bootstrap to preserve autocorrelation seed: int = 42, ): self.n_simulations = n_simulations self.block_size = block_size self.rng = np.random.default_rng(seed) def simulate(self, trade_returns: np.ndarray) -> MonteCarloResult: """Run Monte Carlo simulation on trade returns. Uses block bootstrapping to preserve short-term return autocorrelation. Args: trade_returns: Array of per-trade returns (e.g., +0.02, -0.01). Returns: MonteCarloResult with distribution statistics. """ n_trades = len(trade_returns) if n_trades < 10: logger.warning("Too few trades for Monte Carlo (%d)", n_trades) return MonteCarloResult( n_simulations=0, median_return=0, p5_return=0, p25_return=0, p75_return=0, p95_return=0, prob_positive=0, median_sharpe=0, median_max_dd=0, p95_max_dd=0, simulated_returns=np.array([]), ) n_blocks = max(1, n_trades // self.block_size) sim_total_returns = np.zeros(self.n_simulations) sim_sharpes = np.zeros(self.n_simulations) sim_max_dds = np.zeros(self.n_simulations) for i in range(self.n_simulations): # Block bootstrap block_starts = self.rng.integers( 0, n_trades - self.block_size + 1, size=n_blocks, ) sampled = np.concatenate([ trade_returns[s:s + self.block_size] for s in block_starts ])[:n_trades] # match original length # Compute equity curve equity = np.cumprod(1 + sampled) sim_total_returns[i] = equity[-1] - 1 # Sharpe if np.std(sampled) > 0: sim_sharpes[i] = (np.mean(sampled) / np.std(sampled)) * np.sqrt(252) else: sim_sharpes[i] = 0.0 # Max drawdown peak = np.maximum.accumulate(equity) dd = (peak - equity) / (peak + 1e-9) sim_max_dds[i] = dd.max() return MonteCarloResult( n_simulations=self.n_simulations, median_return=float(np.median(sim_total_returns)), p5_return=float(np.percentile(sim_total_returns, 5)), p25_return=float(np.percentile(sim_total_returns, 25)), p75_return=float(np.percentile(sim_total_returns, 75)), p95_return=float(np.percentile(sim_total_returns, 95)), prob_positive=float(np.mean(sim_total_returns > 0)), median_sharpe=float(np.median(sim_sharpes)), median_max_dd=float(np.median(sim_max_dds)), p95_max_dd=float(np.percentile(sim_max_dds, 95)), simulated_returns=sim_total_returns, ) class StressTester: """Apply stress scenarios to trade results.""" def __init__(self, base_return: float = 0.0, base_sharpe: float = 0.0, base_max_dd: float = 0.0): self.base_return = base_return self.base_sharpe = base_sharpe self.base_max_dd = base_max_dd def stress_higher_costs( self, trade_returns: np.ndarray, cost_multiplier: float = 2.0, base_cost_pct: float = 0.03, ) -> StressTestResult: """Stress: double transaction costs.""" extra_cost = base_cost_pct * (cost_multiplier - 1) / 100 stressed = trade_returns - extra_cost return self._build_result("2x Transaction Costs", trade_returns, stressed) def stress_slippage_spike( self, trade_returns: np.ndarray, slippage_pct: float = 0.05, ) -> StressTestResult: """Stress: add extra slippage to all trades.""" stressed = trade_returns - slippage_pct / 100 return self._build_result(f"Slippage +{slippage_pct}bps", trade_returns, stressed) def stress_execution_delay( self, trade_returns: np.ndarray, delay_penalty_pct: float = 0.02, ) -> StressTestResult: """Stress: execution delay causes adverse price movement.""" penalties = np.abs(trade_returns) * delay_penalty_pct stressed = trade_returns - penalties return self._build_result("Execution Delay", trade_returns, stressed) def stress_tail_events( self, trade_returns: np.ndarray, n_worst: int = 5, amplify_factor: float = 2.0, ) -> StressTestResult: """Stress: amplify the N worst trades (tail risk).""" stressed = trade_returns.copy() sorted_idx = np.argsort(trade_returns) for i in sorted_idx[:n_worst]: if stressed[i] < 0: stressed[i] *= amplify_factor return self._build_result("Tail Risk (2x worst 5)", trade_returns, stressed) def stress_gap_risk( self, trade_returns: np.ndarray, gap_probability: float = 0.02, gap_size_pct: float = 3.0, seed: int = 42, ) -> StressTestResult: """Stress: random overnight gap events.""" rng = np.random.default_rng(seed) stressed = trade_returns.copy() n = len(stressed) n_gaps = max(1, int(n * gap_probability)) gap_idx = rng.choice(n, size=n_gaps, replace=False) direction = rng.choice([-1, 1], size=n_gaps) stressed[gap_idx] += direction * gap_size_pct / 100 return self._build_result("Random Gap Risk", trade_returns, stressed) def _build_result( self, name: str, original: np.ndarray, stressed: np.ndarray, ) -> StressTestResult: o_ret = float(np.sum(original)) s_ret = float(np.sum(stressed)) o_sharpe = self._sharpe(original) s_sharpe = self._sharpe(stressed) o_dd = self._max_dd(original) s_dd = self._max_dd(stressed) deg = ((s_ret - o_ret) / abs(o_ret) * 100) if abs(o_ret) > 1e-9 else 0 return StressTestResult( scenario_name=name, original_return=o_ret, stressed_return=s_ret, original_sharpe=o_sharpe, stressed_sharpe=s_sharpe, original_max_dd=o_dd, stressed_max_dd=s_dd, degradation_pct=deg, ) @staticmethod def _sharpe(returns: np.ndarray) -> float: if len(returns) == 0 or np.std(returns) == 0: return 0.0 return float((np.mean(returns) / np.std(returns)) * np.sqrt(252)) @staticmethod def _max_dd(returns: np.ndarray) -> float: equity = np.cumprod(1 + returns) peak = np.maximum.accumulate(equity) dd = (peak - equity) / (peak + 1e-9) return float(dd.max()) if len(dd) > 0 else 0.0 def run_all(self, trade_returns: np.ndarray) -> list[StressTestResult]: """Run all stress scenarios.""" return [ self.stress_higher_costs(trade_returns), self.stress_slippage_spike(trade_returns), self.stress_execution_delay(trade_returns), self.stress_tail_events(trade_returns), self.stress_gap_risk(trade_returns), ] def compute_robustness_score( mc: MonteCarloResult, stress_tests: list[StressTestResult], ) -> float: """Compute overall robustness score (0-100). Scoring: 40% Monte Carlo: probability of positive returns + consistency 30% Stress resilience: avg degradation across scenarios 30% Drawdown control: worst-case drawdown """ score = 0.0 # MC component (40 pts) if mc.n_simulations > 0: # Prob positive: 0-20 pts score += mc.prob_positive * 20 # Sharpe consistency: 0-10 pts score += min(max(mc.median_sharpe, 0), 3.0) / 3.0 * 10 # P5 not too negative: 0-10 pts p5_score = max(0, min(1.0, (mc.p5_return + 0.20) / 0.20)) score += p5_score * 10 # Stress component (30 pts) if stress_tests: degradations = [abs(st.degradation_pct) for st in stress_tests] avg_deg = np.mean(degradations) # Lower degradation = higher score stress_score = max(0, 30 - avg_deg * 0.3) score += stress_score # Drawdown component (30 pts) if mc.n_simulations > 0: # P95 max DD: scale 0-30 (0% DD = 30, 30% DD = 0) dd_score = max(0, 30 * (1 - mc.p95_max_dd / 0.30)) score += dd_score return float(np.clip(score, 0, 100)) def print_robustness_report(report: RobustnessReport): """Pretty-print the robustness report.""" print("\n" + "=" * 70) print(" ROBUSTNESS TESTING REPORT") print("=" * 70) if report.monte_carlo: mc = report.monte_carlo print(f"\n MONTE CARLO ({mc.n_simulations} simulations)") print(f" {'Metric':<25s} {'Value':>10s}") print(f" {'-'*35}") print(f" {'Median Return':<25s} {mc.median_return:>9.2%}") print(f" {'5th Percentile':<25s} {mc.p5_return:>9.2%}") print(f" {'25th Percentile':<25s} {mc.p25_return:>9.2%}") print(f" {'75th Percentile':<25s} {mc.p75_return:>9.2%}") print(f" {'95th Percentile':<25s} {mc.p95_return:>9.2%}") print(f" {'P(Return > 0)':<25s} {mc.prob_positive:>9.1%}") print(f" {'Median Sharpe':<25s} {mc.median_sharpe:>9.2f}") print(f" {'Median Max DD':<25s} {mc.median_max_dd:>9.2%}") print(f" {'95th pct Max DD':<25s} {mc.p95_max_dd:>9.2%}") if report.stress_tests: print(f"\n STRESS TESTS") print(f" {'Scenario':<25s} {'Orig':>8s} {'Stressed':>8s} {'Degrad':>8s}") print(f" {'-'*49}") for st in report.stress_tests: print( f" {st.scenario_name:<25s}" f" {st.original_return:>7.2%}" f" {st.stressed_return:>7.2%}" f" {st.degradation_pct:>7.1f}%" ) print(f"\n ROBUSTNESS SCORE: {report.overall_score:.0f}/100") if report.overall_score >= 70: print(" Status: ROBUST ✓") elif report.overall_score >= 50: print(" Status: MODERATE — review stress scenarios") else: print(" Status: WEAK — significant improvements needed") print("=" * 70)