From c15be189de05b1d4a2b62fa9aa26954896817e7f Mon Sep 17 00:00:00 2001 From: fall-development-rob Date: Fri, 17 Apr 2026 09:37:05 +0100 Subject: [PATCH] feat: v0.4 stochastic foundation (multi-asset, fat tails, AR(1) inflation, Gompertz longevity, risk metrics) Integrates stochastic samplers into the MC loop while preserving byte-identical output for v0.3-shape scenarios. Adds multi-asset correlated returns (LogNormal, Student-T copula, Bootstrap), AR(1) inflation with calibration presets, Gompertz and cohort longevity models, and institutional risk metrics (VaR/CVaR/Sortino/ drawdown). Sensitivity and backtest modules force-disable stochastic features per ADR-027/031. Version bumped to 0.4.0. 68 smoke tests pass (46 new). Co-Authored-By: claude-flow --- package.json | 8 +- scripts/fix-esm-extensions.mjs | 38 ++ src/__tests__/smoke.test.mjs | 534 +++++++++++++++ src/__tests__/v03-fixture.json | 27 + src/advanced.ts | 5 + src/backtest.ts | 16 +- src/data/life-tables/uk-ons-2020.json | 378 +++++++++++ src/data/life-tables/us-ssa-2020.json | 378 +++++++++++ src/data/shiller-1871-2024.json | 931 ++++++++++++++++++++++++++ src/index.ts | 32 + src/inflation-sampler.ts | 89 +++ src/longevity-sampler.ts | 195 ++++++ src/monte-carlo.ts | 248 ++++++- src/projection.ts | 37 +- src/return-sampler.ts | 672 +++++++++++++++++++ src/risk-metrics.ts | 236 +++++++ src/sensitivity.ts | 22 +- src/types.ts | 180 +++++ 18 files changed, 4003 insertions(+), 23 deletions(-) create mode 100644 scripts/fix-esm-extensions.mjs create mode 100644 src/__tests__/v03-fixture.json create mode 100644 src/data/life-tables/uk-ons-2020.json create mode 100644 src/data/life-tables/us-ssa-2020.json create mode 100644 src/data/shiller-1871-2024.json create mode 100644 src/inflation-sampler.ts create mode 100644 src/longevity-sampler.ts create mode 100644 src/return-sampler.ts create mode 100644 src/risk-metrics.ts diff --git a/package.json b/package.json index e82ca5e..d132aa7 100644 --- a/package.json +++ b/package.json @@ -1,8 +1,9 @@ { "name": "@robotixai/calculator-engine", - "version": "0.3.0", + "version": "0.4.0", "description": "Financial retirement projection engine with Monte Carlo simulation, multi-jurisdiction tax, and withdrawal strategies", "private": false, + "type": "module", "main": "./dist/index.js", "types": "./dist/index.d.ts", "exports": { @@ -14,8 +15,9 @@ }, "files": ["dist", "README.md", "LICENSE"], "scripts": { - "build": "tsc", - "lint": "tsc --noEmit" + "build": "tsc && node scripts/fix-esm-extensions.mjs", + "lint": "tsc --noEmit", + "test": "node src/__tests__/smoke.test.mjs" }, "dependencies": {}, "devDependencies": { diff --git a/scripts/fix-esm-extensions.mjs b/scripts/fix-esm-extensions.mjs new file mode 100644 index 0000000..01a0384 --- /dev/null +++ b/scripts/fix-esm-extensions.mjs @@ -0,0 +1,38 @@ +/** + * Post-build script: adds .js extensions to relative imports in dist/ so that + * Node's native ESM loader can resolve them without --experimental-specifier-resolution. + */ +import { readdir, readFile, writeFile } from 'node:fs/promises'; +import { join } from 'node:path'; + +const DIST = new URL('../dist', import.meta.url).pathname; + +// Match: from './foo' or from '../foo' (no extension) +const RE = /(from\s+['"])(\.\.?\/[^'"]+?)(['"])/g; + +async function walk(dir) { + const entries = await readdir(dir, { withFileTypes: true }); + const files = []; + for (const e of entries) { + const full = join(dir, e.name); + if (e.isDirectory()) files.push(...(await walk(full))); + else if (e.name.endsWith('.js')) files.push(full); + } + return files; +} + +const files = await walk(DIST); +let patched = 0; + +for (const f of files) { + const src = await readFile(f, 'utf8'); + const out = src.replace(RE, (match, pre, path, post) => { + // Skip if already has an extension + if (/\.\w+$/.test(path)) return match; + patched++; + return `${pre}${path}.js${post}`; + }); + if (out !== src) await writeFile(f, out); +} + +console.log(`fix-esm-extensions: patched ${patched} imports in ${files.length} files`); diff --git a/src/__tests__/smoke.test.mjs b/src/__tests__/smoke.test.mjs index 41d23fc..436d21c 100644 --- a/src/__tests__/smoke.test.mjs +++ b/src/__tests__/smoke.test.mjs @@ -10,11 +10,28 @@ import { DEFAULT_SCENARIO, runProjection, + runMonteCarloSimulation, runSensitivityAnalysis, findRequiredSavings, calculateFixedPctWithdrawal, + buildReturnSampler, + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + buildInflationSampler, + buildLongevitySampler, + computeRiskMetrics, + INFLATION_PRESETS, } from '../../dist/index.js'; +import { readFileSync } from 'node:fs'; +import { fileURLToPath } from 'node:url'; +import { dirname, join } from 'node:path'; + +const __dirname = dirname(fileURLToPath(import.meta.url)); +const v03Fixture = JSON.parse( + readFileSync(join(__dirname, 'v03-fixture.json'), 'utf8'), +); + let passed = 0; let failed = 0; @@ -241,6 +258,523 @@ console.log('\nTest 5: Fixed-Pct strategy in runProjection'); ); } +// --------------------------------------------------------------------------- +// Test 6: buildReturnSampler — seed determinism +// --------------------------------------------------------------------------- +console.log('\nTest 6: buildReturnSampler seed determinism'); +{ + const sampler1 = buildReturnSampler( + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + { kind: 'LogNormal' }, + 12345, + ); + const sampler2 = buildReturnSampler( + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + { kind: 'LogNormal' }, + 12345, + ); + const draw1 = sampler1.sample(0); + const draw2 = sampler2.sample(0); + assert( + draw1.us_equity === draw2.us_equity, + `same seed => same us_equity draw (${draw1.us_equity})`, + ); + assert( + draw1.us_bond === draw2.us_bond, + `same seed => same us_bond draw`, + ); + + // Different seed => different draw + const sampler3 = buildReturnSampler( + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + { kind: 'LogNormal' }, + 99999, + ); + const draw3 = sampler3.sample(0); + assert( + draw3.us_equity !== draw1.us_equity, + 'different seed => different draw', + ); +} + +// --------------------------------------------------------------------------- +// Test 7: Cholesky on known PSD matrix +// --------------------------------------------------------------------------- +console.log('\nTest 7: Cholesky basic validity'); +{ + // If buildReturnSampler doesn't throw, Cholesky succeeded on the + // DEFAULT_CORRELATIONS matrix (a 5x5 PSD matrix). + let threw = false; + try { + buildReturnSampler( + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + { kind: 'LogNormal' }, + 42, + ); + } catch { + threw = true; + } + assert(!threw, 'Cholesky succeeds on DEFAULT_CORRELATIONS (PSD)'); + + // 2x2 identity matrix should also work + const twoAssets = [ + { id: 'us_equity', name: 'A', expected_return_pct: 10, return_stdev_pct: 15, weight_pct: 50 }, + { id: 'us_bond', name: 'B', expected_return_pct: 5, return_stdev_pct: 5, weight_pct: 50 }, + ]; + const identityCorr = { + ids: ['us_equity', 'us_bond'], + values: [[1, 0], [0, 1]], + }; + const s = buildReturnSampler(twoAssets, identityCorr, { kind: 'LogNormal' }, 1); + const d = s.sample(0); + assert(typeof d.us_equity === 'number' && isFinite(d.us_equity), 'identity corr produces finite draws'); +} + +// --------------------------------------------------------------------------- +// Test 8: buildInflationSampler — AR(1) mean-reversion +// --------------------------------------------------------------------------- +console.log('\nTest 8: AR(1) inflation mean-reversion'); +{ + const process = { + kind: 'AR1', + long_run_mean_pct: 3.0, + phi: 0.6, + shock_stdev_pct: 1.5, + initial_pct: 3.0, + }; + const sampler = buildInflationSampler(process, 42); + assert(sampler.kind === 'AR1', 'sampler kind is AR1'); + + // Sample 500 years and check mean converges to long_run_mean + let sum = 0; + let prior = process.initial_pct / 100; + const N = 500; + for (let y = 0; y < N; y++) { + const r = sampler.sample(y, prior); + sum += r; + prior = r; + } + const sampleMean = sum / N; + const target = process.long_run_mean_pct / 100; + const pctDiff = Math.abs(sampleMean - target) / target; + assert( + pctDiff < 0.15, + `AR(1) 500-year mean ${(sampleMean * 100).toFixed(2)}% within 15% of ${process.long_run_mean_pct}%`, + ); +} + +// --------------------------------------------------------------------------- +// Test 9: buildInflationSampler — Flat passthrough +// --------------------------------------------------------------------------- +console.log('\nTest 9: Flat inflation passthrough'); +{ + const sampler = buildInflationSampler({ kind: 'Flat', rate_pct: 2.5 }, 42); + assert(sampler.kind === 'Flat', 'sampler kind is Flat'); + assert( + sampler.sample(0, 0) === 0.025, + 'Flat sampler returns configured rate', + ); + assert( + sampler.sample(99, 0.1) === 0.025, + 'Flat sampler ignores prior inflation', + ); +} + +// --------------------------------------------------------------------------- +// Test 10: buildLongevitySampler — Gompertz median +// --------------------------------------------------------------------------- +console.log('\nTest 10: Gompertz longevity median'); +{ + const modal = 88; + const dispersion = 10; + const sampler = buildLongevitySampler( + { kind: 'Gompertz', modal_age: modal, dispersion }, + 42, + ); + assert(sampler.kind === 'Gompertz', 'sampler kind is Gompertz'); + + const med = sampler.median(65); + // Conditional median for a 65yo should be somewhere around 85-92 + assert(med > 80 && med < 100, `Gompertz median for age 65 is ${med.toFixed(1)} (in 80-100)`); + + // Unconditional median formula: modal + dispersion * ln(ln(2)) ≈ 88 + 10*(-0.366) ≈ 84.3 + // But conditional on age 65, it should be higher + assert(med > 84, 'conditional median > unconditional median'); +} + +// --------------------------------------------------------------------------- +// Test 11: buildLongevitySampler — Fixed +// --------------------------------------------------------------------------- +console.log('\nTest 11: Fixed longevity'); +{ + const sampler = buildLongevitySampler({ kind: 'Fixed', end_age: 95 }, 42); + assert(sampler.kind === 'Fixed', 'sampler kind is Fixed'); + assert(sampler.sample(65) === 95, 'Fixed always returns end_age'); + assert(sampler.median(65) === 95, 'Fixed median is end_age'); + assert(sampler.survival(94, 65) === 1, 'survival at 94 is 1'); + assert(sampler.survival(96, 65) === 0, 'survival at 96 is 0'); +} + +// --------------------------------------------------------------------------- +// Test 12: computeRiskMetrics — CVaR <= VaR +// --------------------------------------------------------------------------- +console.log('\nTest 12: computeRiskMetrics invariants'); +{ + // Build synthetic MC data + const terminals = []; + const paths = []; + const returns = []; + for (let i = 0; i < 300; i++) { + const t = 500000 + (i - 150) * 5000 + Math.sin(i) * 100000; + terminals.push(t); + const path = []; + for (let y = 0; y < 30; y++) { + path.push(500000 + (y * 10000) + (i - 150) * 1000); + } + paths.push(path); + returns.push(0.04 + (i - 150) * 0.001); + } + + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 35, + end_age: 65, + current_balance: 500000, + }; + + const metrics = computeRiskMetrics( + { terminal_distribution: terminals, real_balance_paths: paths, annualised_returns: returns }, + scenario, + ); + + assert( + metrics.cvar_95_terminal_real <= metrics.var_95_terminal_real, + `CVaR95 (${metrics.cvar_95_terminal_real.toFixed(0)}) <= VaR95 (${metrics.var_95_terminal_real.toFixed(0)})`, + ); + assert( + metrics.cvar_99_terminal_real <= metrics.var_99_terminal_real, + `CVaR99 (${metrics.cvar_99_terminal_real.toFixed(0)}) <= VaR99 (${metrics.var_99_terminal_real.toFixed(0)})`, + ); + assert( + metrics.max_drawdown_pct >= 0 && metrics.max_drawdown_pct <= 100, + `max_drawdown_pct in [0,100]: ${metrics.max_drawdown_pct.toFixed(2)}`, + ); + + // Per-year P10 <= P50 <= P90 + for (let y = 0; y < metrics.p10_year_by_year_balance_real.length; y++) { + const p10 = metrics.p10_year_by_year_balance_real[y]; + const p50 = metrics.p50_year_by_year_balance_real[y]; + const p90 = metrics.p90_year_by_year_balance_real[y]; + if (p10 > p50 + 0.01 || p50 > p90 + 0.01) { + assert(false, `P10 <= P50 <= P90 violated at year ${y}`); + break; + } + } + assert(true, 'P10 <= P50 <= P90 holds for all years'); +} + +// --------------------------------------------------------------------------- +// Test 13: Backwards compat — v0.3 scenario deterministic projection +// --------------------------------------------------------------------------- +console.log('\nTest 13: v0.3 backwards compat — deterministic projection'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 30, + retirement_age: 65, + end_age: 90, + current_balance: 500_000, + contrib_amount: 10_000, + contrib_cadence: 'Annual', + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: false, + black_swan_enabled: false, + }; + + const { metrics } = runProjection(scenario); + const tolerance = 0.01; // 1 cent + + assert( + Math.abs(metrics.terminal_real - v03Fixture.detRunMetrics.terminal_real) < tolerance, + `det terminal_real matches fixture (${metrics.terminal_real.toFixed(2)} vs ${v03Fixture.detRunMetrics.terminal_real.toFixed(2)})`, + ); + assert( + Math.abs(metrics.terminal_nominal - v03Fixture.detRunMetrics.terminal_nominal) < tolerance, + `det terminal_nominal matches fixture`, + ); +} + +// --------------------------------------------------------------------------- +// Test 14: Backwards compat — v0.3 MC terminal values +// --------------------------------------------------------------------------- +console.log('\nTest 14: v0.3 backwards compat — MC terminal values'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 30, + retirement_age: 65, + end_age: 90, + current_balance: 500_000, + contrib_amount: 10_000, + contrib_cadence: 'Annual', + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: true, + mc_runs: 500, + black_swan_enabled: false, + }; + + const mc = runMonteCarloSimulation(scenario, runProjection, { + runs: 500, + seed: 42, + }); + + assert(mc.runs_completed === 500, `MC completed 500 runs`); + assert( + Math.abs(mc.probability_no_shortfall - v03Fixture.mc.probability_no_shortfall) < 0.01, + `MC prob_no_shortfall matches fixture (${mc.probability_no_shortfall})`, + ); + + // Check first 5 terminal values match fixture + const first5 = mc.terminal_distribution.slice(0, 5); + let allMatch = true; + for (let i = 0; i < 5; i++) { + if (Math.abs(first5[i] - v03Fixture.mc.first5Terminals[i]) > 0.01) { + allMatch = false; + console.error(` terminal[${i}] mismatch: got ${first5[i]}, expected ${v03Fixture.mc.first5Terminals[i]}`); + } + } + assert(allMatch, 'first 5 MC terminals match v0.3 fixture (byte-identical legacy path)'); +} + +// --------------------------------------------------------------------------- +// Test 15: MC with multi-asset returns +// --------------------------------------------------------------------------- +console.log('\nTest 15: MC with multi-asset returns'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 60, + retirement_age: 65, + end_age: 85, + current_balance: 1_000_000, + contrib_amount: 0, + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: true, + mc_runs: 200, + black_swan_enabled: false, + asset_classes: DEFAULT_ASSET_CLASSES, + return_correlation_matrix: DEFAULT_CORRELATIONS, + return_distribution_kind: 'LogNormal', + }; + + const mc = runMonteCarloSimulation(scenario, runProjection, { + runs: 200, + seed: 42, + }); + + assert(mc.runs_completed === 200, 'multi-asset MC completed 200 runs'); + assert(mc.fan_chart.length > 0, 'multi-asset MC produces fan chart'); + assert( + mc.risk_metrics != null, + 'risk_metrics present when runs >= 200 and horizon >= 10', + ); +} + +// --------------------------------------------------------------------------- +// Test 16: MC with AR(1) inflation +// --------------------------------------------------------------------------- +console.log('\nTest 16: MC with AR(1) inflation'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 55, + retirement_age: 65, + end_age: 85, + current_balance: 1_000_000, + contrib_amount: 0, + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: true, + mc_runs: 200, + black_swan_enabled: false, + inflation_model: 'AR1', + inflation_long_run_mean_pct: 3.0, + inflation_ar1_phi: 0.6, + inflation_shock_stdev_pct: 1.5, + inflation_initial_pct: 2.5, + }; + + const mc = runMonteCarloSimulation(scenario, runProjection, { + runs: 200, + seed: 42, + }); + + assert(mc.runs_completed === 200, 'AR(1) MC completed 200 runs'); + assert( + mc.inflation_fan_chart != null && mc.inflation_fan_chart.length > 0, + 'inflation_fan_chart present for AR(1)', + ); +} + +// --------------------------------------------------------------------------- +// Test 17: MC with Gompertz longevity +// --------------------------------------------------------------------------- +console.log('\nTest 17: MC with Gompertz longevity'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 60, + retirement_age: 65, + end_age: 100, + current_balance: 1_000_000, + contrib_amount: 0, + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: true, + mc_runs: 200, + black_swan_enabled: false, + longevity_model: 'Gompertz', + longevity_modal_age: 88, + longevity_dispersion: 10, + }; + + const mc = runMonteCarloSimulation(scenario, runProjection, { + runs: 200, + seed: 42, + }); + + assert(mc.runs_completed === 200, 'Gompertz MC completed 200 runs'); + assert( + mc.lifespan_distribution != null && mc.lifespan_distribution.length === 200, + 'lifespan_distribution has 200 entries', + ); + // Median sampled death age should be plausible + const sorted = [...mc.lifespan_distribution].sort((a, b) => a - b); + const medianDeath = sorted[100]; + assert( + medianDeath >= 80 && medianDeath <= 100, + `median sampled death age ${medianDeath} is plausible (80-100)`, + ); +} + +// --------------------------------------------------------------------------- +// Test 18: DEFAULT_ASSET_CLASSES weights sum to 100 +// --------------------------------------------------------------------------- +console.log('\nTest 18: DEFAULT_ASSET_CLASSES'); +{ + const totalWeight = DEFAULT_ASSET_CLASSES.reduce((s, a) => s + a.weight_pct, 0); + assert( + Math.abs(totalWeight - 100) < 1, + `DEFAULT_ASSET_CLASSES weights sum to ~100 (${totalWeight})`, + ); + assert(DEFAULT_ASSET_CLASSES.length === 5, '5 default asset classes'); +} + +// --------------------------------------------------------------------------- +// Test 19: INFLATION_PRESETS are exported +// --------------------------------------------------------------------------- +console.log('\nTest 19: INFLATION_PRESETS exported'); +{ + assert(INFLATION_PRESETS != null, 'INFLATION_PRESETS is defined'); + assert( + typeof INFLATION_PRESETS['US-CPI'] === 'object', + 'US-CPI preset exists', + ); + assert( + INFLATION_PRESETS['US-CPI'].phi > 0, + 'US-CPI has positive phi', + ); +} + +// --------------------------------------------------------------------------- +// Test 20: Student-T sampler produces finite draws +// --------------------------------------------------------------------------- +console.log('\nTest 20: Student-T sampler'); +{ + const sampler = buildReturnSampler( + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + { kind: 'StudentT', dof: 5 }, + 42, + ); + const draw = sampler.sample(0); + assert(isFinite(draw.us_equity), `Student-T us_equity is finite (${draw.us_equity.toFixed(4)})`); + assert(isFinite(draw.us_bond), `Student-T us_bond is finite (${draw.us_bond.toFixed(4)})`); +} + +// --------------------------------------------------------------------------- +// Test 21: Risk metrics on MC result +// --------------------------------------------------------------------------- +console.log('\nTest 21: risk_metrics on MC result'); +{ + const scenario = { + ...DEFAULT_SCENARIO, + current_age: 35, + retirement_age: 65, + end_age: 90, + current_balance: 500_000, + contrib_amount: 10_000, + contrib_cadence: 'Annual', + nominal_return_pct: 7, + return_stdev_pct: 12, + return_distribution: 'log-normal', + inflation_pct: 2.5, + fee_pct: 0.5, + enable_mc: true, + mc_runs: 200, + black_swan_enabled: false, + }; + + const mc = runMonteCarloSimulation(scenario, runProjection, { + runs: 200, + seed: 42, + }); + + assert(mc.risk_metrics != null, 'risk_metrics attached for 200 runs, 55yr horizon'); + if (mc.risk_metrics) { + assert( + mc.risk_metrics.cvar_95_terminal_real <= mc.risk_metrics.var_95_terminal_real, + 'CVaR95 <= VaR95 on real MC', + ); + assert( + mc.risk_metrics.p10_year_by_year_balance_real.length > 0, + 'year-by-year balance trajectories populated', + ); + } +} + +// --------------------------------------------------------------------------- +// Test 22 (new): Gompertz sampler determinism +// --------------------------------------------------------------------------- +console.log('\nTest 22: Gompertz sampler determinism'); +{ + const s1 = buildLongevitySampler({ kind: 'Gompertz', modal_age: 88, dispersion: 10 }, 42); + const s2 = buildLongevitySampler({ kind: 'Gompertz', modal_age: 88, dispersion: 10 }, 42); + assert(s1.sample(65) === s2.sample(65), 'same seed => same Gompertz sample'); +} + // --------------------------------------------------------------------------- // Summary // --------------------------------------------------------------------------- diff --git a/src/__tests__/v03-fixture.json b/src/__tests__/v03-fixture.json new file mode 100644 index 0000000..81bb4fb --- /dev/null +++ b/src/__tests__/v03-fixture.json @@ -0,0 +1,27 @@ +{ + "detRunMetrics": { + "terminal_nominal": 10560631.870534958, + "terminal_real": 2400258.301747078, + "first_shortfall_age": null, + "readiness_score": 100, + "total_contributions": 350000, + "total_withdrawals": 8156353.487732654, + "total_fees": 1399167.5842914192, + "total_taxes": 0, + "estate_value": 10560631.870534958 + }, + "mc": { + "probability_no_shortfall": 100, + "median_terminal": 1655566.1159770908, + "p10_terminal": 598987.3670926923, + "p90_terminal": 4878362.525935733, + "runs_completed": 500, + "first5Terminals": [ + 850272.1334223117, + 1688225.5733266973, + 7018216.051185282, + 1072436.452260471, + 1316040.4977310624 + ] + } +} \ No newline at end of file diff --git a/src/advanced.ts b/src/advanced.ts index 3f2506f..d9e5928 100644 --- a/src/advanced.ts +++ b/src/advanced.ts @@ -821,6 +821,11 @@ export function runAdvancedProjection( shortfall_withdrawals: shortfallWithdrawals, black_swan_loss: yearBlackSwanLoss, withdrawal_event: withdrawalEvent, + // v0.4 additions: deterministic advanced projection still uses the flat + // legacy inflation rate; asset_returns is null because per-class + // multi-asset mode is not active in this code path. + inflation_this_year: inflation_enabled ? inflation_pct / 100 : 0, + asset_returns: null, }; timeline.push(row); diff --git a/src/backtest.ts b/src/backtest.ts index bc82659..5bd8695 100644 --- a/src/backtest.ts +++ b/src/backtest.ts @@ -202,7 +202,19 @@ export function runHistoricalBacktest( projectionFn: (s: Scenario, returns: number[]) => { metrics: Metrics }, ): BacktestResult { const log = getLogger(); - const span = scenario.end_age - scenario.current_age; + + // v0.4 (ADR-027/031): force-disable stochastic features for backtests. + // Historical windows already contain real-world volatility and inflation; + // layering stochastic samplers would double-count. + const baseScenario: Scenario = { + ...scenario, + inflation_model: 'Flat' as const, + longevity_model: 'Fixed' as const, + return_distribution_kind: 'LogNormal' as const, + asset_classes: [], + }; + + const span = baseScenario.end_age - baseScenario.current_age; // Guard: span must be at least 1 if (span < 1) { @@ -252,7 +264,7 @@ export function runHistoricalBacktest( // ADR-027: backtests already include real historical crashes — layering a // synthetic Black Swan event would double-count. Force-disable it on a // per-window scenario clone. - const periodScenario: Scenario = JSON.parse(JSON.stringify(scenario)); + const periodScenario: Scenario = JSON.parse(JSON.stringify(baseScenario)); periodScenario.black_swan_enabled = false; // Run projection with historical returns diff --git a/src/data/life-tables/uk-ons-2020.json b/src/data/life-tables/uk-ons-2020.json new file mode 100644 index 0000000..8b48bc1 --- /dev/null +++ b/src/data/life-tables/uk-ons-2020.json @@ -0,0 +1,378 @@ +{ + "_source": "synthetic placeholder; Gompertz approximation of UK-ONS-2020-cohort cohort table; replace with authoritative cohort table before production", + "_country": "UK-ONS-2020-cohort", + "_modal_age_M": 86, + "_modal_age_F": 89, + "_dispersion": 10, + "ages": [ + 0, + 1, + 2, + 3, + 4, + 5, + 6, + 7, + 8, + 9, + 10, + 11, + 12, + 13, + 14, + 15, + 16, + 17, + 18, + 19, + 20, + 21, + 22, + 23, + 24, + 25, + 26, + 27, + 28, + 29, + 30, + 31, + 32, + 33, + 34, + 35, + 36, + 37, + 38, + 39, + 40, + 41, + 42, + 43, + 44, + 45, + 46, + 47, + 48, + 49, + 50, + 51, + 52, + 53, + 54, + 55, + 56, + 57, + 58, + 59, + 60, + 61, + 62, + 63, + 64, + 65, + 66, + 67, + 68, + 69, + 70, + 71, + 72, + 73, + 74, + 75, + 76, + 77, + 78, + 79, + 80, + 81, + 82, + 83, + 84, + 85, + 86, + 87, + 88, + 89, + 90, + 91, + 92, + 93, + 94, + 95, + 96, + 97, + 98, + 99, + 100, + 101, + 102, + 103, + 104, + 105, + 106, + 107, + 108, + 109, + 110, + 111, + 112, + 113, + 114, + 115, + 116, + 117, + 118, + 119, + 120 + ], + "survival": { + "M": [ + 1, + 0.999981, + 0.999959, + 0.999936, + 0.999909, + 0.999881, + 0.999849, + 0.999813, + 0.999774, + 0.999731, + 0.999684, + 0.999631, + 0.999573, + 0.999509, + 0.999438, + 0.999359, + 0.999272, + 0.999177, + 0.999071, + 0.998954, + 0.998824, + 0.998682, + 0.998524, + 0.998349, + 0.998156, + 0.997943, + 0.997708, + 0.997448, + 0.997161, + 0.996843, + 0.996492, + 0.996105, + 0.995677, + 0.995204, + 0.994682, + 0.994105, + 0.993468, + 0.992764, + 0.991987, + 0.991128, + 0.990181, + 0.989135, + 0.98798, + 0.986705, + 0.985298, + 0.983745, + 0.982032, + 0.980142, + 0.978058, + 0.975759, + 0.973225, + 0.970433, + 0.967356, + 0.963966, + 0.960234, + 0.956126, + 0.951607, + 0.946637, + 0.941175, + 0.935175, + 0.928589, + 0.921363, + 0.913443, + 0.90477, + 0.89528, + 0.884907, + 0.873584, + 0.861238, + 0.847796, + 0.833185, + 0.81733, + 0.800158, + 0.781599, + 0.761589, + 0.74007, + 0.716995, + 0.692328, + 0.666053, + 0.638174, + 0.608717, + 0.577742, + 0.54534, + 0.511639, + 0.476811, + 0.441072, + 0.404682, + 0.367947, + 0.331215, + 0.294871, + 0.259325, + 0.225003, + 0.192331, + 0.161713, + 0.133511, + 0.108029, + 0.085485, + 0.066, + 0.049589, + 0.036155, + 0.025499, + 0.017335, + 0.011316, + 0.007063, + 0.004195, + 0.002359, + 0.001249, + 0.000618, + 0.000284, + 0.00012, + 0.000047, + 0.000016, + 0.000005, + 0.000001, + 0, + 0, + 0, + 0, + 0, + 0, + 0, + 0 + ], + "F": [ + 1, + 0.999986, + 0.99997, + 0.999952, + 0.999933, + 0.999912, + 0.999888, + 0.999862, + 0.999833, + 0.999801, + 0.999766, + 0.999727, + 0.999684, + 0.999636, + 0.999583, + 0.999525, + 0.999461, + 0.99939, + 0.999312, + 0.999225, + 0.999129, + 0.999023, + 0.998906, + 0.998777, + 0.998634, + 0.998476, + 0.998302, + 0.998109, + 0.997896, + 0.99766, + 0.9974, + 0.997113, + 0.996796, + 0.996445, + 0.996057, + 0.995629, + 0.995157, + 0.994634, + 0.994057, + 0.99342, + 0.992716, + 0.991939, + 0.991081, + 0.990134, + 0.989087, + 0.987932, + 0.986658, + 0.985251, + 0.983698, + 0.981985, + 0.980095, + 0.978011, + 0.975713, + 0.973179, + 0.970386, + 0.967309, + 0.96392, + 0.960188, + 0.956081, + 0.951562, + 0.946592, + 0.94113, + 0.935131, + 0.928544, + 0.921319, + 0.9134, + 0.904727, + 0.895237, + 0.884865, + 0.873542, + 0.861197, + 0.847756, + 0.833145, + 0.817291, + 0.80012, + 0.781562, + 0.761553, + 0.740035, + 0.71696, + 0.692295, + 0.666022, + 0.638143, + 0.608688, + 0.577715, + 0.545314, + 0.511615, + 0.476789, + 0.441051, + 0.404663, + 0.36793, + 0.331199, + 0.294857, + 0.259312, + 0.224992, + 0.192322, + 0.161705, + 0.133505, + 0.108024, + 0.085481, + 0.065997, + 0.049587, + 0.036154, + 0.025498, + 0.017334, + 0.011316, + 0.007063, + 0.004195, + 0.002359, + 0.001249, + 0.000618, + 0.000284, + 0.00012, + 0.000047, + 0.000016, + 0.000005, + 0.000001, + 0, + 0, + 0, + 0, + 0 + ] + } +} \ No newline at end of file diff --git a/src/data/life-tables/us-ssa-2020.json b/src/data/life-tables/us-ssa-2020.json new file mode 100644 index 0000000..e868914 --- /dev/null +++ b/src/data/life-tables/us-ssa-2020.json @@ -0,0 +1,378 @@ +{ + "_source": "synthetic placeholder; Gompertz approximation of US-SSA-2020-cohort cohort table; replace with authoritative cohort table before production", + "_country": "US-SSA-2020-cohort", + "_modal_age_M": 85, + "_modal_age_F": 88, + "_dispersion": 10, + "ages": [ + 0, + 1, + 2, + 3, + 4, + 5, + 6, + 7, + 8, + 9, + 10, + 11, + 12, + 13, + 14, + 15, + 16, + 17, + 18, + 19, + 20, + 21, + 22, + 23, + 24, + 25, + 26, + 27, + 28, + 29, + 30, + 31, + 32, + 33, + 34, + 35, + 36, + 37, + 38, + 39, + 40, + 41, + 42, + 43, + 44, + 45, + 46, + 47, + 48, + 49, + 50, + 51, + 52, + 53, + 54, + 55, + 56, + 57, + 58, + 59, + 60, + 61, + 62, + 63, + 64, + 65, + 66, + 67, + 68, + 69, + 70, + 71, + 72, + 73, + 74, + 75, + 76, + 77, + 78, + 79, + 80, + 81, + 82, + 83, + 84, + 85, + 86, + 87, + 88, + 89, + 90, + 91, + 92, + 93, + 94, + 95, + 96, + 97, + 98, + 99, + 100, + 101, + 102, + 103, + 104, + 105, + 106, + 107, + 108, + 109, + 110, + 111, + 112, + 113, + 114, + 115, + 116, + 117, + 118, + 119, + 120 + ], + "survival": { + "M": [ + 1, + 0.999979, + 0.999955, + 0.999929, + 0.9999, + 0.999868, + 0.999833, + 0.999794, + 0.999751, + 0.999703, + 0.99965, + 0.999592, + 0.999528, + 0.999457, + 0.999379, + 0.999292, + 0.999196, + 0.99909, + 0.998973, + 0.998844, + 0.998701, + 0.998543, + 0.998368, + 0.998176, + 0.997963, + 0.997727, + 0.997467, + 0.99718, + 0.996862, + 0.996512, + 0.996124, + 0.995696, + 0.995223, + 0.994701, + 0.994124, + 0.993487, + 0.992783, + 0.992006, + 0.991148, + 0.9902, + 0.989154, + 0.987999, + 0.986724, + 0.985317, + 0.983764, + 0.982051, + 0.980161, + 0.978077, + 0.975778, + 0.973244, + 0.970451, + 0.967374, + 0.963985, + 0.960253, + 0.956145, + 0.951626, + 0.946656, + 0.941193, + 0.935193, + 0.928607, + 0.921381, + 0.913461, + 0.904787, + 0.895297, + 0.884924, + 0.873601, + 0.861255, + 0.847813, + 0.833201, + 0.817346, + 0.800174, + 0.781615, + 0.761604, + 0.740085, + 0.717008, + 0.692341, + 0.666066, + 0.638186, + 0.608729, + 0.577753, + 0.54535, + 0.511649, + 0.476821, + 0.441081, + 0.40469, + 0.367954, + 0.331222, + 0.294876, + 0.25933, + 0.225008, + 0.192335, + 0.161716, + 0.133514, + 0.108031, + 0.085486, + 0.066001, + 0.04959, + 0.036156, + 0.0255, + 0.017336, + 0.011317, + 0.007063, + 0.004195, + 0.002359, + 0.001249, + 0.000618, + 0.000284, + 0.00012, + 0.000047, + 0.000016, + 0.000005, + 0.000001, + 0, + 0, + 0, + 0, + 0, + 0, + 0, + 0, + 0 + ], + "F": [ + 1, + 0.999984, + 0.999967, + 0.999947, + 0.999926, + 0.999902, + 0.999876, + 0.999847, + 0.999815, + 0.99978, + 0.999741, + 0.999698, + 0.99965, + 0.999598, + 0.99954, + 0.999475, + 0.999404, + 0.999326, + 0.999239, + 0.999143, + 0.999037, + 0.99892, + 0.998791, + 0.998648, + 0.99849, + 0.998316, + 0.998123, + 0.99791, + 0.997675, + 0.997415, + 0.997127, + 0.99681, + 0.996459, + 0.996072, + 0.995644, + 0.995171, + 0.994649, + 0.994072, + 0.993434, + 0.992731, + 0.991954, + 0.991095, + 0.990148, + 0.989102, + 0.987947, + 0.986672, + 0.985265, + 0.983712, + 0.981999, + 0.980109, + 0.978025, + 0.975727, + 0.973193, + 0.9704, + 0.967323, + 0.963934, + 0.960202, + 0.956095, + 0.951575, + 0.946606, + 0.941144, + 0.935144, + 0.928558, + 0.921333, + 0.913413, + 0.90474, + 0.89525, + 0.884878, + 0.873555, + 0.861209, + 0.847768, + 0.833157, + 0.817303, + 0.800131, + 0.781573, + 0.761564, + 0.740046, + 0.716971, + 0.692305, + 0.666031, + 0.638152, + 0.608697, + 0.577723, + 0.545321, + 0.511622, + 0.476796, + 0.441058, + 0.404669, + 0.367935, + 0.331204, + 0.294861, + 0.259316, + 0.224996, + 0.192325, + 0.161707, + 0.133507, + 0.108025, + 0.085482, + 0.065998, + 0.049588, + 0.036154, + 0.025498, + 0.017335, + 0.011316, + 0.007063, + 0.004195, + 0.002359, + 0.001249, + 0.000618, + 0.000284, + 0.00012, + 0.000047, + 0.000016, + 0.000005, + 0.000001, + 0, + 0, + 0, + 0, + 0, + 0 + ] + } +} \ No newline at end of file diff --git a/src/data/shiller-1871-2024.json b/src/data/shiller-1871-2024.json new file mode 100644 index 0000000..357ca2f --- /dev/null +++ b/src/data/shiller-1871-2024.json @@ -0,0 +1,931 @@ +{ + "_source": "us_equity from Shiller real total stock return; us_bond and us_cpi are synthetic placeholders generated with a seeded RNG to preserve realistic mean/stdev/correlation. Replace with authoritative Ibbotson + BLS CPI data before production.", + "_first_year": 1871, + "_last_year": 2024, + "rows": [ + { + "year": 1871, + "us_equity": 0.1478, + "us_bond": 0.0494, + "us_cpi": 0.0183 + }, + { + "year": 1872, + "us_equity": 0.1076, + "us_bond": 0.0498, + "us_cpi": 0.0097 + }, + { + "year": 1873, + "us_equity": -0.0415, + "us_bond": 0.005, + "us_cpi": -0.0033 + }, + { + "year": 1874, + "us_equity": 0.1249, + "us_bond": 0.0511, + "us_cpi": 0.0187 + }, + { + "year": 1875, + "us_equity": 0.0555, + "us_bond": 0.0965, + "us_cpi": 0.0399 + }, + { + "year": 1876, + "us_equity": -0.029, + "us_bond": 0.0513, + "us_cpi": 0.0423 + }, + { + "year": 1877, + "us_equity": 0.0187, + "us_bond": 0.0143, + "us_cpi": 0.0203 + }, + { + "year": 1878, + "us_equity": 0.1797, + "us_bond": 0.0276, + "us_cpi": 0.0335 + }, + { + "year": 1879, + "us_equity": 0.2102, + "us_bond": -0.0174, + "us_cpi": 0.0382 + }, + { + "year": 1880, + "us_equity": 0.2361, + "us_bond": 0.021, + "us_cpi": 0.0299 + }, + { + "year": 1881, + "us_equity": 0.0118, + "us_bond": 0.0246, + "us_cpi": 0.0269 + }, + { + "year": 1882, + "us_equity": 0.0043, + "us_bond": 0.0261, + "us_cpi": 0.0309 + }, + { + "year": 1883, + "us_equity": -0.0041, + "us_bond": 0.0313, + "us_cpi": 0.0304 + }, + { + "year": 1884, + "us_equity": -0.0785, + "us_bond": 0.0569, + "us_cpi": 0.0451 + }, + { + "year": 1885, + "us_equity": 0.3034, + "us_bond": 0.0324, + "us_cpi": 0.0225 + }, + { + "year": 1886, + "us_equity": 0.135, + "us_bond": -0.0152, + "us_cpi": 0.0404 + }, + { + "year": 1887, + "us_equity": -0.0122, + "us_bond": 0.0881, + "us_cpi": 0.0346 + }, + { + "year": 1888, + "us_equity": 0.0447, + "us_bond": 0.0863, + "us_cpi": 0.0193 + }, + { + "year": 1889, + "us_equity": 0.0511, + "us_bond": 0.0272, + "us_cpi": 0.0418 + }, + { + "year": 1890, + "us_equity": -0.0678, + "us_bond": 0.0405, + "us_cpi": 0.0424 + }, + { + "year": 1891, + "us_equity": 0.192, + "us_bond": 0.067, + "us_cpi": 0.019 + }, + { + "year": 1892, + "us_equity": 0.0486, + "us_bond": -0.0007, + "us_cpi": 0.0058 + }, + { + "year": 1893, + "us_equity": -0.1716, + "us_bond": 0.0502, + "us_cpi": -0.0033 + }, + { + "year": 1894, + "us_equity": 0.0311, + "us_bond": 0.0342, + "us_cpi": 0.0101 + }, + { + "year": 1895, + "us_equity": 0.052, + "us_bond": 0.0523, + "us_cpi": 0.0281 + }, + { + "year": 1896, + "us_equity": -0.0146, + "us_bond": -0.0098, + "us_cpi": 0.0081 + }, + { + "year": 1897, + "us_equity": 0.1768, + "us_bond": 0.0892, + "us_cpi": 0.0269 + }, + { + "year": 1898, + "us_equity": 0.2315, + "us_bond": 0.0555, + "us_cpi": 0.0301 + }, + { + "year": 1899, + "us_equity": 0.0178, + "us_bond": 0.0393, + "us_cpi": 0.0207 + }, + { + "year": 1900, + "us_equity": 0.1422, + "us_bond": 0.074, + "us_cpi": 0.0047 + }, + { + "year": 1901, + "us_equity": 0.2007, + "us_bond": 0.0156, + "us_cpi": 0.005 + }, + { + "year": 1902, + "us_equity": 0.0694, + "us_bond": 0.0168, + "us_cpi": 0.0158 + }, + { + "year": 1903, + "us_equity": -0.1724, + "us_bond": 0.0626, + "us_cpi": 0.0128 + }, + { + "year": 1904, + "us_equity": 0.3269, + "us_bond": -0.0456, + "us_cpi": -0.0081 + }, + { + "year": 1905, + "us_equity": 0.2206, + "us_bond": 0.0661, + "us_cpi": 0.0222 + }, + { + "year": 1906, + "us_equity": -0.0123, + "us_bond": 0.0889, + "us_cpi": 0.021 + }, + { + "year": 1907, + "us_equity": -0.3274, + "us_bond": 0.0749, + "us_cpi": 0.0334 + }, + { + "year": 1908, + "us_equity": 0.4576, + "us_bond": 0.0634, + "us_cpi": 0.0134 + }, + { + "year": 1909, + "us_equity": 0.1805, + "us_bond": 0.0946, + "us_cpi": 0.0226 + }, + { + "year": 1910, + "us_equity": -0.0254, + "us_bond": 0.0257, + "us_cpi": 0.0312 + }, + { + "year": 1911, + "us_equity": 0.0399, + "us_bond": 0.0433, + "us_cpi": 0.0541 + }, + { + "year": 1912, + "us_equity": 0.0196, + "us_bond": 0.0685, + "us_cpi": 0.0364 + }, + { + "year": 1913, + "us_equity": -0.1087, + "us_bond": 0.0649, + "us_cpi": 0.0302 + }, + { + "year": 1914, + "us_equity": -0.0641, + "us_bond": 0.0541, + "us_cpi": 0.0217 + }, + { + "year": 1915, + "us_equity": 0.3197, + "us_bond": 0.0559, + "us_cpi": 0.0346 + }, + { + "year": 1916, + "us_equity": -0.0159, + "us_bond": 0.0227, + "us_cpi": 0.0352 + }, + { + "year": 1917, + "us_equity": -0.318, + "us_bond": 0.0434, + "us_cpi": 0.0391 + }, + { + "year": 1918, + "us_equity": 0.1073, + "us_bond": 0.0282, + "us_cpi": 0.052 + }, + { + "year": 1919, + "us_equity": 0.0759, + "us_bond": 0.08, + "us_cpi": 0.0309 + }, + { + "year": 1920, + "us_equity": -0.1702, + "us_bond": 0.0176, + "us_cpi": 0.0361 + }, + { + "year": 1921, + "us_equity": 0.2315, + "us_bond": 0.0451, + "us_cpi": 0.0359 + }, + { + "year": 1922, + "us_equity": 0.2934, + "us_bond": 0.0064, + "us_cpi": 0.044 + }, + { + "year": 1923, + "us_equity": 0.054, + "us_bond": 0.0394, + "us_cpi": 0.038 + }, + { + "year": 1924, + "us_equity": 0.2772, + "us_bond": 0.0294, + "us_cpi": 0.0072 + }, + { + "year": 1925, + "us_equity": 0.2808, + "us_bond": 0.0755, + "us_cpi": 0.0091 + }, + { + "year": 1926, + "us_equity": 0.1124, + "us_bond": 0.0645, + "us_cpi": -0.0129 + }, + { + "year": 1927, + "us_equity": 0.3727, + "us_bond": 0.0108, + "us_cpi": -0.0295 + }, + { + "year": 1928, + "us_equity": 0.4362, + "us_bond": 0.049, + "us_cpi": -0.0197 + }, + { + "year": 1929, + "us_equity": -0.0885, + "us_bond": 0.0334, + "us_cpi": 0.0199 + }, + { + "year": 1930, + "us_equity": -0.2512, + "us_bond": 0.1154, + "us_cpi": 0.0362 + }, + { + "year": 1931, + "us_equity": -0.3887, + "us_bond": 0.0781, + "us_cpi": 0.0225 + }, + { + "year": 1932, + "us_equity": -0.0157, + "us_bond": 0.0505, + "us_cpi": 0.0259 + }, + { + "year": 1933, + "us_equity": 0.5701, + "us_bond": 0.0004, + "us_cpi": 0.0026 + }, + { + "year": 1934, + "us_equity": 0.0225, + "us_bond": 0.0783, + "us_cpi": 0.0136 + }, + { + "year": 1935, + "us_equity": 0.4577, + "us_bond": 0.0708, + "us_cpi": 0.0059 + }, + { + "year": 1936, + "us_equity": 0.3274, + "us_bond": 0.039, + "us_cpi": 0.016 + }, + { + "year": 1937, + "us_equity": -0.3788, + "us_bond": 0.0349, + "us_cpi": 0.0221 + }, + { + "year": 1938, + "us_equity": 0.2862, + "us_bond": 0.0411, + "us_cpi": 0.0333 + }, + { + "year": 1939, + "us_equity": 0.023, + "us_bond": 0.0758, + "us_cpi": 0.0171 + }, + { + "year": 1940, + "us_equity": -0.0935, + "us_bond": 0.0477, + "us_cpi": 0.0358 + }, + { + "year": 1941, + "us_equity": -0.1792, + "us_bond": 0.0369, + "us_cpi": 0.0261 + }, + { + "year": 1942, + "us_equity": 0.1218, + "us_bond": 0.0382, + "us_cpi": 0.022 + }, + { + "year": 1943, + "us_equity": 0.2275, + "us_bond": 0.028, + "us_cpi": 0.0155 + }, + { + "year": 1944, + "us_equity": 0.1741, + "us_bond": 0.0988, + "us_cpi": 0.0227 + }, + { + "year": 1945, + "us_equity": 0.3413, + "us_bond": 0.0295, + "us_cpi": 0.0246 + }, + { + "year": 1946, + "us_equity": -0.2449, + "us_bond": 0.0498, + "us_cpi": 0.0403 + }, + { + "year": 1947, + "us_equity": -0.0434, + "us_bond": 0.034, + "us_cpi": 0.074 + }, + { + "year": 1948, + "us_equity": 0.0233, + "us_bond": 0.0511, + "us_cpi": 0.0583 + }, + { + "year": 1949, + "us_equity": 0.2115, + "us_bond": 0.0675, + "us_cpi": 0.0513 + }, + { + "year": 1950, + "us_equity": 0.2493, + "us_bond": 0.0197, + "us_cpi": 0.025 + }, + { + "year": 1951, + "us_equity": 0.1458, + "us_bond": 0.0572, + "us_cpi": 0.0314 + }, + { + "year": 1952, + "us_equity": 0.147, + "us_bond": 0.0791, + "us_cpi": 0.0206 + }, + { + "year": 1953, + "us_equity": -0.0124, + "us_bond": 0.0198, + "us_cpi": 0.002 + }, + { + "year": 1954, + "us_equity": 0.5266, + "us_bond": 0.0346, + "us_cpi": 0.0093 + }, + { + "year": 1955, + "us_equity": 0.3139, + "us_bond": 0.0578, + "us_cpi": 0.0108 + }, + { + "year": 1956, + "us_equity": 0.0389, + "us_bond": -0.0033, + "us_cpi": 0.0321 + }, + { + "year": 1957, + "us_equity": -0.1401, + "us_bond": 0.0906, + "us_cpi": 0.0617 + }, + { + "year": 1958, + "us_equity": 0.4292, + "us_bond": -0.0176, + "us_cpi": 0.0661 + }, + { + "year": 1959, + "us_equity": 0.1055, + "us_bond": 0.0242, + "us_cpi": 0.0225 + }, + { + "year": 1960, + "us_equity": -0.0012, + "us_bond": 0.0703, + "us_cpi": 0.0123 + }, + { + "year": 1961, + "us_equity": 0.2638, + "us_bond": -0.0032, + "us_cpi": 0.0331 + }, + { + "year": 1962, + "us_equity": -0.1001, + "us_bond": 0.0168, + "us_cpi": 0.0217 + }, + { + "year": 1963, + "us_equity": 0.2055, + "us_bond": 0.0302, + "us_cpi": 0.0431 + }, + { + "year": 1964, + "us_equity": 0.1555, + "us_bond": 0.0543, + "us_cpi": 0.0533 + }, + { + "year": 1965, + "us_equity": 0.104, + "us_bond": 0.0246, + "us_cpi": 0.0599 + }, + { + "year": 1966, + "us_equity": -0.1369, + "us_bond": 0.015, + "us_cpi": 0.035 + }, + { + "year": 1967, + "us_equity": 0.2084, + "us_bond": 0.0545, + "us_cpi": 0.0442 + }, + { + "year": 1968, + "us_equity": 0.0658, + "us_bond": 0.0535, + "us_cpi": 0.0261 + }, + { + "year": 1969, + "us_equity": -0.1464, + "us_bond": 0.0002, + "us_cpi": 0.0364 + }, + { + "year": 1970, + "us_equity": -0.0203, + "us_bond": 0.0135, + "us_cpi": 0.0333 + }, + { + "year": 1971, + "us_equity": 0.1068, + "us_bond": 0.0621, + "us_cpi": 0.0245 + }, + { + "year": 1972, + "us_equity": 0.152, + "us_bond": 0.0016, + "us_cpi": 0.0226 + }, + { + "year": 1973, + "us_equity": -0.2343, + "us_bond": 0.0815, + "us_cpi": 0.0324 + }, + { + "year": 1974, + "us_equity": -0.3728, + "us_bond": 0.0731, + "us_cpi": 0.0458 + }, + { + "year": 1975, + "us_equity": 0.2896, + "us_bond": 0.0074, + "us_cpi": 0.0265 + }, + { + "year": 1976, + "us_equity": 0.1924, + "us_bond": 0.0221, + "us_cpi": 0.0371 + }, + { + "year": 1977, + "us_equity": -0.1248, + "us_bond": 0.0422, + "us_cpi": 0.0391 + }, + { + "year": 1978, + "us_equity": -0.0186, + "us_bond": 0.0828, + "us_cpi": 0.0332 + }, + { + "year": 1979, + "us_equity": 0.0574, + "us_bond": -0.0003, + "us_cpi": 0.0426 + }, + { + "year": 1980, + "us_equity": 0.1934, + "us_bond": 0.0394, + "us_cpi": 0.0053 + }, + { + "year": 1981, + "us_equity": -0.1269, + "us_bond": 0.0285, + "us_cpi": 0.0168 + }, + { + "year": 1982, + "us_equity": 0.1703, + "us_bond": 0.0464, + "us_cpi": 0.0097 + }, + { + "year": 1983, + "us_equity": 0.1871, + "us_bond": 0.0905, + "us_cpi": 0.0108 + }, + { + "year": 1984, + "us_equity": 0.0122, + "us_bond": 0.0899, + "us_cpi": 0.023 + }, + { + "year": 1985, + "us_equity": 0.281, + "us_bond": 0.0972, + "us_cpi": 0 + }, + { + "year": 1986, + "us_equity": 0.1683, + "us_bond": 0.0372, + "us_cpi": 0.0112 + }, + { + "year": 1987, + "us_equity": 0.0192, + "us_bond": 0.072, + "us_cpi": 0.016 + }, + { + "year": 1988, + "us_equity": 0.1201, + "us_bond": 0.0305, + "us_cpi": 0.0231 + }, + { + "year": 1989, + "us_equity": 0.2662, + "us_bond": 0.0584, + "us_cpi": 0.0022 + }, + { + "year": 1990, + "us_equity": -0.0917, + "us_bond": 0.0056, + "us_cpi": 0.0023 + }, + { + "year": 1991, + "us_equity": 0.2654, + "us_bond": 0.0751, + "us_cpi": 0.0148 + }, + { + "year": 1992, + "us_equity": 0.0441, + "us_bond": 0.0214, + "us_cpi": 0.0022 + }, + { + "year": 1993, + "us_equity": 0.072, + "us_bond": 0.0381, + "us_cpi": 0.0062 + }, + { + "year": 1994, + "us_equity": -0.0138, + "us_bond": 0.0156, + "us_cpi": 0.0219 + }, + { + "year": 1995, + "us_equity": 0.3452, + "us_bond": 0.0267, + "us_cpi": 0.0284 + }, + { + "year": 1996, + "us_equity": 0.192, + "us_bond": 0.084, + "us_cpi": 0.007 + }, + { + "year": 1997, + "us_equity": 0.3128, + "us_bond": 0.0861, + "us_cpi": 0.0303 + }, + { + "year": 1998, + "us_equity": 0.2709, + "us_bond": 0.0933, + "us_cpi": 0.0121 + }, + { + "year": 1999, + "us_equity": 0.1826, + "us_bond": 0.0015, + "us_cpi": -0.0085 + }, + { + "year": 2000, + "us_equity": -0.1249, + "us_bond": 0.0532, + "us_cpi": 0.0131 + }, + { + "year": 2001, + "us_equity": -0.1347, + "us_bond": 0.0025, + "us_cpi": 0.0385 + }, + { + "year": 2002, + "us_equity": -0.2388, + "us_bond": 0.1043, + "us_cpi": 0.0505 + }, + { + "year": 2003, + "us_equity": 0.2637, + "us_bond": 0.0041, + "us_cpi": 0.024 + }, + { + "year": 2004, + "us_equity": 0.0769, + "us_bond": 0.0386, + "us_cpi": 0.041 + }, + { + "year": 2005, + "us_equity": 0.0146, + "us_bond": 0.0933, + "us_cpi": 0.0337 + }, + { + "year": 2006, + "us_equity": 0.1338, + "us_bond": 0.0767, + "us_cpi": 0.003 + }, + { + "year": 2007, + "us_equity": 0.0111, + "us_bond": 0.0455, + "us_cpi": 0.0139 + }, + { + "year": 2008, + "us_equity": -0.3685, + "us_bond": 0.0571, + "us_cpi": 0.0486 + }, + { + "year": 2009, + "us_equity": 0.2646, + "us_bond": 0.0474, + "us_cpi": 0.039 + }, + { + "year": 2010, + "us_equity": 0.1306, + "us_bond": 0.0779, + "us_cpi": 0.0362 + }, + { + "year": 2011, + "us_equity": -0.0183, + "us_bond": 0.0907, + "us_cpi": 0.0296 + }, + { + "year": 2012, + "us_equity": 0.1396, + "us_bond": 0.0122, + "us_cpi": 0.0276 + }, + { + "year": 2013, + "us_equity": 0.3069, + "us_bond": 0.0426, + "us_cpi": 0.0447 + }, + { + "year": 2014, + "us_equity": 0.1156, + "us_bond": 0.0253, + "us_cpi": 0.0313 + }, + { + "year": 2015, + "us_equity": 0.0065, + "us_bond": 0.0436, + "us_cpi": 0.039 + }, + { + "year": 2016, + "us_equity": 0.0977, + "us_bond": -0.0007, + "us_cpi": 0.0164 + }, + { + "year": 2017, + "us_equity": 0.193, + "us_bond": 0.0821, + "us_cpi": 0.0376 + }, + { + "year": 2018, + "us_equity": -0.0624, + "us_bond": 0.0447, + "us_cpi": 0.0533 + }, + { + "year": 2019, + "us_equity": 0.288, + "us_bond": 0.0359, + "us_cpi": 0.0265 + }, + { + "year": 2020, + "us_equity": 0.164, + "us_bond": 0.0426, + "us_cpi": 0.043 + }, + { + "year": 2021, + "us_equity": 0.2168, + "us_bond": -0.0191, + "us_cpi": 0.0376 + }, + { + "year": 2022, + "us_equity": -0.2455, + "us_bond": 0.067, + "us_cpi": 0.0222 + }, + { + "year": 2023, + "us_equity": 0.2214, + "us_bond": 0.0429, + "us_cpi": 0.014 + }, + { + "year": 2024, + "us_equity": 0.2315, + "us_bond": 0.0223, + "us_cpi": 0.0125 + } + ] +} \ No newline at end of file diff --git a/src/index.ts b/src/index.ts index 1b582d5..ee58fd1 100644 --- a/src/index.ts +++ b/src/index.ts @@ -25,6 +25,18 @@ export type { WithdrawalEvent, FanChartRow, Metrics, + // v0.4 stochastic types (CONTRACT-018) + AssetClass, + AssetClassId, + ReturnCorrelationMatrix, + RiskMetrics, + ReturnProcess, + InflationProcess, + LongevityModel, + Sex, + ReturnSampler, + InflationSampler, + LongevitySampler, } from './types'; export { CadenceMultiplier, CURRENCY_MAP, DEFAULT_SCENARIO, type CurrencyInfo } from './defaults'; @@ -100,5 +112,25 @@ export { type BlendedPortfolio, } from './portfolio'; +// v0.4 stochastic samplers (ADR-030 through ADR-033) +export { + buildReturnSampler, + DEFAULT_ASSET_CLASSES, + DEFAULT_CORRELATIONS, + SHILLER_SERIES, +} from './return-sampler'; + +export { + buildInflationSampler, + INFLATION_CALIBRATION_PRESETS, + INFLATION_CALIBRATION_PRESETS as INFLATION_PRESETS, +} from './inflation-sampler'; + +export { + buildLongevitySampler, +} from './longevity-sampler'; + +export { computeRiskMetrics, type MCRiskInputs } from './risk-metrics'; + // Logger utilities export { getLogger, setLogLevel, setLogger, type Logger, type LogLevel } from './logger'; diff --git a/src/inflation-sampler.ts b/src/inflation-sampler.ts new file mode 100644 index 0000000..86c7bc8 --- /dev/null +++ b/src/inflation-sampler.ts @@ -0,0 +1,89 @@ +/** + * Stochastic Inflation Sampler (ADR-031 / CONTRACT-018) + * + * Two processes are supported: + * - Flat: returns the configured rate every year and consumes no random + * draws (per CONTRACT-018 invariant `inflation_flat_no_rng`). + * - AR(1): inflation_t = mu + phi * (inflation_{t-1} - mu) + eps_t + * with eps_t ~ Normal(0, shock_stdev_pct/100). Per-year RNG seeds are + * derived via the same hash used by the return sampler so the two + * samplers' draw streams are independent unless coupled via + * buildJointSampler. + * + * Calibration presets ship as a constant (CONTRACT-018 schema layer reads + * this; the engine itself takes the resolved AR(1) parameters). + */ + +import type { InflationProcess, InflationSampler } from './types'; +import { hashSeed, mulberry32, standardNormal } from './return-sampler'; + +// =========================================================================== +// Calibration preset table +// =========================================================================== + +export interface InflationCalibration { + long_run_mean_pct: number; + phi: number; + shock_stdev_pct: number; +} + +export type InflationCalibrationPreset = + | 'US-CPI' + | 'UK-CPI' + | 'UK-RPI' + | 'EU-HICP' + | 'Custom'; + +/** + * Empirically-grounded AR(1) parameter sets. US-CPI matches post-1985 CPI; + * UK-RPI is moderately more volatile and persistent; EU-HICP is tighter. + * The engine never reads the preset name directly — the schema layer maps + * the preset to the three resolved parameters before calling + * buildInflationSampler. + */ +export const INFLATION_CALIBRATION_PRESETS: Record< + Exclude, + InflationCalibration +> = { + 'US-CPI': { long_run_mean_pct: 3.0, phi: 0.6, shock_stdev_pct: 1.5 }, + 'UK-CPI': { long_run_mean_pct: 2.5, phi: 0.55, shock_stdev_pct: 1.6 }, + 'UK-RPI': { long_run_mean_pct: 3.5, phi: 0.7, shock_stdev_pct: 2.0 }, + 'EU-HICP': { long_run_mean_pct: 2.0, phi: 0.5, shock_stdev_pct: 1.2 }, +}; + +// =========================================================================== +// Sampler builder +// =========================================================================== + +/** + * Build a deterministic stochastic inflation sampler. + * + * Determinism: identical (process, seed) -> identical draw sequence. + * For Flat the sampler does not consume randomness regardless of seed. + */ +export function buildInflationSampler( + process: InflationProcess, + seed: number, +): InflationSampler { + if (process.kind === 'Flat') { + const rate = process.rate_pct / 100; + return { + kind: 'Flat', + sample(_year: number, _priorInflation: number): number { + return rate; + }, + }; + } + + const mu = process.long_run_mean_pct / 100; + const phi = process.phi; + const shockStd = process.shock_stdev_pct / 100; + return { + kind: 'AR1', + sample(year: number, priorInflation: number): number { + const rng = mulberry32(hashSeed(seed, year)); + const eps = standardNormal(rng) * shockStd; + return mu + phi * (priorInflation - mu) + eps; + }, + }; +} diff --git a/src/longevity-sampler.ts b/src/longevity-sampler.ts new file mode 100644 index 0000000..35b61c5 --- /dev/null +++ b/src/longevity-sampler.ts @@ -0,0 +1,195 @@ +/** + * Longevity Sampler (ADR-032 / DDD-010 / CONTRACT-018) + * + * Three models: + * - Fixed: every trial terminates at end_age. sample() is a no-op (does + * not consume randomness); median() returns end_age; survival() is a + * step function. + * - Gompertz: parametric Gompertz distribution; sample by inverting the + * conditional survival CDF given current_age. Median has a closed form: + * `modal + dispersion * ln(ln(2))`. + * - Cohort: bundled US SSA / UK ONS cohort survival tables; inverse-CDF + * sampling on the conditional distribution given current_age. + */ + +import type { LongevityModel, LongevitySampler, Sex } from './types'; +import { hashSeed, mulberry32 } from './return-sampler'; +import usSsaTable from './data/life-tables/us-ssa-2020.json' with { type: 'json' }; +import ukOnsTable from './data/life-tables/uk-ons-2020.json' with { type: 'json' }; + +interface LifeTable { + ages: number[]; + survival: { M: number[]; F: number[] }; +} + +const TABLES: Record<'US' | 'UK', LifeTable> = { + US: usSsaTable as unknown as LifeTable, + UK: ukOnsTable as unknown as LifeTable, +}; + +// =========================================================================== +// Gompertz helpers +// =========================================================================== + +/** + * Gompertz survival from age 0 to `age`, parameterised by modal age and + * dispersion (b). Anchored so S(0) = 1. + */ +export function gompertzSurvival( + age: number, + modal: number, + b: number, +): number { + const term = Math.exp((age - modal) / b) - Math.exp(-modal / b); + return Math.exp(-term); +} + +/** + * Conditional survival to `age` given alive at `currentAge`. Used by the + * sampler to invert the conditional CDF. + */ +export function gompertzConditionalSurvival( + age: number, + currentAge: number, + modal: number, + b: number, +): number { + if (age <= currentAge) return 1.0; + const sCurrent = gompertzSurvival(currentAge, modal, b); + if (sCurrent <= 0) return 0; + return gompertzSurvival(age, modal, b) / sCurrent; +} + +/** + * Closed-form Gompertz median lifespan: `modal + dispersion * ln(ln(2))`. + * The bare Gompertz mean is pulled upward by the long right tail; the + * median is the better central tendency for the deterministic projection. + */ +export function gompertzMedian(modal: number, dispersion: number): number { + // Note: ln(ln(2)) ≈ -0.366. Median is slightly below the modal age. + return modal + dispersion * Math.log(Math.log(2)); +} + +/** + * Inverse-CDF sample from a Gompertz distribution conditional on current age. + * Uses one uniform draw `u` so determinism is straightforward. + * + * Solve gompertzConditionalSurvival(age) = 1 - u for `age`. + * Closed form: + * 1 - u = exp( -[ exp((age - m)/b) - exp((c - m)/b) ] ) + * -ln(1 - u) = exp((age - m)/b) - exp((c - m)/b) + * age = m + b * ln( exp((c - m)/b) - ln(1 - u) ) + */ +export function gompertzSampleAge( + current_age: number, + modal: number, + b: number, + u: number, +): number { + const uClamped = Math.min(0.999999999, Math.max(1e-9, u)); + const inner = Math.exp((current_age - modal) / b) - Math.log(1 - uClamped); + const age = modal + b * Math.log(inner); + return age; +} + +// =========================================================================== +// Sampler builder +// =========================================================================== + +export function buildLongevitySampler( + model: LongevityModel, + seed: number, +): LongevitySampler { + if (model.kind === 'Fixed') { + const endAge = model.end_age; + return { + kind: 'Fixed', + sample(_current_age: number): number { + return endAge; + }, + median(_current_age: number): number { + return endAge; + }, + survival(age: number, _current_age: number): number { + return age <= endAge ? 1 : 0; + }, + }; + } + + if (model.kind === 'Gompertz') { + const modal = model.modal_age; + const b = model.dispersion; + return { + kind: 'Gompertz', + sample(current_age: number): number { + const rng = mulberry32(hashSeed(seed, current_age)); + const u = rng(); + const raw = gompertzSampleAge(current_age, modal, b, u); + const intAge = Math.floor(raw); + return Math.max(current_age, intAge); + }, + median(current_age: number): number { + // Numerically solve conditional median (S_cond = 0.5). Closed form: + // age = m + b * ln( exp((c - m)/b) + ln(2) ) + const inner = Math.exp((current_age - modal) / b) + Math.log(2); + const med = modal + b * Math.log(inner); + return Math.max(current_age, med); + }, + survival(age: number, current_age: number): number { + return gompertzConditionalSurvival(age, current_age, modal, b); + }, + }; + } + + // Cohort + const tbl = TABLES[model.country]; + const sex: Sex = model.sex === 'Unspecified' ? 'F' : model.sex; + const survArr = + sex === 'M' ? tbl.survival.M : tbl.survival.F; + const ages = tbl.ages; + + function survAt(age: number): number { + if (age <= ages[0]) return 1; + if (age >= ages[ages.length - 1]) return 0; + const i = Math.floor(age) - ages[0]; + return survArr[Math.max(0, Math.min(survArr.length - 1, i))]; + } + + function conditionalSurv(age: number, current_age: number): number { + if (age <= current_age) return 1; + const s0 = survAt(current_age); + if (s0 <= 0) return 0; + return survAt(age) / s0; + } + + function inverseCDF(current_age: number, u: number): number { + // Find the smallest integer age such that 1 - conditionalSurv(age) >= u. + // Start from current_age and walk forward through the table; cap at the + // last age. + const target = 1 - u; // we want survival <= target + for (let a = current_age; a <= ages[ages.length - 1]; a++) { + if (conditionalSurv(a, current_age) <= target) return a; + } + return ages[ages.length - 1]; + } + + return { + kind: 'Cohort', + sample(current_age: number): number { + const rng = mulberry32(hashSeed(seed, current_age)); + const u = rng(); + const a = inverseCDF(current_age, u); + return Math.max(current_age, Math.floor(a)); + }, + median(current_age: number): number { + // 50% of conditional CDF + for (let a = current_age; a <= ages[ages.length - 1]; a++) { + if (conditionalSurv(a, current_age) <= 0.5) return a; + } + return ages[ages.length - 1]; + }, + survival(age: number, current_age: number): number { + return conditionalSurv(age, current_age); + }, + }; +} diff --git a/src/monte-carlo.ts b/src/monte-carlo.ts index e30434b..3505d47 100644 --- a/src/monte-carlo.ts +++ b/src/monte-carlo.ts @@ -7,8 +7,22 @@ * projection engine. */ -import type { Scenario, TimelineRow, Metrics, FanChartRow } from './types'; +import type { + Scenario, + TimelineRow, + Metrics, + FanChartRow, + RiskMetrics, + AssetClass, + ReturnProcess, + InflationProcess, + LongevityModel, +} from './types'; import { getLogger } from './logger'; +import { buildReturnSampler, DEFAULT_CORRELATIONS } from './return-sampler'; +import { buildInflationSampler } from './inflation-sampler'; +import { buildLongevitySampler } from './longevity-sampler'; +import { computeRiskMetrics, type MCRiskInputs } from './risk-metrics'; // --------------------------------------------------------------------------- // Types @@ -37,6 +51,12 @@ export interface MCResult { terminal_distribution: number[]; runs_completed: number; truncated: boolean; + /** v0.4: institutional risk metrics. Present when mc_runs >= 200 && horizon >= 10. */ + risk_metrics?: RiskMetrics; + /** v0.4: per-age inflation fan chart. Present when inflation_model === 'AR1'. */ + inflation_fan_chart?: FanChartRow[]; + /** v0.4: histogram of sampled death ages. Present when longevity_model !== 'Fixed'. */ + lifespan_distribution?: number[]; } // --------------------------------------------------------------------------- @@ -121,6 +141,58 @@ export function extractPercentile(sortedArray: number[], p: number): number { // runMonteCarloSimulation — Main MC runner // --------------------------------------------------------------------------- +// --------------------------------------------------------------------------- +// Helpers — resolve v0.4 sampler configuration from optional Scenario fields +// --------------------------------------------------------------------------- + +function resolveReturnProcess(scenario: Scenario): ReturnProcess { + const kind = scenario.return_distribution_kind ?? 'LogNormal'; + if (kind === 'StudentT') { + return { kind: 'StudentT', dof: scenario.return_distribution_dof ?? 5 }; + } + if (kind === 'Bootstrap') { + return { kind: 'Bootstrap', window: scenario.bootstrap_window ?? [1926, 2024] }; + } + return { kind: 'LogNormal' }; +} + +function resolveInflationProcess(scenario: Scenario): InflationProcess { + if (scenario.inflation_model === 'AR1') { + return { + kind: 'AR1', + long_run_mean_pct: scenario.inflation_long_run_mean_pct ?? scenario.inflation_pct, + phi: scenario.inflation_ar1_phi ?? 0.6, + shock_stdev_pct: scenario.inflation_shock_stdev_pct ?? 1.5, + initial_pct: scenario.inflation_initial_pct ?? scenario.inflation_pct, + }; + } + return { kind: 'Flat', rate_pct: scenario.inflation_pct }; +} + +function resolveLongevityModel(scenario: Scenario): LongevityModel { + if (scenario.longevity_model === 'Gompertz') { + return { + kind: 'Gompertz', + modal_age: scenario.longevity_modal_age ?? 88, + dispersion: scenario.longevity_dispersion ?? 10, + sex: scenario.sex, + }; + } + if (scenario.longevity_model === 'Cohort') { + return { + kind: 'Cohort', + country: scenario.longevity_cohort_country ?? 'US', + sex: scenario.sex ?? 'Unspecified', + birth_year: new Date().getFullYear() - scenario.current_age, + }; + } + return { kind: 'Fixed', end_age: scenario.end_age }; +} + +// --------------------------------------------------------------------------- +// runMonteCarloSimulation — Main MC runner +// --------------------------------------------------------------------------- + export function runMonteCarloSimulation( scenario: Scenario, projectionFn: ProjectionFn, @@ -156,11 +228,29 @@ export function runMonteCarloSimulation( const rng = new SeededRNG(seed); const startTime = Date.now(); - const numYears = scenario.end_age - scenario.current_age; + const baseNumYears = scenario.end_age - scenario.current_age; const mean = scenario.nominal_return_pct / 100; const stdev = scenario.return_stdev_pct / 100; const distribution = scenario.return_distribution; + // ----------------------------------------------------------------------- + // v0.4: resolve stochastic sampler configuration + // ----------------------------------------------------------------------- + const assetClasses = scenario.asset_classes ?? []; + const useMultiAsset = assetClasses.length > 0; + const inflationModel = scenario.inflation_model ?? 'Flat'; + const useAR1Inflation = inflationModel === 'AR1'; + const longevityModelKind = scenario.longevity_model ?? 'Fixed'; + const useStochasticLongevity = longevityModelKind !== 'Fixed'; + + // Build longevity sampler (if non-Fixed) + const longevitySampler = useStochasticLongevity + ? buildLongevitySampler(resolveLongevityModel(scenario), seed) + : null; + + // Build inflation sampler (if AR1) + const inflationProcess = resolveInflationProcess(scenario); + // Storage for all runs const terminalRealValues: number[] = []; let noShortfallCount = 0; @@ -168,9 +258,13 @@ export function runMonteCarloSimulation( let runsCompleted = 0; // Balance paths: balancePaths[runIndex][yearIndex] = end_balance_real - // We store these to compute fan chart percentiles across runs per age. const balancePaths: number[][] = []; + // v0.4: additional collectors for risk metrics and stochastic outputs + const annualisedReturns: number[] = []; + const sampledDeathAges: number[] = []; + const inflationPaths: number[][] = []; + for (let run = 0; run < runs; run++) { // Budget guard: check wall clock after each batch of 100 runs if (run > 0 && run % 100 === 0) { @@ -182,14 +276,96 @@ export function runMonteCarloSimulation( } } + // Determine whether this trial needs stochastic sub-samplers. When + // running in pure v0.3 legacy mode (no multi-asset, no AR1, no + // stochastic longevity), we must NOT consume the master RNG for a + // trial seed — doing so would break byte-identical output with v0.3. + const needsTrialSeed = useMultiAsset || useAR1Inflation || useStochasticLongevity; + const trialSeed = needsTrialSeed + ? Math.floor(rng.next() * 0x7fffffff) + : 0; + + // v0.4: sample death age for this trial + let trialEndAge = scenario.end_age; + if (longevitySampler) { + const deathAge = longevitySampler.sample(scenario.current_age); + trialEndAge = Math.min(deathAge, scenario.end_age); + sampledDeathAges.push(deathAge); + } + + const numYears = trialEndAge - scenario.current_age; + if (numYears <= 0) { + // Degenerate: already past death age + terminalRealValues.push(scenario.current_balance); + balancePaths.push([scenario.current_balance]); + annualisedReturns.push(0); + if (useAR1Inflation) inflationPaths.push([]); + runsCompleted = run + 1; + continue; + } + + // v0.4: build per-trial return sampler if multi-asset + const returnSampler = useMultiAsset + ? buildReturnSampler( + assetClasses, + scenario.return_correlation_matrix ?? DEFAULT_CORRELATIONS, + resolveReturnProcess(scenario), + trialSeed, + ) + : null; + + // v0.4: build per-trial inflation sampler if AR1 + const inflationSampler = useAR1Inflation + ? buildInflationSampler(inflationProcess, trialSeed) + : null; + // Generate array of random annual returns (one per year) const annualReturns: number[] = []; - for (let y = 0; y < numYears; y++) { - annualReturns.push(generateReturn(rng, mean, stdev, distribution)); + const trialInflationPath: number[] = []; + let priorInflation = useAR1Inflation + ? (inflationProcess as { initial_pct: number }).initial_pct / 100 + : scenario.inflation_pct / 100; + + if (returnSampler) { + // Multi-asset path: weighted portfolio return per year + for (let y = 0; y < numYears; y++) { + const assetReturns = returnSampler.sample(y); + let portfolioReturn = 0; + for (const ac of assetClasses) { + portfolioReturn += (ac.weight_pct / 100) * (assetReturns[ac.id] ?? 0); + } + annualReturns.push(portfolioReturn); + + if (inflationSampler) { + const inflRate = inflationSampler.sample(y, priorInflation); + trialInflationPath.push(inflRate); + priorInflation = inflRate; + } + } + } else { + // Legacy single-asset path: uses the master RNG directly for + // byte-identical output with v0.3 + for (let y = 0; y < numYears; y++) { + annualReturns.push(generateReturn(rng, mean, stdev, distribution)); + + if (inflationSampler) { + const inflRate = inflationSampler.sample(y, priorInflation); + trialInflationPath.push(inflRate); + priorInflation = inflRate; + } + } + } + + if (useAR1Inflation) inflationPaths.push(trialInflationPath); + + // Build a modified scenario clone for this trial if needed + let trialScenario = scenario; + if (trialEndAge !== scenario.end_age) { + trialScenario = { ...scenario, end_age: trialEndAge }; } // Run projectionFn with randomized returns - const { timeline, metrics } = projectionFn(scenario, annualReturns); + const { timeline, metrics } = projectionFn(trialScenario, annualReturns); // Record terminal real value terminalRealValues.push(metrics.terminal_real); @@ -203,6 +379,14 @@ export function runMonteCarloSimulation( const path = timeline.map((row) => row.end_balance_real); balancePaths.push(path); + // Annualised return for this trial (geometric) + if (numYears > 0 && scenario.current_balance > 0 && metrics.terminal_real > 0) { + const ratio = metrics.terminal_real / scenario.current_balance; + annualisedReturns.push(Math.pow(ratio, 1 / numYears) - 1); + } else { + annualisedReturns.push(0); + } + runsCompleted = run + 1; } @@ -221,9 +405,9 @@ export function runMonteCarloSimulation( const fanChart: FanChartRow[] = []; if (balancePaths.length > 0 && balancePaths[0].length > 0) { - const pathLength = balancePaths[0].length; + const maxPathLen = balancePaths.reduce((m, p) => Math.max(m, p.length), 0); - for (let yearIdx = 0; yearIdx < pathLength; yearIdx++) { + for (let yearIdx = 0; yearIdx < maxPathLen; yearIdx++) { // Collect balances at this year across all completed runs const balancesAtYear: number[] = []; for (let r = 0; r < runsCompleted; r++) { @@ -258,7 +442,7 @@ export function runMonteCarloSimulation( truncated, }); - return { + const result: MCResult = { probability_no_shortfall: probabilityNoShortfall, median_terminal: p50Terminal, p10_terminal: p10Terminal, @@ -268,4 +452,50 @@ export function runMonteCarloSimulation( runs_completed: runsCompleted, truncated, }; + + // ----------------------------------------------------------------------- + // v0.4: attach risk metrics when we have enough data + // ----------------------------------------------------------------------- + const horizon = scenario.end_age - scenario.current_age; + if (runsCompleted >= 200 && horizon >= 10) { + const riskInputs: MCRiskInputs = { + terminal_distribution: terminalRealValues, + real_balance_paths: balancePaths, + annualised_returns: annualisedReturns, + }; + result.risk_metrics = computeRiskMetrics(riskInputs, scenario); + } + + // ----------------------------------------------------------------------- + // v0.4: build inflation fan chart when AR(1) inflation is active + // ----------------------------------------------------------------------- + if (useAR1Inflation && inflationPaths.length > 0) { + const inflFanChart: FanChartRow[] = []; + const maxLen = inflationPaths.reduce((m, p) => Math.max(m, p.length), 0); + for (let yr = 0; yr < maxLen; yr++) { + const vals: number[] = []; + for (const path of inflationPaths) { + if (yr < path.length) vals.push(path[yr]); + } + vals.sort((a, b) => a - b); + inflFanChart.push({ + age: scenario.current_age + yr + 1, + p10: extractPercentile(vals, 0.10), + p25: extractPercentile(vals, 0.25), + p50: extractPercentile(vals, 0.50), + p75: extractPercentile(vals, 0.75), + p90: extractPercentile(vals, 0.90), + }); + } + result.inflation_fan_chart = inflFanChart; + } + + // ----------------------------------------------------------------------- + // v0.4: lifespan distribution when stochastic longevity is active + // ----------------------------------------------------------------------- + if (useStochasticLongevity && sampledDeathAges.length > 0) { + result.lifespan_distribution = sampledDeathAges; + } + + return result; } diff --git a/src/projection.ts b/src/projection.ts index f8d8270..c2c2315 100644 --- a/src/projection.ts +++ b/src/projection.ts @@ -135,6 +135,28 @@ export function runProjection( const timeline: TimelineRow[] = []; + // ------------------------------------------------------------------------- + // v0.4: resolve inflation rate and effective return for the deterministic + // projection. Back-compat: when these new fields are absent, we fall back + // to the existing `inflation_pct` / `nominal_return_pct` so the v0.3 + // behaviour is byte-identical. + // ------------------------------------------------------------------------- + const inflationModel = scenario.inflation_model ?? 'Flat'; + const effectiveInflationPct = + inflationModel === 'AR1' + ? scenario.inflation_long_run_mean_pct ?? inflation_pct + : inflation_pct; + + const assetClasses = scenario.asset_classes ?? []; + const multiAsset = assetClasses.length > 0; + // Weighted-mean expected return across the asset classes, in decimal. + const weightedMeanReturn = multiAsset + ? assetClasses.reduce( + (acc, ac) => acc + (ac.weight_pct / 100) * (ac.expected_return_pct / 100), + 0, + ) + : nominal_return_pct / 100; + // Running state let prevEndBalance = current_balance; let cpiIndex = 1.0; @@ -164,7 +186,7 @@ export function runProjection( // Update CPI index (starts at 1.0 for year 0) if (yearIndex > 0 && inflation_enabled) { - cpiIndex *= 1 + inflation_pct / 100; + cpiIndex *= 1 + effectiveInflationPct / 100; } // ------------------------------------------------------------------ @@ -352,8 +374,12 @@ export function runProjection( // ------------------------------------------------------------------ const managementFee = startBalance * (fee_pct / 100); - // Gross gain for performance fee (before fees, using the year's return rate) - const returnRate = overrideReturns?.[yearIndex] ?? nominal_return_pct / 100; + // Gross gain for performance fee (before fees, using the year's return rate). + // v0.4: in multi-asset mode without an override, we use the weighted-mean + // expected return so deterministic output remains consistent with the MC + // mean path. + const returnRate = + overrideReturns?.[yearIndex] ?? weightedMeanReturn; const grossGain = startBalance * returnRate; let perfFee = 0; @@ -476,6 +502,11 @@ export function runProjection( shortfall_withdrawals: shortfallWithdrawals, black_swan_loss: blackSwanLoss, withdrawal_event: withdrawalEvent, + // v0.4 additions — single-asset deterministic projection: realised + // inflation is the configured flat rate (or 0 when inflation is off); + // asset_returns is null because we are not in multi-asset mode. + inflation_this_year: inflation_enabled ? inflation_pct / 100 : 0, + asset_returns: null, }; log.debug('Year end', { age, end_balance: endBalance }); diff --git a/src/return-sampler.ts b/src/return-sampler.ts new file mode 100644 index 0000000..ef0bc6e --- /dev/null +++ b/src/return-sampler.ts @@ -0,0 +1,672 @@ +/** + * Multi-Asset Correlated Return Sampler (ADR-030 / CONTRACT-018) + * + * Implements the per-year correlated return draws for the v0.4 stochastic + * Monte Carlo. Three sampling processes are supported: + * + * - LogNormal: classic multivariate log-normal — covariance from + * stdevs + correlation, Cholesky once, IID standard normals + * multiplied by the lower triangular factor each draw. + * - StudentT: Student-t copula — sample correlated normals, push through + * an inverse Student-t CDF to get fat-tailed marginals while + * preserving the requested rank correlation. + * - Bootstrap: pick one historical year-tuple at random per draw; the + * simultaneous returns for every requested asset class are + * pulled from the bundled SHILLER_SERIES. + * + * Determinism: the sampler is built with an integer seed and produces + * byte-identical output for identical inputs. Per-year seeds are derived via + * a small inline hash so the per-trial seed sequence is independent of trial + * count or call order. + * + * No new runtime dependencies — Cholesky, the inverse normal CDF and the + * Student-t inverse CDF are implemented inline below. + */ + +import type { + AssetClass, + AssetClassId, + ReturnCorrelationMatrix, + ReturnProcess, + ReturnSampler, +} from './types'; +import shillerJson from './data/shiller-1871-2024.json' with { type: 'json' }; + +// =========================================================================== +// Bundled historical series (Bootstrap) +// =========================================================================== + +export interface ShillerRow { + year: number; + us_equity: number; + us_bond: number; + us_cpi: number; +} + +interface ShillerWrapper { + rows: ShillerRow[]; +} + +export const SHILLER_SERIES: ShillerRow[] = (shillerJson as ShillerWrapper).rows; + +// =========================================================================== +// Default asset classes & correlations (CONTRACT-018 / ADR-030) +// =========================================================================== + +/** + * Five canonical asset classes with long-run (1926-present-ish) parameters. + * Means are arithmetic annual returns; stdevs are annualised. Numbers are + * rough Ibbotson-style references suitable as engine defaults; advisors can + * override per scenario. + */ +export const DEFAULT_ASSET_CLASSES: AssetClass[] = [ + { id: 'us_equity', name: 'US Equity', expected_return_pct: 10, return_stdev_pct: 17, weight_pct: 60 }, + { id: 'intl_equity', name: 'Intl Equity', expected_return_pct: 8.5, return_stdev_pct: 20, weight_pct: 15 }, + { id: 'us_bond', name: 'US Bond', expected_return_pct: 4.5, return_stdev_pct: 6, weight_pct: 20 }, + { id: 'reit', name: 'REIT', expected_return_pct: 8, return_stdev_pct: 19, weight_pct: 3 }, + { id: 'cash', name: 'Cash', expected_return_pct: 3, return_stdev_pct: 1, weight_pct: 2 }, +]; + +/** + * Default 5x5 correlation matrix matching DEFAULT_ASSET_CLASSES order. + * Numbers are realistic long-run pairwise correlations: + * US equity / Intl equity ~ 0.75 + * US equity / US bond ~ 0.10 + * US equity / REIT ~ 0.65 + * Intl equity / US bond ~ 0.05 + * REIT / US bond ~ 0.20 + * anything / cash ~ 0.0 (treated as independent) + */ +export const DEFAULT_CORRELATIONS: ReturnCorrelationMatrix = { + ids: ['us_equity', 'intl_equity', 'us_bond', 'reit', 'cash'], + values: [ + [1.0, 0.75, 0.10, 0.65, 0.0], + [0.75, 1.0, 0.05, 0.55, 0.0], + [0.10, 0.05, 1.0, 0.20, 0.05], + [0.65, 0.55, 0.20, 1.0, 0.0], + [0.0, 0.0, 0.05, 0.0, 1.0], + ], +}; + +// =========================================================================== +// Linear algebra primitives — written inline (no new deps) +// =========================================================================== + +/** + * Lower-triangular Cholesky factor of a symmetric positive-semi-definite + * matrix. Throws an Error with `code: 'NON_PSD_CORRELATION'` (and an + * `offending_eigenvalue` hint) if the matrix is not PSD. + * + * Accepts a small numerical tolerance on the diagonal (a tiny negative pivot + * within `tol` is treated as zero) so that degenerate-but-valid matrices do + * not spuriously fail. + */ +export function cholesky(matrix: number[][], tol: number = 1e-10): number[][] { + const n = matrix.length; + if (n === 0) return []; + for (let i = 0; i < n; i++) { + if (matrix[i].length !== n) { + const err = new Error('cholesky: matrix is not square') as Error & { code?: string }; + err.code = 'NON_PSD_CORRELATION'; + throw err; + } + } + + const L: number[][] = Array.from({ length: n }, () => new Array(n).fill(0)); + + for (let i = 0; i < n; i++) { + for (let j = 0; j <= i; j++) { + let sum = 0; + for (let k = 0; k < j; k++) sum += L[i][k] * L[j][k]; + + if (i === j) { + const pivot = matrix[i][i] - sum; + if (pivot < -tol) { + const err = new Error( + `Correlation matrix is not positive semi-definite (pivot ${pivot.toExponential(3)} at row ${i}).`, + ) as Error & { code?: string; offending_eigenvalue?: number }; + err.code = 'NON_PSD_CORRELATION'; + err.offending_eigenvalue = pivot; + throw err; + } + L[i][j] = Math.sqrt(Math.max(0, pivot)); + } else { + const denom = L[j][j]; + if (denom === 0) { + // Zero pivot on a strictly-PSD matrix can occur for a redundant + // dimension; we propagate a 0 row, which produces a zero draw for + // that dimension — acceptable degenerate behaviour. + L[i][j] = 0; + } else { + L[i][j] = (matrix[i][j] - sum) / denom; + } + } + } + } + return L; +} + +/** + * Build a covariance matrix from per-asset stdevs (in percent) and a + * correlation matrix indexed by the same ids. Values are scaled to decimals. + */ +export function buildCovariance( + assetClasses: AssetClass[], + correlation: ReturnCorrelationMatrix, +): number[][] { + const n = assetClasses.length; + const idIndex = new Map(correlation.ids.map((id, i) => [id, i])); + const cov: number[][] = Array.from({ length: n }, () => new Array(n).fill(0)); + for (let i = 0; i < n; i++) { + const ai = idIndex.get(assetClasses[i].id); + if (ai === undefined) { + throw new Error( + `correlation matrix is missing id "${assetClasses[i].id}"`, + ); + } + for (let j = 0; j < n; j++) { + const aj = idIndex.get(assetClasses[j].id); + if (aj === undefined) { + throw new Error( + `correlation matrix is missing id "${assetClasses[j].id}"`, + ); + } + const sigmaI = assetClasses[i].return_stdev_pct / 100; + const sigmaJ = assetClasses[j].return_stdev_pct / 100; + cov[i][j] = correlation.values[ai][aj] * sigmaI * sigmaJ; + } + } + return cov; +} + +// =========================================================================== +// RNG primitives — splitmix32 hash + mulberry32 stream +// =========================================================================== + +/** Stable 32-bit hash combining a seed with a year-index. Used to derive + * per-year RNGs so the draw sequence is independent of trial count. */ +export function hashSeed(seed: number, year: number): number { + let x = (seed | 0) ^ Math.imul(year + 1, 0x9e3779b9); + x = Math.imul(x ^ (x >>> 16), 0x85ebca6b); + x = Math.imul(x ^ (x >>> 13), 0xc2b2ae35); + return (x ^ (x >>> 16)) >>> 0; +} + +/** Mulberry32 PRNG. Matches the engine-wide convention from monte-carlo.ts. */ +export function mulberry32(seed: number): () => number { + let s = seed | 0; + return function next(): number { + s = (s + 0x6d2b79f5) | 0; + let t = Math.imul(s ^ (s >>> 15), 1 | s); + t = (t + Math.imul(t ^ (t >>> 7), 61 | t)) ^ t; + return ((t ^ (t >>> 14)) >>> 0) / 4294967296; + }; +} + +/** Standard-normal draw via Box-Muller; matches monte-carlo.ts gaussian(). */ +export function standardNormal(rng: () => number): number { + const u1 = Math.max(rng(), 1e-10); + const u2 = rng(); + return Math.sqrt(-2 * Math.log(u1)) * Math.cos(2 * Math.PI * u2); +} + +// =========================================================================== +// Inverse normal & inverse Student-t CDFs (for Student-t copula path) +// =========================================================================== + +/** + * Beasley-Springer-Moro approximation of the inverse standard-normal CDF. + * Adequate for our copula use (~6 sig figs in the tails). + */ +export function inverseNormalCDF(p: number): number { + if (p <= 0) return -Infinity; + if (p >= 1) return Infinity; + const a = [ + -3.969683028665376e1, 2.209460984245205e2, -2.759285104469687e2, + 1.38357751867269e2, -3.066479806614716e1, 2.506628277459239, + ]; + const b = [ + -5.447609879822406e1, 1.615858368580409e2, -1.556989798598866e2, + 6.680131188771972e1, -1.328068155288572e1, + ]; + const c = [ + -7.784894002430293e-3, -3.223964580411365e-1, -2.400758277161838, + -2.549732539343734, 4.374664141464968, 2.938163982698783, + ]; + const d = [ + 7.784695709041462e-3, 3.224671290700398e-1, 2.445134137142996, + 3.754408661907416, + ]; + const pLow = 0.02425; + const pHigh = 1 - pLow; + let q: number; + let r: number; + if (p < pLow) { + q = Math.sqrt(-2 * Math.log(p)); + return (((((c[0] * q + c[1]) * q + c[2]) * q + c[3]) * q + c[4]) * q + c[5]) / + ((((d[0] * q + d[1]) * q + d[2]) * q + d[3]) * q + 1); + } + if (p > pHigh) { + q = Math.sqrt(-2 * Math.log(1 - p)); + return -(((((c[0] * q + c[1]) * q + c[2]) * q + c[3]) * q + c[4]) * q + c[5]) / + ((((d[0] * q + d[1]) * q + d[2]) * q + d[3]) * q + 1); + } + q = p - 0.5; + r = q * q; + return (((((a[0] * r + a[1]) * r + a[2]) * r + a[3]) * r + a[4]) * r + a[5]) * q / + (((((b[0] * r + b[1]) * r + b[2]) * r + b[3]) * r + b[4]) * r + 1); +} + +/** Standard-normal CDF (Abramowitz & Stegun 7.1.26 — ~7 sig figs). */ +export function standardNormalCDF(x: number): number { + const sign = x < 0 ? -1 : 1; + const a1 = 0.254829592; + const a2 = -0.284496736; + const a3 = 1.421413741; + const a4 = -1.453152027; + const a5 = 1.061405429; + const p = 0.3275911; + const ax = Math.abs(x) / Math.SQRT2; + const t = 1 / (1 + p * ax); + const y = + 1 - (((((a5 * t + a4) * t) + a3) * t + a2) * t + a1) * t * Math.exp(-ax * ax); + return 0.5 * (1 + sign * y); +} + +/** Lanczos log-gamma — used by the Student-t inverse CDF. */ +function logGamma(x: number): number { + const c = [ + 76.18009172947146, -86.50532032941677, 24.01409824083091, + -1.231739572450155, 0.1208650973866179e-2, -0.5395239384953e-5, + ]; + let y = x; + let tmp = x + 5.5; + tmp -= (x + 0.5) * Math.log(tmp); + let ser = 1.000000000190015; + for (let j = 0; j < 6; j++) { + y += 1; + ser += c[j] / y; + } + return -tmp + Math.log((2.5066282746310005 * ser) / x); +} + +/** Regularised incomplete beta function — continued-fraction approach. */ +function betaIncomplete(a: number, b: number, x: number): number { + if (x <= 0) return 0; + if (x >= 1) return 1; + const lbeta = logGamma(a + b) - logGamma(a) - logGamma(b); + const front = Math.exp(Math.log(x) * a + Math.log(1 - x) * b + lbeta) / a; + // Lentz's algorithm for the continued fraction + let f = 1.0; + let c = 1.0; + let d = 0.0; + for (let i = 0; i < 200; i++) { + const m = i / 2; + let numerator: number; + if (i === 0) numerator = 1; + else if (i % 2 === 0) + numerator = (m * (b - m) * x) / ((a + 2 * m - 1) * (a + 2 * m)); + else + numerator = -(((a + m) * (a + b + m) * x) / ((a + 2 * m) * (a + 2 * m + 1))); + d = 1 + numerator * d; + if (Math.abs(d) < 1e-30) d = 1e-30; + d = 1 / d; + c = 1 + numerator / c; + if (Math.abs(c) < 1e-30) c = 1e-30; + const delta = c * d; + f *= delta; + if (Math.abs(delta - 1) < 1e-12) break; + } + return front * (f - 1); +} + +/** Student-t CDF for x with `dof` degrees of freedom. */ +export function studentTCDF(x: number, dof: number): number { + const v = dof; + const xx = (v) / (v + x * x); + const ib = betaIncomplete(v / 2, 0.5, xx); + if (x >= 0) return 1 - 0.5 * ib; + return 0.5 * ib; +} + +/** + * Inverse Student-t CDF via bisection. Good to ~5 sig figs which is plenty + * for the copula transformation step. + */ +export function inverseStudentT(p: number, dof: number): number { + if (p <= 0) return -Infinity; + if (p >= 1) return Infinity; + let lo = -50; + let hi = 50; + for (let i = 0; i < 60; i++) { + const mid = 0.5 * (lo + hi); + const cdf = studentTCDF(mid, dof); + if (cdf < p) lo = mid; + else hi = mid; + if (hi - lo < 1e-7) break; + } + return 0.5 * (lo + hi); +} + +// =========================================================================== +// Sampler builder +// =========================================================================== + +/** + * Build a deterministic correlated multi-asset return sampler for one MC run. + * + * @param assetClasses Active asset classes — at least one entry. + * @param correlation NxN correlation matrix (validated PSD via Cholesky). + * @param process Distribution choice (LogNormal | StudentT | Bootstrap). + * @param seed Run-level seed; per-year seeds are derived deterministically. + */ +export function buildReturnSampler( + assetClasses: AssetClass[], + correlation: ReturnCorrelationMatrix, + process: ReturnProcess, + seed: number, +): ReturnSampler { + if (assetClasses.length === 0) { + throw new Error('buildReturnSampler: assetClasses must be non-empty'); + } + if (process.kind === 'StudentT' && process.dof <= 2) { + const err = new Error( + `Student-T degrees of freedom must be > 2 for finite variance (got ${process.dof}).`, + ) as Error & { code?: string }; + err.code = 'STUDENT_T_DOF_TOO_LOW'; + throw err; + } + + const n = assetClasses.length; + const means = assetClasses.map((a) => a.expected_return_pct / 100); + const stdevs = assetClasses.map((a) => a.return_stdev_pct / 100); + + if (process.kind === 'Bootstrap') { + // Filter the historical series to the requested window once. + const [yLo, yHi] = process.window; + if (yLo >= yHi) { + const err = new Error( + `bootstrap window invalid: ${yLo}-${yHi}`, + ) as Error & { code?: string }; + err.code = 'BOOTSTRAP_WINDOW_INVALID'; + throw err; + } + const window = SHILLER_SERIES.filter( + (r) => r.year >= yLo && r.year <= yHi, + ); + if (window.length === 0) { + const err = new Error( + `bootstrap window contains no rows (${yLo}-${yHi})`, + ) as Error & { code?: string }; + err.code = 'BOOTSTRAP_WINDOW_INVALID'; + throw err; + } + return { + sample(year: number): Record { + const rng = mulberry32(hashSeed(seed, year)); + const idx = Math.floor(rng() * window.length) % window.length; + const row = window[idx]; + // Map every requested asset class to a column from the historical + // row. Unknown ids fall back to the asset's expected return — we + // always return the same shape regardless of what's bundled. + const out: Record = {}; + for (let i = 0; i < n; i++) { + const id = assetClasses[i].id; + if (id === 'us_equity' || id === 'intl_equity' || id === 'reit') { + out[id] = row.us_equity; + } else if (id === 'us_bond' || id === 'intl_bond') { + out[id] = row.us_bond; + } else if (id === 'cash') { + // Cash is largely unaffected by a bootstrap sample; use the mean. + out[id] = means[i]; + } else { + out[id] = means[i]; + } + } + return out; + }, + }; + } + + // LogNormal & StudentT both use the Cholesky factor of the covariance + // matrix. We compute it once at build time; per-year work is one matrix- + // vector multiply plus inverse-CDF transforms. + const cov = buildCovariance(assetClasses, correlation); + const L = cholesky(cov); + const dof = process.kind === 'StudentT' ? process.dof : 0; + const studentScale = + dof > 2 ? Math.sqrt(dof / (dof - 2)) : 1; + + return { + sample(year: number): Record { + const rng = mulberry32(hashSeed(seed, year)); + // Step 1: n IID standard-normal draws. + const z: number[] = new Array(n); + for (let i = 0; i < n; i++) z[i] = standardNormal(rng); + // Step 2: correlated normals via L*z + const correlatedZ: number[] = new Array(n).fill(0); + for (let i = 0; i < n; i++) { + for (let j = 0; j <= i; j++) { + correlatedZ[i] += L[i][j] * z[j]; + } + } + // Step 3: convert to per-asset returns. + const out: Record = {}; + for (let i = 0; i < n; i++) { + const id = assetClasses[i].id; + if (process.kind === 'LogNormal') { + // Standard log-normal: returnNominal = exp(mu_ln + sigma_ln * z) - 1 + // where the (mu, sigma) for the log space are derived from arith + // mean + stdev per-class. We re-derive sigma_ln per class instead + // of using the correlated covariance directly — the Cholesky- + // driven correlatedZ supplies the standard-normal shock with the + // requested cross-correlation, but each class still uses its own + // parameters in log-space so the marginal mean/stdev are right. + const ratio = stdevs[i] / (1 + means[i]); + const sigmaLn = Math.sqrt(Math.log(1 + ratio * ratio)); + const muLn = Math.log(1 + means[i]) - (sigmaLn * sigmaLn) / 2; + // Normalise the correlated draw to unit stdev then apply per-class + // log-normal parameters. Without this, the log-normal mean drifts + // because we'd be exponentiating a draw with the wrong scale. + const unitScale = stdevs[i] > 0 ? correlatedZ[i] / stdevs[i] : 0; + out[id] = Math.exp(muLn + sigmaLn * unitScale) - 1; + } else { + // Student-t copula: convert correlated normal to a uniform via + // standard-normal CDF, then push through inverse Student-t CDF + // with the configured dof. Final scaling preserves the requested + // arithmetic mean and stdev. + const unitScale = stdevs[i] > 0 ? correlatedZ[i] / stdevs[i] : 0; + const u = standardNormalCDF(unitScale); + // Clamp away from 0/1 to avoid +/-Infinity from the inverse CDF. + const uClamped = Math.min(0.9999999, Math.max(1e-7, u)); + const tDraw = inverseStudentT(uClamped, dof); + // Student-t with dof has variance dof/(dof-2); divide to standardise. + const standardised = tDraw / studentScale; + out[id] = means[i] + stdevs[i] * standardised; + } + // Floor return at -1 (cannot lose more than 100%). + if (out[id] < -1) out[id] = -1; + } + return out; + }, + }; +} + +// =========================================================================== +// Joint return-inflation sampler (Section 4 of the v0.4 task list) +// =========================================================================== + +/** + * Couple a return sampler with an inflation shock dimension via an augmented + * covariance matrix. Equity classes use `return_inflation_correlation`; bond + * classes use `bond_inflation_correlation`; cash is treated as orthogonal to + * the inflation shock. + * + * The augmented matrix is PSD-validated via Cholesky at build time; failure + * throws an Error with code 'NON_PSD_CORRELATION'. + * + * Per year the sampler produces: + * - a Record of nominal returns + * - a single inflation rate (decimal) — the AR(1) recurrence applied if + * the supplied `InflationProcess` is AR(1); flat passthrough otherwise. + */ +export interface JointSample { + returns: Record; + inflation: number; +} + +export interface JointSampler { + sample(year: number, priorInflation: number): JointSample; +} + +const EQUITY_LIKE = new Set([ + 'us_equity', + 'intl_equity', + 'reit', + 'commodities', +]); +const BOND_LIKE = new Set(['us_bond', 'intl_bond']); + +function inflationCorrelationFor( + id: AssetClassId, + returnInflationCorr: number, + bondInflationCorr: number, +): number { + if (BOND_LIKE.has(id)) return bondInflationCorr; + if (EQUITY_LIKE.has(id)) return returnInflationCorr; + // Cash and unknown classes: treated as uncorrelated with inflation shocks. + return 0; +} + +/** + * Build an augmented (returns + inflation shock) sampler. Falls back to two + * independent samplers when the supplied inflation process is Flat (since + * Flat draws no randomness and therefore has nothing to correlate). + */ +export function buildJointSampler( + assetClasses: AssetClass[], + correlation: ReturnCorrelationMatrix, + process: ReturnProcess, + inflationProcess: + | { kind: 'Flat'; rate_pct: number } + | { kind: 'AR1'; long_run_mean_pct: number; phi: number; shock_stdev_pct: number; initial_pct: number }, + returnInflationCorr: number, + bondInflationCorr: number, + seed: number, +): JointSampler { + // Flat inflation: no augmentation needed; delegate to the plain return + // sampler and constant inflation. + if (inflationProcess.kind === 'Flat') { + const rs = buildReturnSampler(assetClasses, correlation, process, seed); + const flatRate = inflationProcess.rate_pct / 100; + return { + sample(year: number, _priorInflation: number): JointSample { + return { returns: rs.sample(year), inflation: flatRate }; + }, + }; + } + + // Bootstrap path: historical years carry their own inflation, so the + // "joint" draw simply pairs a sampled history-year's CPI with its returns. + if (process.kind === 'Bootstrap') { + const [yLo, yHi] = process.window; + const window = SHILLER_SERIES.filter( + (r) => r.year >= yLo && r.year <= yHi, + ); + if (window.length === 0) { + const err = new Error( + `bootstrap window contains no rows (${yLo}-${yHi})`, + ) as Error & { code?: string }; + err.code = 'BOOTSTRAP_WINDOW_INVALID'; + throw err; + } + const n = assetClasses.length; + const means = assetClasses.map((a) => a.expected_return_pct / 100); + return { + sample(year: number, _priorInflation: number): JointSample { + const rng = mulberry32(hashSeed(seed, year)); + const idx = Math.floor(rng() * window.length) % window.length; + const row = window[idx]; + const ret: Record = {}; + for (let i = 0; i < n; i++) { + const id = assetClasses[i].id; + if (EQUITY_LIKE.has(id)) ret[id] = row.us_equity; + else if (BOND_LIKE.has(id)) ret[id] = row.us_bond; + else if (id === 'cash') ret[id] = means[i]; + else ret[id] = means[i]; + } + return { returns: ret, inflation: row.us_cpi }; + }, + }; + } + + // LogNormal / StudentT path: build the augmented covariance and Cholesky- + // factorise once. The extra dimension is the inflation shock. + const n = assetClasses.length; + const means = assetClasses.map((a) => a.expected_return_pct / 100); + const stdevs = assetClasses.map((a) => a.return_stdev_pct / 100); + const inflStd = inflationProcess.shock_stdev_pct / 100; + const inflMean = inflationProcess.long_run_mean_pct / 100; + const phi = inflationProcess.phi; + + // Build (n+1)x(n+1) covariance matrix. Top-left n x n block is the + // existing return covariance; the extra row/column couples each asset to + // the inflation shock. + const baseCov = buildCovariance(assetClasses, correlation); + const aug: number[][] = baseCov.map((row) => [...row, 0]); + aug.push(new Array(n + 1).fill(0)); + for (let i = 0; i < n; i++) { + const corr = inflationCorrelationFor( + assetClasses[i].id, + returnInflationCorr, + bondInflationCorr, + ); + const cov = corr * stdevs[i] * inflStd; + aug[i][n] = cov; + aug[n][i] = cov; + } + aug[n][n] = inflStd * inflStd; + + const L = cholesky(aug); // throws NON_PSD_CORRELATION if augmented matrix invalid + const dof = process.kind === 'StudentT' ? process.dof : 0; + const studentScale = dof > 2 ? Math.sqrt(dof / (dof - 2)) : 1; + + return { + sample(year: number, priorInflation: number): JointSample { + const rng = mulberry32(hashSeed(seed, year)); + const z: number[] = new Array(n + 1); + for (let i = 0; i < n + 1; i++) z[i] = standardNormal(rng); + const correlatedZ: number[] = new Array(n + 1).fill(0); + for (let i = 0; i < n + 1; i++) { + for (let j = 0; j <= i; j++) correlatedZ[i] += L[i][j] * z[j]; + } + // Returns + const ret: Record = {}; + for (let i = 0; i < n; i++) { + const id = assetClasses[i].id; + const unit = stdevs[i] > 0 ? correlatedZ[i] / stdevs[i] : 0; + let r: number; + if (process.kind === 'LogNormal') { + const ratio = stdevs[i] / (1 + means[i]); + const sigmaLn = Math.sqrt(Math.log(1 + ratio * ratio)); + const muLn = Math.log(1 + means[i]) - (sigmaLn * sigmaLn) / 2; + r = Math.exp(muLn + sigmaLn * unit) - 1; + } else { + const u = standardNormalCDF(unit); + const uClamped = Math.min(0.9999999, Math.max(1e-7, u)); + const tDraw = inverseStudentT(uClamped, dof); + r = means[i] + stdevs[i] * (tDraw / studentScale); + } + if (r < -1) r = -1; + ret[id] = r; + } + // Inflation: AR(1) using the augmented draw as epsilon (already scaled + // by inflStd via the Cholesky factor — divide back out so we have a + // unit-stdev shock, then multiply by inflStd for the recurrence). + const epsUnit = inflStd > 0 ? correlatedZ[n] / inflStd : 0; + const epsilon = inflStd * epsUnit; // == correlatedZ[n] but kept explicit + const nextInflation = inflMean + phi * (priorInflation - inflMean) + epsilon; + return { returns: ret, inflation: nextInflation }; + }, + }; +} diff --git a/src/risk-metrics.ts b/src/risk-metrics.ts new file mode 100644 index 0000000..f8213c3 --- /dev/null +++ b/src/risk-metrics.ts @@ -0,0 +1,236 @@ +/** + * Risk Metrics Module (ADR-033 / CONTRACT-018) + * + * Pure function computing VaR, CVaR, Sortino, drawdown and per-year + * quantile trajectories from a completed MCResult. + * + * Invariants (from CONTRACT-018): + * - Pure: identical inputs -> identical outputs. + * - cvar_95 <= var_95 (and the analogous monotonicity at 99%). The MC + * convention expresses VaR as "quantile of the terminal-balance + * distribution", so "breach" means "below"; CVaR is the conditional + * mean *below* VaR, which is always <= VaR. + * - max_drawdown_pct in [0, 100]. + * - p10 <= p50 <= p90 per year. + * + * The module does not import MCResult from monte-carlo.ts to avoid a + * circular dependency; we describe the minimum shape we need via an inline + * type and let monte-carlo.ts pass the richer value in. + */ + +import type { Scenario, RiskMetrics } from './types'; + +// =========================================================================== +// Minimal MCResult-shape we rely on (decoupled to avoid circular imports) +// =========================================================================== + +export interface MCRiskInputs { + /** Real terminal balance per trial. */ + terminal_distribution: number[]; + /** + * Per-trial year-by-year real balance paths. Shape: paths[trialIdx][yearIdx]. + * Paths may be different lengths (stochastic longevity terminates some + * trials early); the utility functions tolerate ragged input. + */ + real_balance_paths: number[][]; + /** Per-trial annualised realised portfolio return, decimal. */ + annualised_returns: number[]; +} + +// =========================================================================== +// Quantile helpers +// =========================================================================== + +/** Sort-free safe percentile: expects `sorted` in ascending order. */ +export function quantile(sorted: number[], p: number): number { + if (sorted.length === 0) return 0; + if (p <= 0) return sorted[0]; + if (p >= 1) return sorted[sorted.length - 1]; + const idx = p * (sorted.length - 1); + const lo = Math.floor(idx); + const hi = Math.ceil(idx); + if (lo === hi) return sorted[lo]; + const frac = idx - lo; + return sorted[lo] * (1 - frac) + sorted[hi] * frac; +} + +// =========================================================================== +// Drawdown utilities +// =========================================================================== + +export interface DrawdownStats { + /** Maximum peak-to-trough drop as a percentage of peak (0-100). */ + maxDrawdownPct: number; + /** Years from trough back to prior peak; Infinity if not recovered. */ + recoveryYears: number; +} + +/** + * Peak-to-trough drawdown on a real balance path. The "worst drawdown" is the + * deepest percentage drop from any prior peak. + */ +export function computeDrawdown(path: number[]): DrawdownStats { + if (path.length === 0) return { maxDrawdownPct: 0, recoveryYears: Infinity }; + let peak = path[0]; + let peakIdx = 0; + let worstDDPct = 0; + let worstTroughIdx = 0; + let worstPeakAtTrough = peak; + for (let i = 0; i < path.length; i++) { + if (path[i] > peak) { + peak = path[i]; + peakIdx = i; + } + if (peak <= 0) continue; + const ddPct = ((peak - path[i]) / peak) * 100; + if (ddPct > worstDDPct) { + worstDDPct = ddPct; + worstTroughIdx = i; + worstPeakAtTrough = peak; + } + } + void peakIdx; // not reported directly + let recoveryYears: number = Infinity; + for (let j = worstTroughIdx + 1; j < path.length; j++) { + if (path[j] >= worstPeakAtTrough) { + recoveryYears = j - worstTroughIdx; + break; + } + } + return { maxDrawdownPct: worstDDPct, recoveryYears }; +} + +// =========================================================================== +// Sortino +// =========================================================================== + +/** + * Sortino ratio on realised annualised portfolio returns. Downside deviation + * is the root-mean-square of shortfalls below the minimum-acceptable return + * (MAR). We return NaN when downside deviation is zero (signalled by the + * caller as "N/A"); ADR-033 permits either Infinity or NaN. + */ +export function computeSortino( + returns: number[], + mar: number, +): number { + if (returns.length === 0) return NaN; + const mean = returns.reduce((a, b) => a + b, 0) / returns.length; + let downsideSq = 0; + let downsideCount = 0; + for (const r of returns) { + if (r < mar) { + downsideSq += (mar - r) * (mar - r); + downsideCount++; + } + } + if (downsideCount === 0) return NaN; + const downsideDev = Math.sqrt(downsideSq / returns.length); + if (downsideDev === 0) return NaN; + return (mean - mar) / downsideDev; +} + +// =========================================================================== +// Main entry point +// =========================================================================== + +/** + * Compute the full RiskMetrics record from a completed Monte Carlo result. + * + * Scenario provides: + * - risk_free_rate_pct (default 3.5) — used as both the Sortino MAR and + * the risk-free baseline for the VaR/CVaR return-space expression. + * - current_age / end_age — span used to compound VaR from terminal space + * back into an annual return. + */ +export function computeRiskMetrics( + inputs: MCRiskInputs, + scenario: Scenario, +): RiskMetrics { + const { terminal_distribution, real_balance_paths, annualised_returns } = inputs; + const mar = (scenario.risk_free_rate_pct ?? 3.5) / 100; + + const sortedTerminals = [...terminal_distribution].sort((a, b) => a - b); + const var95 = quantile(sortedTerminals, 0.05); + const var99 = quantile(sortedTerminals, 0.01); + + // CVaR: mean of the worst-k% tail. For CVaR-95 the "worst 5%" is the + // bottom 5% of the sorted terminal distribution. + const n = sortedTerminals.length; + const tail95Count = Math.max(1, Math.ceil(n * 0.05)); + const tail99Count = Math.max(1, Math.ceil(n * 0.01)); + const cvar95 = + sortedTerminals.slice(0, tail95Count).reduce((a, b) => a + b, 0) / tail95Count; + const cvar99 = + sortedTerminals.slice(0, tail99Count).reduce((a, b) => a + b, 0) / tail99Count; + + // VaR expressed as an annualised real return, starting from current_balance. + // If we can't sensibly invert (non-positive terminal or starting balance), + // leave the field as 0 — callers should prefer the terminal-balance form + // in that case. Horizon uses years (end_age - current_age). + const years = Math.max(1, scenario.end_age - scenario.current_age); + const startingBalance = Math.max(1, scenario.current_balance); + let var95ReturnPct = 0; + if (var95 > 0 && startingBalance > 0) { + const ratio = var95 / startingBalance; + var95ReturnPct = (Math.pow(ratio, 1 / years) - 1) * 100; + } + + // Sortino + const sortino = computeSortino(annualised_returns, mar); + + // Drawdowns per trial + const drawdowns = real_balance_paths.map(computeDrawdown); + const maxDD = drawdowns.reduce( + (acc, dd) => (dd.maxDrawdownPct > acc.maxDrawdownPct ? dd : acc), + { maxDrawdownPct: 0, recoveryYears: Infinity }, + ); + const sortedDDPct = drawdowns.map((dd) => dd.maxDrawdownPct).sort((a, b) => a - b); + const medianDDPct = quantile(sortedDDPct, 0.5); + // Recovery: median across trials (Infinity treated as end-of-list — if + // >=50% of trials never recovered, the median is Infinity). + const recoveryYears = drawdowns + .map((dd) => dd.recoveryYears) + .sort((a, b) => a - b); + const medianRecovery = + recoveryYears.length === 0 + ? Infinity + : recoveryYears[Math.floor(recoveryYears.length / 2)]; + + // Per-year P10 / P50 / P90 balance trajectories. We collect all trials' + // values at each year index; trials that terminated early contribute only + // to the years they actually spanned, which is the right reporting shape. + const maxLen = real_balance_paths.reduce( + (m, p) => Math.max(m, p.length), + 0, + ); + const p10Path: number[] = []; + const p50Path: number[] = []; + const p90Path: number[] = []; + for (let yr = 0; yr < maxLen; yr++) { + const ys: number[] = []; + for (const p of real_balance_paths) { + if (yr < p.length) ys.push(p[yr]); + } + ys.sort((a, b) => a - b); + p10Path.push(quantile(ys, 0.1)); + p50Path.push(quantile(ys, 0.5)); + p90Path.push(quantile(ys, 0.9)); + } + + return { + var_95_terminal_real: var95, + var_99_terminal_real: var99, + cvar_95_terminal_real: cvar95, + cvar_99_terminal_real: cvar99, + var_95_return_pct: var95ReturnPct, + sortino_ratio: sortino, + max_drawdown_pct: maxDD.maxDrawdownPct, + max_drawdown_recovery_years: maxDD.recoveryYears, + median_drawdown_pct: medianDDPct, + median_drawdown_recovery_years: medianRecovery, + p10_year_by_year_balance_real: p10Path, + p50_year_by_year_balance_real: p50Path, + p90_year_by_year_balance_real: p90Path, + }; +} diff --git a/src/sensitivity.ts b/src/sensitivity.ts index b4abec0..85920ca 100644 --- a/src/sensitivity.ts +++ b/src/sensitivity.ts @@ -70,15 +70,25 @@ export function runSensitivityAnalysis( const log = getLogger(); log.info('Starting sensitivity analysis', { parameterCount: PARAMETERS.length }); + // v0.4 (ADR-027/031): force-disable stochastic features so the tornado + // chart reflects single-parameter sensitivity in isolation. + const baseScenario: Scenario = { + ...scenario, + inflation_model: 'Flat' as const, + longevity_model: 'Fixed' as const, + return_distribution_kind: 'LogNormal' as const, + asset_classes: [], + }; + const factors: SensitivityFactor[] = []; for (const param of PARAMETERS) { // Skip withdrawal_pct when strategy is Age-Banded (not applicable) - if (param.name === 'withdrawal_pct' && scenario.withdrawal_strategy === 'Age-Banded') { + if (param.name === 'withdrawal_pct' && baseScenario.withdrawal_strategy === 'Age-Banded') { continue; } - const baselineValue = (scenario as Record)[param.name] as number; + const baselineValue = (baseScenario as Record)[param.name] as number; // Compute absolute delta const absDelta = param.deltaIsPct @@ -91,8 +101,8 @@ export function runSensitivityAnalysis( // --- Clamping guards --- if (param.name === 'retirement_age') { // retirement_age must stay in (current_age, end_age) - lowValue = clamp(Math.round(lowValue), scenario.current_age + 1, scenario.end_age - 1); - highValue = clamp(Math.round(highValue), scenario.current_age + 1, scenario.end_age - 1); + lowValue = clamp(Math.round(lowValue), baseScenario.current_age + 1, baseScenario.end_age - 1); + highValue = clamp(Math.round(highValue), baseScenario.current_age + 1, baseScenario.end_age - 1); } else if ( param.name === 'nominal_return_pct' || param.name === 'inflation_pct' || @@ -109,7 +119,7 @@ export function runSensitivityAnalysis( } // Run projection with low value - const lowScenario = cloneScenario(scenario); + const lowScenario = cloneScenario(baseScenario); (lowScenario as Record)[param.name] = lowValue; // ADR-027: sensitivity must always exclude the Black Swan stress event so // bands are interpretable in isolation. @@ -117,7 +127,7 @@ export function runSensitivityAnalysis( const lowResult = projectionFn(lowScenario); // Run projection with high value - const highScenario = cloneScenario(scenario); + const highScenario = cloneScenario(baseScenario); (highScenario as Record)[param.name] = highValue; highScenario.black_swan_enabled = false; const highResult = projectionFn(highScenario); diff --git a/src/types.ts b/src/types.ts index 20066d9..b29765e 100644 --- a/src/types.ts +++ b/src/types.ts @@ -333,6 +333,48 @@ export interface Scenario { // Withdrawal order withdrawal_order: string; + + // -------------------------------------------------------------------------- + // v0.4 stochastic foundation (CONTRACT-018) — all optional, all defaulted + // downstream. A v0.3-shape Scenario continues to typecheck and compute + // identical outputs. + // -------------------------------------------------------------------------- + + // Multi-asset (ADR-030) + asset_classes?: AssetClass[]; + return_correlation_matrix?: ReturnCorrelationMatrix | null; + return_distribution_kind?: 'LogNormal' | 'StudentT' | 'Bootstrap'; + return_distribution_dof?: number; + bootstrap_window?: [number, number]; + + // Inflation (ADR-031) + inflation_model?: 'Flat' | 'AR1'; + inflation_long_run_mean_pct?: number; + inflation_ar1_phi?: number; + inflation_shock_stdev_pct?: number; + inflation_initial_pct?: number; + inflation_calibration_preset?: + | 'US-CPI' + | 'UK-CPI' + | 'UK-RPI' + | 'EU-HICP' + | 'Custom'; + return_inflation_correlation?: number; + bond_inflation_correlation?: number; + + // Longevity (ADR-032) + longevity_model?: 'Fixed' | 'Gompertz' | 'Cohort'; + longevity_modal_age?: number; + longevity_dispersion?: number; + longevity_cohort_country?: 'US' | 'UK'; + sex?: Sex; + longevity_partner_modal_age?: number | null; + longevity_partner_dispersion?: number; + longevity_partner_cohort_country?: 'US' | 'UK'; + partner_sex?: Sex; + + // Risk metrics (ADR-033) + risk_free_rate_pct?: number; } // --------------------------------------------------------------------------- @@ -382,6 +424,15 @@ export interface TimelineRow { black_swan_loss: number; /** Tag indicating which strategy event produced this year's withdrawal. */ withdrawal_event: WithdrawalEvent; + + // CONTRACT-018 / ADR-030..033 additions (v0.4). All defaulted; existing + // consumers that do not inspect these fields are unaffected. + /** Realised inflation rate for this year (decimal, e.g. 0.025 == 2.5%). */ + inflation_this_year: number; + /** Per-asset realised return by AssetClassId. Null in single-asset mode. */ + asset_returns: Record | null; + /** True on the final row of an MC trial that was terminated by a stochastic longevity draw (not at end_age). */ + death_sampled_this_trial?: boolean; } export interface FanChartRow { @@ -404,3 +455,132 @@ export interface Metrics { total_taxes: number; estate_value: number; } + +// --------------------------------------------------------------------------- +// v0.4 — Stochastic Foundation (ADR-029 through ADR-033, CONTRACT-018) +// --------------------------------------------------------------------------- + +/** + * A stable identifier for an asset class within a scenario. Common canonical + * values are listed first; arbitrary user-defined ids are also permitted. + */ +export type AssetClassId = + | 'us_equity' + | 'intl_equity' + | 'us_bond' + | 'intl_bond' + | 'reit' + | 'commodities' + | 'cash' + | string; + +/** + * One row of the multi-asset portfolio. + * Weights across all AssetClasses in a Scenario must sum to 100 +/- 1. + */ +export interface AssetClass { + id: AssetClassId; + name: string; + expected_return_pct: number; + return_stdev_pct: number; + weight_pct: number; +} + +/** + * Square symmetric positive-semi-definite matrix indexed by AssetClassId. + * Validated at the engine boundary via attempted Cholesky factorisation. + */ +export interface ReturnCorrelationMatrix { + ids: AssetClassId[]; + values: number[][]; +} + +/** + * Sex designation used by the Gompertz and cohort longevity models. + */ +export type Sex = 'M' | 'F' | 'Unspecified'; + +/** + * Discriminated union describing the sampling process for annual returns. + * - LogNormal: current v0.3 behaviour (multivariate log-normal). + * - StudentT: Student-t copula with configurable degrees of freedom. + * - Bootstrap: resamples simultaneous tuples from a bundled historical series. + */ +export type ReturnProcess = + | { kind: 'LogNormal' } + | { kind: 'StudentT'; dof: number } + | { kind: 'Bootstrap'; window: [number, number] }; + +/** + * Stochastic inflation process. Flat is a no-op passthrough; AR(1) follows the + * recurrence defined in ADR-031. + */ +export type InflationProcess = + | { kind: 'Flat'; rate_pct: number } + | { + kind: 'AR1'; + long_run_mean_pct: number; + phi: number; + shock_stdev_pct: number; + initial_pct: number; + }; + +/** + * Longevity sampling model (ADR-032 / DDD-010). + */ +export type LongevityModel = + | { kind: 'Fixed'; end_age: number } + | { kind: 'Gompertz'; modal_age: number; dispersion: number; sex?: Sex } + | { + kind: 'Cohort'; + country: 'US' | 'UK'; + sex: Sex; + birth_year: number; + }; + +/** + * Produces a correlated tuple of asset returns for a given year index. + * Must be deterministic given the build-time seed. + */ +export interface ReturnSampler { + sample(year: number): Record; +} + +/** + * Produces an inflation rate (decimal, not percent) for a given year. The + * AR(1) variant consumes one draw per call; Flat is a no-op. + */ +export interface InflationSampler { + sample(year: number, priorInflation: number): number; + readonly kind: 'Flat' | 'AR1'; +} + +/** + * Samples a death age from the configured longevity model. `median` and + * `survival` are deterministic and do not consume randomness. + */ +export interface LongevitySampler { + sample(current_age: number): number; + median(current_age: number): number; + survival(age: number, current_age: number): number; + readonly kind: 'Fixed' | 'Gompertz' | 'Cohort'; +} + +/** + * Institutional risk metrics attached to a Monte Carlo result (ADR-033). + */ +export interface RiskMetrics { + var_95_terminal_real: number; + var_99_terminal_real: number; + cvar_95_terminal_real: number; + cvar_99_terminal_real: number; + var_95_return_pct: number; + sortino_ratio: number; + max_drawdown_pct: number; + max_drawdown_recovery_years: number; + median_drawdown_pct: number; + median_drawdown_recovery_years: number; + p10_year_by_year_balance_real: number[]; + p50_year_by_year_balance_real: number[]; + p90_year_by_year_balance_real: number[]; +}