"""Service-level checks: the scenario optimiser reproduces the exact Gaussian (mean-variance-equivalent) optimum when the scenarios are moment-matched Sobol' points, the drift convention is arithmetic, and importance sampling / row generation leave the optimum unchanged.""" import numpy as np import pandas as pd import pytest from quantcore.optimization.problem import MeanCVaRProblem from quantcore.optimization.service import OptimizationRequest, OptimizationService from tests.conftest import SYMBOLS BASE = dict(symbols=SYMBOLS, cash=500_000.0, cvar_limit=12_000.0, solver="highs", horizon="week", costs={"proportional_buy": 0.0, "proportional_sell": 0.0}, importance_sampling="off", constraints={"max_weight": 0.15, "max_sector_weight": 0.4}) @pytest.fixture(scope="module") def svc(settings, md): return OptimizationService(settings, md) def gaussian_cvar_of(svc, req, result): """Closed-form CVaR of a proposal under the Gaussian model of the request.""" symbols = svc.candidate_symbols(req) mu_ann, _, _ = svc.expected_returns(req, symbols) spec, _ = svc.build_gaussian_spec(req, symbols, mu_ann) shares = {p["symbol"]: p["shares"] for p in result.positions} value = np.array([shares.get(s, 0) for s in symbols]) * spec.prices return MeanCVaRProblem(spec).gaussian_risk(value, result.total_cost) def test_sobol_moment_matched_scenarios_reproduce_the_exact_optimum(svc): exact = svc.optimize(OptimizationRequest(scenario_model="gaussian_exact", **BASE)) assert exact.status == "Optimal" and exact.reproducible and exact.solver_info["outer_approximation"]["converged"] scen = svc.optimize(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **BASE)) assert scen.status == "Optimal" and scen.model["scenarios"]["moment_matched"] # The GBM scenarios are lognormal in simple returns while the exact mode is normal: at weekly # volatilities of 4-5 % the lognormal left tail is ~2-3 % thinner, so the scenario optimum takes # ~3 % more exposure. Within that model difference the two agree (the sharp test of the sampling # itself, with Gaussian scenarios used linearly, is in test_rowgen_gaussian.py). assert abs(scen.expected_pnl - exact.expected_pnl) <= 0.05 * exact.expected_pnl _, cvar_true, _ = gaussian_cvar_of(svc, OptimizationRequest(scenario_model="gaussian_exact", **BASE), scen) assert cvar_true <= BASE["cvar_limit"] * 1.05 # the expected horizon return reported per position is the model's (arithmetic) expectation mu_ann, _, _ = svc.expected_returns(OptimizationRequest(**BASE), SYMBOLS) for p in scen.positions: assert abs(p["expected_return_horizon"] - mu_ann[p["symbol"]] / 252 * 5) < 1e-12 # plain Monte Carlo without moment matching: the sample mean disagrees with the model mean mc = svc.optimize(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, sampling="mc", moment_match=False, **BASE)) assert mc.status == "Optimal" dev_mc = abs(mc.solver_info["objective"] - mc.expected_pnl) # objective uses the model mean either way assert dev_mc < 1e-6 # (the objective is always the model expectation) ... scen_mean_mc = svc.build_spec(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, sampling="mc", moment_match=False, **BASE), SYMBOLS, svc.scenarios(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, sampling="mc", moment_match=False, **BASE), SYMBOLS, mu_ann)) # ... but the sample mean of a plain MC set is visibly off the model mean, a moment-matched one is not err_mc = np.abs(scen_mean_mc.scenarios.mean(axis=0) - scen_mean_mc.expected_returns).max() sobol_spec = svc.build_spec(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **BASE), SYMBOLS, svc.scenarios(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **BASE), SYMBOLS, mu_ann)) err_sobol = np.abs(sobol_spec.scenarios.mean(axis=0) - sobol_spec.expected_returns).max() assert err_sobol < 0.1 * err_mc def test_importance_sampling_and_row_generation_keep_the_optimum(svc): plain = svc.optimize(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **BASE)) tilted = svc.optimize(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **{**BASE, "importance_sampling": "market"})) assert tilted.status == "Optimal" and tilted.model["scenarios"]["weighted"] assert tilted.model["scenarios"]["effective_sample_size"] > 2000 assert abs(tilted.expected_pnl - plain.expected_pnl) <= 0.02 * plain.expected_pnl exact_req = OptimizationRequest(scenario_model="gaussian_exact", **BASE) assert gaussian_cvar_of(svc, exact_req, tilted)[1] <= BASE["cvar_limit"] * 1.05 rg = svc.optimize(OptimizationRequest(scenario_model="gbm", n_scenarios=4096, **{**BASE, "row_generation": "on"})) assert rg.status == "Optimal" and rg.solver_info["row_generation"]["converged"] assert rg.scenario_hash == plain.scenario_hash # same scenario set ... assert abs(rg.expected_pnl - plain.expected_pnl) <= 1e-4 * plain.expected_pnl # ... same optimum assert max(rg.solver_info["row_generation"]["working_set_sizes"]) < 4096 def test_empirical_model_options(svc): req = OptimizationRequest(scenario_model="empirical", n_scenarios=2048, **{**BASE, "importance_sampling": "market"}) res = svc.optimize(req) assert res.status == "Optimal" and res.model["scenarios"]["model"] == "empirical" assert all(t > 0 for t in res.model["scenarios"]["importance_sampling"]["theta"]) assert res.model["scenarios"]["importance_sampling"]["direction"] == "market" assert res.cvar <= BASE["cvar_limit"] * (1 + 1e-6) bal = svc.optimize(OptimizationRequest(scenario_model="empirical", n_scenarios=1000, sampling="balanced", **BASE)) assert bal.status == "Optimal" and bal.model["scenarios"]["n_scenarios"] == 1000 with pytest.raises(ValueError): OptimizationRequest(scenario_model="bogus", **BASE) with pytest.raises(ValueError): OptimizationRequest(**{**BASE, "importance_sampling": "maybe"}) def test_one_scenario_set_for_a_family_of_runs(svc): """Risk level, constraints, costs, objective, cash and holdings must not change the scenarios; a different expected-return model must change only their drift, not their shocks.""" svc.clear_scenario_cache() base = dict(BASE, importance_sampling="market", scenario_model="gbm", n_scenarios=2048) ref = svc.optimize(OptimizationRequest(**base)) variants = { "other CVaR limit": dict(base, cvar_limit=6_000.0), "limit as a fraction": dict(base, cvar_limit=None, cvar_limit_fraction=0.02), "confidence 99 %": dict(base, confidence=0.99, cvar_limit=20_000.0), "tighter caps": dict(base, constraints={"max_weight": 0.10, "max_sector_weight": 0.3}), "with costs": dict(base, costs={"proportional_buy": 0.002, "proportional_sell": 0.002}), "min-CVaR objective": dict(base, objective="min_cvar", cvar_limit=None, target_return=0.001), "utility objective": dict(base, objective="utility", cvar_limit=None, risk_aversion=0.2), "more cash": dict(base, cash=2_000_000.0), "rebalancing": dict(base, mode="rebalance", cash=50_000.0, holdings={p["symbol"]: p["shares"] for p in ref.positions[:3]}), } for label, kwargs in variants.items(): res = svc.optimize(OptimizationRequest(**kwargs)) assert res.status == "Optimal", label assert res.scenario_key["id"] == ref.scenario_key["id"], label assert res.scenario_hash == ref.scenario_hash, label assert res.solver_info["scenario_cache"]["reused"], label # and it was not regenerated # a different expected-return model shifts the scenarios but leaves the shocks alone symbols = svc.candidate_symbols(OptimizationRequest(**base)) req = OptimizationRequest(**base) mu_a, _, _ = svc.expected_returns(req, symbols) scen_a = svc.scenarios(req, symbols, mu_a) req_views = OptimizationRequest(**dict(base, views=[{"symbol": SYMBOLS[0], "sentiment": 0.9, "confidence": 0.8}])) mu_b, _, _ = svc.expected_returns(req_views, symbols) scen_b = svc.scenarios(req_views, symbols, mu_b) assert scen_b.hash != scen_a.hash and scen_b.shock_hash == scen_a.shock_hash np.testing.assert_allclose(scen_b.shocks, scen_a.shocks, atol=1e-15) assert (scen_b.expected_simple[0] > scen_a.expected_simple[0]) # the view moved the drift # what *must* change the set: symbols, horizon, covariance model, seed, sampling, scenario count for label, kwargs in { "fewer symbols": dict(base, symbols=SYMBOLS[:8]), "monthly horizon": dict(base, horizon="month"), "other covariance": dict(base, covariance={"method": "rolling", "window": 250}), "other seed": dict(base, seed=7), "plain Monte Carlo": dict(base, sampling="mc"), "more scenarios": dict(base, n_scenarios=4096), "no tilt": dict(base, importance_sampling="off"), }.items(): other = OptimizationRequest(**kwargs) assert svc.scenario_key(other, svc.candidate_symbols(other))["id"] != ref.scenario_key["id"], label # the adaptive tilt is honest about depending on the risk parameters ad = OptimizationRequest(**dict(base, importance_sampling="adaptive")) ad_other = OptimizationRequest(**dict(base, importance_sampling="adaptive", cvar_limit=6_000.0)) assert svc.scenario_key(ad, symbols)["id"] != svc.scenario_key(ad_other, symbols)["id"] assert "adaptive_tilt" in svc.scenario_key(ad, symbols) def test_tilt_is_driven_by_the_covariance_and_reports_its_alignment(svc): """The tilt shifts each asset by its covariance with the reference portfolio (not by its own volatility), leaves the model correlations exactly intact, and the result says how well the solution's loss lines up with the tilted direction.""" from scipy.stats import norm from quantcore.models import sampling as smp from quantcore.models.covariance import CovarianceModel from quantcore.models.simulation import GBMSimulator, TiltSpec from quantcore.optimization.problem import OptimizationConstraints, ProblemSpec, TransactionCosts req = OptimizationRequest(**dict(BASE, scenario_model="gbm", n_scenarios=2048, importance_sampling="market")) symbols = svc.candidate_symbols(req) r, cov_model = svc._returns_and_cov_model(req, symbols) h = req.horizon_periods cov_h = cov_model.estimate(r).values * h root = smp.pca_root(cov_h) reference = svc.market_exposure(symbols) * 1e6 # 1. the shift of a tilted block is exactly -z_beta * Sigma v / sigma: proportional to each # asset's covariance with the reference, i.e. to its beta, not to its own volatility z = norm.ppf(0.95) direction = smp.loss_direction(root, reference) shift = z * (root @ direction) sigma = np.sqrt(reference @ cov_h @ reference) np.testing.assert_allclose(shift, -z * (cov_h @ reference) / sigma, atol=1e-12) beta_i = (cov_h @ reference) / (reference @ cov_h @ reference) assert np.corrcoef(shift, beta_i)[0, 1] < -0.999 # perfectly ordered by beta if np.std(np.sqrt(np.diag(cov_h))) > 0: # ... and only loosely by own volatility assert np.corrcoef(shift, np.sqrt(np.diag(cov_h)))[0, 1] > -0.999 # 2. the weighted tilted sample has exactly the model covariance (and correlations) mu_ann, _, _ = svc.expected_returns(req, symbols) sim = GBMSimulator(cov_model) tilted = sim.simulate(r, n_scenarios=2048, horizon_periods=h, mu=mu_ann / 252, seed=1, tilt=TiltSpec(exposure=reference)) w = tilted.weights centred = tilted.log_returns - w @ tilted.log_returns cov_t = (centred * w[:, None]).T @ centred np.testing.assert_allclose(cov_t, cov_h, rtol=1e-8, atol=1e-14) unweighted = np.cov(tilted.log_returns, rowvar=False) assert np.abs(unweighted - cov_h).max() > 1e-6 # the raw tilted sample is *not* the model # 3. alignment: 1 for the reference itself, negative for a portfolio built to oppose it spec = ProblemSpec(scenarios=tilted.simple_returns(), probabilities=tilted.probabilities, expected_returns=tilted.expected_simple, **svc._spec_common(req, symbols)) assert abs(svc.tilt_alignment(spec, reference, svc.market_exposure(symbols)) - 1) < 1e-9 # itself opposed = np.zeros(len(symbols)) opposed[np.argsort(beta_i)[:2]] = 1e5 # the two assets least exposed to the reference assert svc.tilt_alignment(spec, opposed, svc.market_exposure(symbols)) < 0.9 # a real run reports the figure and stays well above the warning threshold res = svc.optimize(req) assert res.solver_info["tilt_alignment"] > 0.5 and not res.warnings def test_risk_drift_zero_keeps_the_limit_binding_when_the_forecast_is_large(md, settings): """Over a long horizon a big expected return shifts the whole tail into profit, and a CVaR limit measured on the drifted scenarios stops constraining anything. The zero-drift convention measures risk on the centred scenarios -- so the risk figure does not depend on the return forecast at all -- while the objective still uses it.""" svc = OptimizationService(settings, md) base = dict(symbols=SYMBOLS[:8], cash=1_000_000.0, horizon=252, freq="daily", confidence=0.99, cvar_limit_fraction=0.12, scenario_model="gbm", n_scenarios=2048, solver="highs", expected_returns="zero", row_generation="off", importance_sampling="off") symbols = svc.candidate_symbols(OptimizationRequest(**base)) v = np.full(len(symbols), 125_000.0) # an equally weighted portfolio of the whole wealth risk = {} for mu_val in (0.8, 1.5): # 80 % and 150 % p.a. forecasts for mode in ("model", "zero"): req = OptimizationRequest(**{**base, "risk_drift": mode}) spec = svc.build_spec(req, symbols, svc.scenarios(req, symbols, pd.Series(mu_val, index=symbols))) risk[(mu_val, mode)] = OptimizationService.portfolio_risk_on_scenarios(spec, v / spec.prices, 0.99) if mode == "zero": np.testing.assert_allclose(spec.mean_returns, svc.build_spec( OptimizationRequest(**{**base, "risk_drift": "model"}), symbols, svc.scenarios(OptimizationRequest(**base), symbols, pd.Series(mu_val, index=symbols)) ).mean_returns) # the objective sees the forecast either way # the zero-drift risk is the same whatever the forecast; the drifted one is not assert risk[(0.8, "zero")]["cvar"] == pytest.approx(risk[(1.5, "zero")]["cvar"], rel=1e-12) assert risk[(1.5, "model")]["cvar"] < 0 < risk[(1.5, "zero")]["cvar"] # the tail became a gain assert risk[(0.8, "model")]["cvar"] < risk[(0.8, "zero")]["cvar"] # the drift enters multiplicatively (e^{r+d} - 1), so the gap is well below mu * V gap = risk[(1.5, "zero")]["cvar"] - risk[(1.5, "model")]["cvar"] assert 0.3 * 1.5e6 < gap < 1.5e6 res = svc.optimize(OptimizationRequest(**{**base, "expected_returns": "historical"})) if res.cvar <= 0: # the historical forecast is large enough to make the limit vacuous assert any("constrains nothing" in w for w in res.warnings) def test_an_all_cash_utility_solution_says_why(md, settings): """objective='utility' with a large lambda is optimal in cash -- which must not look like a bug.""" from quantcore.optimization.service import OptimizationRequest, OptimizationService svc = OptimizationService(settings, md) common = dict(cash=1_000_000.0, max_assets=6, horizon="week", scenario_model="gaussian_exact", solver="highs", time_limit=120, objective="utility") big = svc.optimize(OptimizationRequest(risk_aversion=50.0, **common)) assert big.status in ("Optimal", "FeasibleFound") and big.invested == 0 assert any("nothing was invested" in w and "risk_aversion" in w for w in big.warnings) small = svc.optimize(OptimizationRequest(risk_aversion=0.0, **common)) assert small.invested > 0 and not any("nothing was invested" in w for w in small.warnings) def test_a_weighted_cvar_objective_is_routed_to_highs_and_says_so(md, settings): """cuOpt aborts when the CVaR is the objective *and* the scenarios carry probabilities. The routing is that narrow on purpose: the same objective with uniform weights solves on the GPU, so only the weighted case may be diverted, and only under solver='auto'. """ from quantcore.optimization.service import OptimizationRequest, OptimizationService svc = OptimizationService(settings, md) common = dict(cash=1_000_000.0, max_assets=6, horizon="week", scenario_model="gbm", n_scenarios=1024, solver="auto", time_limit=120, risk_drift="zero") weighted = svc.optimize(OptimizationRequest(objective="min_cvar", target_return=0.0005, importance_sampling="market", **common)) assert weighted.solver == "highs" assert any("probability-weighted" in w for w in weighted.warnings) # uniform weights: no diversion (cuOpt is unreachable in the tests, so it still lands on HiGHS -- # through the fallback, which says something different) plain = svc.optimize(OptimizationRequest(objective="min_cvar", target_return=0.0005, importance_sampling="off", **common)) assert not any("probability-weighted" in w for w in plain.warnings) # the constrained form is never diverted either other = svc.optimize(OptimizationRequest(objective="max_return", cvar_limit_fraction=0.05, importance_sampling="market", **common)) assert not any("probability-weighted" in w for w in other.warnings)