Coverage for src/jquantstats/_stats/_montecarlo.py: 100%
71 statements
« prev ^ index » next coverage.py v7.17.1, created at 2026-11-03 15:52 +0000
« prev ^ index » next coverage.py v7.17.1, created at 2026-11-03 15:52 +0000
1"""Monte Carlo simulation statistics."""
3from __future__ import annotations
5import math
6from typing import TYPE_CHECKING
8import numpy as np
9import polars as pl
11from ..exceptions import NonPositivePeriodsPerYearError, NonPositiveWindowError
13if TYPE_CHECKING:
14 from ..data import Data
17class _MonteCarloStatsMixin:
18 """Mixin providing Monte Carlo simulation statistics for return series."""
20 _data: Data
21 all: pl.DataFrame
23 if TYPE_CHECKING:
24 from .._protocol import DataLike
26 data: DataLike
28 @staticmethod
29 def _validate_positive_integer(name: str, value: int) -> int:
30 """Validate that *value* is a positive integer."""
31 if isinstance(value, bool) or not isinstance(value, int) or value <= 0:
32 raise NonPositiveWindowError(name)
33 return value
35 @staticmethod
36 def _block_bootstrap_paths(values: np.ndarray, n: int, period: int, block_size: int) -> np.ndarray:
37 """Generate *n* return paths via circular block bootstrap, returning an ``(n, period)`` array.
39 Block starts are drawn uniformly over all observations and a block that
40 runs past the last observation wraps around to the first (Politis and
41 Romano, 1992). Every observation is therefore drawn with equal
42 probability. A non-circular moving-block bootstrap would draw the first
43 and last observations at ``1 / block_size`` of the rate of interior
44 ones, biasing any statistic of a series whose behaviour changes through
45 the sample.
47 All *n* paths are sampled simultaneously using fully vectorised numpy
48 indexing — no Python-level loop is required.
50 Args:
51 values: 1-D array of historical returns.
52 n: Number of simulation paths.
53 period: Length of each simulated path (number of return observations).
54 block_size: Block length for the block bootstrap resampling.
56 Returns:
57 np.ndarray: Float64 array of shape ``(n, period)``.
59 """
60 n_obs = values.size
61 if n_obs == 0: # pragma: no cover
62 return np.full((n, period), np.nan, dtype=np.float64)
64 n_blocks = math.ceil(period / block_size)
65 # Draw all block starting indices at once: shape (n, n_blocks)
66 starts = np.random.randint(0, n_obs, size=(n, n_blocks))
67 # Build full index array: (n, n_blocks, block_size), wrapping past the end
68 idx = (starts[:, :, np.newaxis] + np.arange(block_size)[np.newaxis, np.newaxis, :]) % n_obs
69 # Flatten and trim to the requested period length
70 return np.asarray(values[idx].reshape(n, -1)[:, :period])
72 def _simulate_distribution(self, n: int, period: int) -> dict[str, np.ndarray]:
73 """Prepare validated inputs and sample *n* block-bootstrap paths per asset.
75 Returns a dict mapping each asset column name to an ``(n, period)``
76 float64 array of simulated return paths. Assets with no usable data
77 map to an ``(n, period)`` array of NaN.
79 Args:
80 n: Number of simulation paths (must be a positive integer).
81 period: Path length in return observations (must be a positive integer).
83 Returns:
84 dict[str, np.ndarray]: Asset name → ``(n, period)`` paths array.
86 """
87 n = self._validate_positive_integer("n", n)
88 period = self._validate_positive_integer("period", period)
89 block_size = max(1, round(period**0.5))
91 paths: dict[str, np.ndarray] = {}
92 for col, series in self._data.items():
93 clean = series.cast(pl.Float64).drop_nulls().drop_nans()
94 values = np.asarray(clean.to_numpy(), dtype=np.float64)
95 if values.size == 0:
96 paths[col] = np.full((n, period), np.nan, dtype=np.float64)
97 else:
98 block = min(block_size, values.size)
99 paths[col] = self._block_bootstrap_paths(values, n, period, block)
100 return paths
102 def montecarlo(self, n: int = 1000, period: int = 252) -> pl.DataFrame:
103 """Simulate cumulative returns across *n* block-bootstrap paths.
105 Args:
106 n: Number of Monte Carlo paths. Defaults to 1000.
107 period: Simulation horizon in return observations. Defaults to 252.
109 Returns:
110 pl.DataFrame: Shape ``(n, n_assets)`` — one simulated terminal
111 cumulative return per path and asset.
114 Returns NaN when:
115 Entries are NaN for assets with no usable (non-null, non-NaN)
116 observations.
117 """
118 paths = self._simulate_distribution(n=n, period=period)
119 result = {col: np.prod(1.0 + arr, axis=1) - 1.0 for col, arr in paths.items()}
120 return pl.DataFrame(result)
122 def montecarlo_sharpe(
123 self,
124 n: int = 1000,
125 period: int = 252,
126 periods_per_year: int | float | None = None,
127 ) -> pl.DataFrame:
128 """Simulate the Sharpe-ratio distribution across block-bootstrap paths.
130 Args:
131 n: Number of Monte Carlo paths. Defaults to 1000.
132 period: Simulation horizon in return observations. Defaults to 252.
133 periods_per_year: Annualisation factor. Defaults to the value
134 inferred from the data.
136 Returns:
137 pl.DataFrame: Shape ``(n, n_assets)`` — one simulated annualised
138 Sharpe ratio per path and asset.
141 Returns NaN when:
142 Entries are NaN when a path's standard deviation is zero or the asset
143 has no usable observations.
144 """
145 ppy = self._data._periods_per_year if periods_per_year is None else periods_per_year
146 if ppy <= 0:
147 raise NonPositivePeriodsPerYearError
148 scale = math.sqrt(ppy)
149 paths = self._simulate_distribution(n=n, period=period)
150 result: dict[str, np.ndarray] = {}
151 for col, arr in paths.items():
152 means = arr.mean(axis=1)
153 stds = arr.std(axis=1, ddof=1)
154 with np.errstate(invalid="ignore", divide="ignore"):
155 result[col] = np.where(stds == 0.0, np.nan, means / stds * scale)
156 return pl.DataFrame(result)
158 def montecarlo_drawdown(self, n: int = 1000, period: int = 252) -> pl.DataFrame:
159 """Simulate the maximum-drawdown distribution across block-bootstrap paths.
161 Args:
162 n: Number of Monte Carlo paths. Defaults to 1000.
163 period: Simulation horizon in return observations. Defaults to 252.
165 Returns:
166 pl.DataFrame: Shape ``(n, n_assets)`` — one simulated maximum
167 drawdown per path and asset (values in ``[-1, 0]``).
170 Returns NaN when:
171 Entries are NaN for assets with no usable (non-null, non-NaN)
172 observations.
173 """
174 paths = self._simulate_distribution(n=n, period=period)
175 result: dict[str, np.ndarray] = {}
176 for col, arr in paths.items():
177 nav = np.cumprod(1.0 + arr, axis=1)
178 hwm = np.maximum.accumulate(nav, axis=1)
179 result[col] = np.min(nav / hwm - 1.0, axis=1)
180 return pl.DataFrame(result)
182 def montecarlo_cagr(
183 self,
184 n: int = 1000,
185 period: int = 252,
186 periods_per_year: int | float | None = None,
187 ) -> pl.DataFrame:
188 """Simulate the CAGR distribution across block-bootstrap paths.
190 Args:
191 n: Number of Monte Carlo paths. Defaults to 1000.
192 period: Simulation horizon in return observations. Defaults to 252.
193 periods_per_year: Annualisation factor. Defaults to the value
194 inferred from the data.
196 Returns:
197 pl.DataFrame: Shape ``(n, n_assets)`` — one simulated annualised
198 CAGR per path and asset.
201 Returns NaN when:
202 Entries are NaN when a path's terminal compound return is non-positive
203 or the asset has no usable observations.
204 """
205 ppy = self._data._periods_per_year if periods_per_year is None else periods_per_year
206 if ppy <= 0:
207 raise NonPositivePeriodsPerYearError
208 years = period / ppy
209 paths = self._simulate_distribution(n=n, period=period)
210 result: dict[str, np.ndarray] = {}
211 for col, arr in paths.items():
212 totals = np.prod(1.0 + arr, axis=1)
213 with np.errstate(invalid="ignore"):
214 result[col] = np.where(totals > 0, totals ** (1.0 / years) - 1.0, np.nan)
215 return pl.DataFrame(result)