Coverage for src/cvxcla/types.py: 100%
113 statements
« prev ^ index » next coverage.py v7.16.2, created at 2026-09-29 05:24 +0000
1"""Type definitions and classes for the Critical Line Algorithm.
3This module defines the core data structures used in the Critical Line Algorithm:
4- FrontierPoint: Represents a point on the efficient frontier.
5- TurningPoint: Represents a turning point on the efficient frontier.
6- Frontier: Represents the entire efficient frontier.
8It also defines type aliases for commonly used types.
9"""
11from __future__ import annotations
13from collections.abc import Iterator
14from dataclasses import dataclass, field
15from typing import TYPE_CHECKING
17if TYPE_CHECKING:
18 import plotly.graph_objects as go
20import numpy as np
21from numpy.typing import NDArray
23from .operators import CovarianceOperator
26def _covariance_matvec(
27 covariance: NDArray[np.float64] | CovarianceOperator, x: NDArray[np.float64]
28) -> NDArray[np.float64]:
29 """Compute ``Sigma @ x`` for a dense matrix or a ``CovarianceOperator`` backend."""
30 if isinstance(covariance, CovarianceOperator):
31 return covariance.matvec(x)
32 return covariance @ x
35@dataclass(frozen=True)
36class FrontierPoint:
37 """A point on the efficient frontier.
39 This class represents a portfolio on the efficient frontier, defined by its weights.
40 It provides methods to compute the expected return and variance of the portfolio.
42 Attributes:
43 weights: Vector of portfolio weights for each asset.
45 """
47 weights: NDArray[np.float64]
49 def mean(self, mean: NDArray[np.float64]) -> float:
50 """Compute the expected return of the portfolio.
52 Args:
53 mean: Vector of expected returns for each asset.
55 Returns:
56 The expected return of the portfolio.
58 Examples:
59 >>> import numpy as np
60 >>> fp = FrontierPoint(weights=np.array([0.5, 0.5]))
61 >>> fp.mean(np.array([0.1, 0.2]))
62 0.15000000000000002
64 """
65 return float(mean.T @ self.weights)
67 def variance(self, covariance: NDArray[np.float64] | CovarianceOperator) -> float:
68 """Compute the expected variance of the portfolio.
70 Args:
71 covariance: Covariance matrix of asset returns, either as a dense
72 matrix or as a ``CovarianceOperator`` backend.
74 Returns:
75 The expected variance of the portfolio.
77 Examples:
78 >>> import numpy as np
79 >>> fp = FrontierPoint(weights=np.array([0.5, 0.5]))
80 >>> covariance = np.array([[0.04, 0.0], [0.0, 0.09]])
81 >>> fp.variance(covariance)
82 0.0325
84 A ``CovarianceOperator`` backend gives the same number without
85 forming the matrix:
87 >>> from cvxcla import DenseCovariance
88 >>> fp.variance(DenseCovariance(covariance))
89 0.0325
91 """
92 return float(self.weights.T @ _covariance_matvec(covariance, self.weights))
95@dataclass(frozen=True)
96class TurningPoint(FrontierPoint):
97 """Turning point.
99 A turning point is a vector of weights, a lambda value, and a boolean vector
100 indicating which assets are free. All assets that are not free are blocked.
102 For problems with general inequality constraints ``G w <= h`` the turning
103 point also records which inequality *rows* are active (held at equality
104 ``g_i w = h_i``) via ``active_ineq``. These are the row analogue of the box
105 active set: a free weight on a bound is a per-variable active constraint, an
106 active inequality row is a per-row one. The default is an empty mask, so
107 box-and-equality problems (and the LASSO path) are unaffected.
109 Examples:
110 >>> import numpy as np
111 >>> tp = TurningPoint(
112 ... weights=np.array([0.5, 0.5, 0.0]),
113 ... free=np.array([True, True, False]),
114 ... lamb=1.25,
115 ... )
116 >>> tp.free_indices
117 array([0, 1])
118 >>> tp.blocked_indices
119 array([2])
121 ``free`` and ``blocked`` partition the assets, so the two index sets are
122 always complementary:
124 >>> len(tp.free_indices) + len(tp.blocked_indices) == len(tp.weights)
125 True
126 """
128 free: NDArray[np.bool_]
129 lamb: float = np.inf
130 active_ineq: NDArray[np.bool_] = field(default_factory=lambda: np.zeros(0, dtype=bool))
132 @property
133 def free_indices(self) -> np.ndarray:
134 """Returns the indices of the free assets."""
135 return np.where(self.free)[0]
137 @property
138 def blocked_indices(self) -> np.ndarray:
139 """Returns the indices of the blocked assets."""
140 return np.where(~self.free)[0]
143@dataclass(frozen=True)
144class Frontier:
145 """A frontier is a list of frontier points. Some of them might be turning points.
147 The frontier owns the ``mean`` and ``covariance`` the points are scored
148 against, so the risk/return properties below (``returns``, ``variance``,
149 ``volatility``, ``sharpe_ratio``) are vectors with one entry per point, in
150 frontier order -- from the maximum-return end towards minimum variance.
152 Examples:
153 >>> import numpy as np
154 >>> mean = np.array([0.1, 0.2])
155 >>> covariance = np.array([[0.04, 0.0], [0.0, 0.09]])
156 >>> frontier = Frontier(
157 ... mean=mean,
158 ... covariance=covariance,
159 ... frontier=[
160 ... FrontierPoint(weights=np.array([0.0, 1.0])),
161 ... FrontierPoint(weights=np.array([1.0, 0.0])),
162 ... ],
163 ... )
164 >>> len(frontier)
165 2
166 >>> frontier.returns
167 array([0.2, 0.1])
168 >>> frontier.volatility
169 array([0.3, 0.2])
171 Usually you get one from a solved :class:`~cvxcla.cla.CLA` rather than
172 by hand:
174 >>> from cvxcla import CLA
175 >>> solved = CLA.problem(mean, covariance).long_only().budget().trace()
176 >>> isinstance(solved.frontier, Frontier)
177 True
178 """
180 mean: NDArray[np.float64]
181 covariance: NDArray[np.float64] | CovarianceOperator
182 frontier: list[FrontierPoint] = field(default_factory=list)
184 def interpolate(self, num: int = 100) -> Frontier:
185 """Interpolate the frontier with additional points between existing points.
187 This method creates a new Frontier object with additional points interpolated
188 between the existing points. This is useful for creating a smoother representation
189 of the efficient frontier for visualization or analysis.
191 Args:
192 num: The number of points to use in the interpolation. The method will create
193 num-1 new points between each pair of adjacent existing points.
195 Returns:
196 A new Frontier object with the interpolated points.
198 Examples:
199 >>> import numpy as np
200 >>> frontier = Frontier(
201 ... mean=np.array([0.1, 0.2]),
202 ... covariance=np.array([[0.04, 0.0], [0.0, 0.09]]),
203 ... frontier=[
204 ... FrontierPoint(weights=np.array([0.0, 1.0])),
205 ... FrontierPoint(weights=np.array([1.0, 0.0])),
206 ... ],
207 ... )
208 >>> len(frontier.interpolate(num=10))
209 9
211 One adjacent pair yields ``num - 1`` points, and every interpolated
212 point is a convex combination of its neighbours -- so a budget
213 constraint satisfied at the turning points still holds between them:
215 >>> dense = frontier.interpolate(num=10)
216 >>> bool(np.allclose(dense.weights.sum(axis=1), 1.0))
217 True
219 """
221 def _interpolate() -> Iterator[FrontierPoint]:
222 """Yield interpolated frontier points between each adjacent pair."""
223 for w_right, w_left in zip(self.weights[0:-1], self.weights[1:], strict=False): # pragma: no mutate
224 for lamb in np.linspace(0, 1, num):
225 if lamb > 0:
226 yield FrontierPoint(weights=lamb * w_left + (1 - lamb) * w_right)
228 points = list(_interpolate())
229 return Frontier(frontier=points, mean=self.mean, covariance=self.covariance)
231 def __iter__(self) -> Iterator[FrontierPoint]:
232 """Iterate over all frontier points."""
233 yield from self.frontier
235 def __len__(self) -> int:
236 """Give number of frontier points."""
237 return len(self.frontier)
239 @property
240 def weights(self) -> np.ndarray:
241 """Matrix of weights. One row per point."""
242 return np.array([point.weights for point in self])
244 @property
245 def returns(self) -> np.ndarray:
246 """Vector of expected returns."""
247 return np.array([point.mean(self.mean) for point in self])
249 @property
250 def variance(self) -> np.ndarray:
251 """Vector of expected variances."""
252 return np.array([point.variance(self.covariance) for point in self])
254 @property
255 def sharpe_ratio(self) -> np.ndarray:
256 """Vector of expected Sharpe ratios."""
257 ratios: np.ndarray = self.returns / self.volatility
258 return ratios
260 @property
261 def volatility(self) -> np.ndarray:
262 """Vector of expected volatilities."""
263 vol: np.ndarray = np.sqrt(self.variance)
264 return vol
266 @property
267 def max_sharpe(self) -> tuple[float, np.ndarray]:
268 """Maximal Sharpe ratio on the frontier.
270 The maximiser lies on one of the two affine segments adjacent to the
271 turning point of largest discrete Sharpe ratio. On each segment the Sharpe
272 ratio has a closed-form maximiser (see :meth:`_segment_max_sharpe`), so the
273 result is exact rather than the product of a numerical line search.
275 Returns:
276 Tuple of maximal Sharpe ratio and the weights to achieve it
278 Examples:
279 >>> import numpy as np
280 >>> frontier = Frontier(
281 ... mean=np.array([0.1, 0.2]),
282 ... covariance=np.array([[0.04, 0.0], [0.0, 0.09]]),
283 ... frontier=[
284 ... FrontierPoint(weights=np.array([0.0, 1.0])),
285 ... FrontierPoint(weights=np.array([1.0, 0.0])),
286 ... ],
287 ... )
288 >>> sharpe, weights = frontier.max_sharpe
289 >>> float(np.round(sharpe, 6))
290 0.833333
291 >>> bool(np.allclose(weights, [0.529412, 0.470588], atol=1e-6))
292 True
294 The continuous optimum sits *between* two turning points, so it
295 beats every discrete point on the frontier -- which is the reason
296 this is not simply ``argmax(sharpe_ratio)``:
298 >>> bool(sharpe > frontier.sharpe_ratio.max())
299 True
301 """
302 weights = self.weights
303 sharpe_ratios = self.sharpe_ratio
305 # The discrete maximum brackets the continuous one: the optimum sits on a
306 # segment touching the turning point of largest Sharpe ratio.
307 sr_position_max = int(np.argmax(sharpe_ratios))
308 right = min(sr_position_max + 1, len(self) - 1)
309 left = max(0, sr_position_max - 1)
311 # Look to the left and to the right of the discrete maximum.
312 if right > sr_position_max:
313 sharpe_ratio_right, w_right = self._segment_max_sharpe(weights[sr_position_max], weights[right])
314 else:
315 w_right = weights[sr_position_max]
316 sharpe_ratio_right = sharpe_ratios[sr_position_max]
318 if left < sr_position_max:
319 sharpe_ratio_left, w_left = self._segment_max_sharpe(weights[left], weights[sr_position_max])
320 else:
321 w_left = weights[sr_position_max]
322 sharpe_ratio_left = sharpe_ratios[sr_position_max]
324 if sharpe_ratio_left > sharpe_ratio_right:
325 return sharpe_ratio_left, w_left
327 return sharpe_ratio_right, w_right
329 def _segment_max_sharpe(self, w0: np.ndarray, w1: np.ndarray) -> tuple[float, np.ndarray]:
330 """Closed-form maximum Sharpe ratio on the affine segment between two points.
332 Parametrise the segment as ``w(t) = (1 - t) w0 + t w1`` for ``t`` in
333 ``[0, 1]``. The expected return is affine and the variance quadratic in
334 ``t``, so the Sharpe ratio is::
336 S(t) = (a0 + a1 t) / sqrt(c0 + c1 t + c2 t**2)
338 Its derivative has a *linear* numerator (the ``t**2`` terms cancel), so
339 there is a single stationary point
340 ``t* = (a0 c1 - 2 a1 c0) / (a1 c1 - 2 a0 c2)``. The maximiser over the
341 segment is therefore whichever of ``{0, 1, clamp(t*)}`` yields the largest
342 Sharpe ratio, evaluated in closed form rather than by a bounded line search.
344 Args:
345 w0: Weights at the ``t = 0`` end of the segment.
346 w1: Weights at the ``t = 1`` end of the segment.
348 Returns:
349 Tuple of the maximal Sharpe ratio on the segment and its weights.
351 """
352 delta = w1 - w0
353 sigma_w0 = _covariance_matvec(self.covariance, w0)
354 sigma_delta = _covariance_matvec(self.covariance, delta)
355 a0 = float(self.mean @ w0)
356 a1 = float(self.mean @ delta)
357 c0 = float(w0 @ sigma_w0)
358 c1 = 2.0 * float(w0 @ sigma_delta)
359 c2 = float(delta @ sigma_delta)
361 def sharpe_at(t: float) -> tuple[float, np.ndarray]:
362 """Sharpe ratio and weights at position ``t`` along the segment."""
363 weight = w0 + t * delta
364 sharpe = (a0 + a1 * t) / np.sqrt(c0 + c1 * t + c2 * t * t)
365 return float(sharpe), weight
367 # Candidate positions: the two endpoints and the interior stationary point
368 # (only when it falls strictly inside the segment).
369 candidates = [0.0, 1.0]
370 denominator = a1 * c1 - 2.0 * a0 * c2
371 if denominator != 0.0:
372 t_star = (a0 * c1 - 2.0 * a1 * c0) / denominator
373 if 0.0 < t_star < 1.0:
374 candidates.append(t_star)
376 return max((sharpe_at(t) for t in candidates), key=lambda item: item[0])
378 def plot(self, volatility: bool = False, markers: bool = True) -> go.Figure:
379 """Plot the efficient frontier.
381 This function creates a line plot of the efficient frontier, with expected return
382 on the y-axis and either variance or volatility on the x-axis.
384 Args:
385 volatility: If True, plot volatility (standard deviation) on the x-axis.
386 If False, plot variance on the x-axis.
387 markers: If True, show markers at each point on the frontier.
389 Returns:
390 A plotly Figure object that can be displayed or saved.
392 """
393 try:
394 import plotly.graph_objects as go
395 except ImportError as e:
396 msg = "Plotting requires plotly. Install it with: pip install cvxcla[plot]"
397 raise ImportError(msg) from e
399 fig = go.Figure()
401 x = self.volatility if volatility else self.variance
402 axis_title = "Expected volatility" if volatility else "Expected variance"
404 fig.add_trace(
405 go.Scatter(x=x, y=self.returns, mode="lines+markers" if markers else "lines", name="Efficient Frontier")
406 )
408 fig.update_layout(
409 xaxis_title=axis_title,
410 yaxis_title="Expected Return",
411 )
413 return fig