Coverage for src/pyhrp/algos.py: 100%
69 statements
« prev ^ index » next coverage.py v7.14.3, created at 2026-11-03 05:36 +0000
1"""Portfolio optimization algorithms for hierarchical risk parity.
3This module implements various portfolio optimization algorithms:
4- risk_parity: The main hierarchical risk parity algorithm
5- schur_risk_parity: Schur Complementary Allocation (Cotton, arXiv:2411.05807)
6- one_over_n: A simple equal-weight allocation strategy
8Allocator contract
9------------------
10All three allocators take the same inputs — a ``Cluster`` tree (``root``) plus
11the asset names — and none of them changes the *shape* of that tree: no node is
12added, removed or re-parented. They differ in where the weights land.
14``risk_parity`` and ``schur_risk_parity`` share the recursive ``_allocate_with``
15scaffolding and write **into the tree they are given**: every node's
16``portfolio`` is replaced rather than accumulated into, and the same root object
17is returned. Replacing makes them idempotent — re-running on an already weighted
18tree, with a different covariance matrix or gamma, gives the same answer as
19running on a fresh one — but the caller's tree *is* modified. Note in particular
20that ``Dendrogram`` is a frozen dataclass while the ``Cluster`` it holds is not:
21passing ``dendrogram.root`` to either allocator rewrites the portfolios inside
22that dendrogram. Pass ``copy.deepcopy(dendrogram.root)`` to keep the original
23tree unweighted.
25``one_over_n`` differs on both counts. It accumulates into a local buffer, so it
26leaves the input tree untouched entirely, and its *output* is a generator
27yielding the equal-weight portfolio one tree level at a time (see its
28docstring), because its purpose is to expose the allocation as it deepens rather
29than a single final result.
30"""
32from __future__ import annotations
34from collections.abc import Callable, Generator, Sequence
35from copy import deepcopy
37import numpy as np
38import polars as pl
40from .cluster import Cluster, Portfolio, _portfolio
41from .operators import CovarianceOperator, as_operator, solve_block
43__all__ = ["one_over_n", "risk_parity", "schur_risk_parity"]
46def risk_parity(root: Cluster, cov: pl.DataFrame | CovarianceOperator, assets: Sequence[str] | None = None) -> Cluster:
47 """Compute hierarchical risk parity weights for a cluster tree.
49 This is the main algorithm for hierarchical risk parity. It recursively
50 traverses the cluster tree and assigns weights to each node based on
51 the risk parity principle.
53 Note:
54 The tree is modified in place: the portfolio of every node is rebuilt
55 from scratch and ``root`` itself is returned, so the function is
56 idempotent and a tree can be reused with a different covariance matrix.
57 Deep-copy the tree first if you need to keep it unweighted.
59 The covariance is reached only through block products, one per child of
60 every node, so a factor or Gram :class:`~pyhrp.operators.CovarianceOperator`
61 allocates without forming the ``n x n`` matrix.
63 Args:
64 root (Cluster): The root node of the cluster tree
65 cov (pl.DataFrame | CovarianceOperator): Covariance matrix of asset returns,
66 as a DataFrame or as an operator
67 assets (Sequence[str], optional): Asset names for an operator ``cov``;
68 defaults to ``"0", "1", ...``. A DataFrame names its own assets.
70 Returns:
71 Cluster: The same root node, with portfolio weights assigned
73 Examples:
74 >>> import polars as pl
75 >>> from pyhrp.cluster import Cluster
76 >>> from pyhrp.algos import risk_parity
77 >>> cov = pl.DataFrame({"A": [4.0, 0.0], "B": [0.0, 1.0]})
78 >>> root = Cluster(2, left=Cluster(0), right=Cluster(1))
79 >>> cluster = risk_parity(root=root, cov=cov)
80 >>> round(cluster.portfolio["B"], 1)
81 0.8
82 """
84 def node_variances(left: Block, right: Block, op: CovarianceOperator) -> tuple[float, float]:
85 """Plain block variance of each child sub-portfolio."""
86 return _block_variance(left, op), _block_variance(right, op)
88 return _allocate_with(root, cov, node_variances, assets)
91def schur_risk_parity(
92 root: Cluster,
93 cov: pl.DataFrame | CovarianceOperator,
94 gamma: float = 0.5,
95 assets: Sequence[str] | None = None,
96) -> Cluster:
97 """Compute Schur Complementary Allocation weights for a cluster tree.
99 An extension of HRP introduced by Peter Cotton (arXiv:2411.05807) that augments
100 sub-covariance matrices with off-diagonal block information via Schur complements.
101 At gamma=0 this recovers standard HRP; at gamma=1 it recovers the minimum-variance
102 portfolio through the same recursive structure.
104 Note:
105 The tree is modified in place: the portfolio of every node is rebuilt
106 from scratch and ``root`` itself is returned, so the function is
107 idempotent and a tree can be reused with a different covariance matrix
108 or gamma. Deep-copy the tree first if you need to keep it unweighted.
110 The Schur terms need, per node, the two cross-block products
111 ``Sigma[L, R] w_R`` and ``Sigma[R, L] w_L`` and one solve against each
112 child's principal block, which a factor operator does by Woodbury in
113 ``O(|C| k**2)`` rather than ``O(|C|**3)``. The block quadratic forms
114 ``B D^{-1} B^T`` are never formed: only their action on the child weights is.
116 Args:
117 root (Cluster): The root node of the cluster tree
118 cov (pl.DataFrame | CovarianceOperator): Covariance matrix of asset returns,
119 as a DataFrame or as an operator
120 gamma (float): Interpolation parameter in [0, 1]. 0 = HRP, 1 = minimum variance.
121 assets (Sequence[str], optional): Asset names for an operator ``cov``;
122 defaults to ``"0", "1", ...``. A DataFrame names its own assets.
124 Returns:
125 Cluster: The same root node, with portfolio weights assigned
127 Raises:
128 ValueError: If gamma is outside the interval [0, 1].
130 Examples:
131 >>> import polars as pl
132 >>> from pyhrp.cluster import Cluster
133 >>> from pyhrp.algos import schur_risk_parity
134 >>> cov = pl.DataFrame({"A": [4.0, 0.0], "B": [0.0, 1.0]})
135 >>> root = Cluster(2, left=Cluster(0), right=Cluster(1))
136 >>> cluster = schur_risk_parity(root=root, cov=cov, gamma=0.5)
137 >>> round(cluster.portfolio["B"], 1)
138 0.8
139 """
140 if not 0.0 <= gamma <= 1.0:
141 msg = f"gamma must be in [0, 1], got {gamma}"
142 raise ValueError(msg)
144 def node_variances(left: Block, right: Block, op: CovarianceOperator) -> tuple[float, float]:
145 """Schur-augmented block variance of each child, conditioned on the other."""
146 (li, w_left), (ri, w_right) = left, right
148 v_left = _block_variance(left, op)
149 v_right = _block_variance(right, op)
150 if gamma == 0.0:
151 return v_left, v_right
153 # Schur-augmented blocks A - gamma B D^{-1} B^T and D - gamma B^T A^{-1} B
154 # (A = Sigma_LL, B = Sigma_LR, D = Sigma_RR), applied to the child weights
155 # through one cross product and one block solve each.
156 bt_w = op.block_matvec(ri, li, w_left)
157 b_w = op.block_matvec(li, ri, w_right)
158 return (
159 _schur_variance(v_left, gamma * float(bt_w @ solve_block(op, ri, bt_w))),
160 _schur_variance(v_right, gamma * float(b_w @ solve_block(op, li, b_w))),
161 )
163 return _allocate_with(root, cov, node_variances, assets)
166# A child sub-portfolio as (positions into the covariance, weights), in leaf order.
167Block = tuple[np.ndarray, np.ndarray]
169# Given a node's two children and the covariance operator, return the
170# (v_left, v_right) risk pair used to split the node.
171NodeVariances = Callable[[Block, Block, CovarianceOperator], tuple[float, float]]
174def _allocate_with(
175 root: Cluster,
176 cov: pl.DataFrame | CovarianceOperator,
177 node_variances: NodeVariances,
178 assets: Sequence[str] | None = None,
179) -> Cluster:
180 """Shared scaffolding for the recursive risk-based allocators.
182 Wraps the covariance as an operator, then walks the tree bottom-up, splitting
183 each node's weight between its children inversely to the ``(v_left, v_right)``
184 pair supplied by ``node_variances``. The only thing that distinguishes
185 ``risk_parity`` from ``schur_risk_parity`` is that per-node variance rule;
186 everything else lives here.
188 The walk holds one weight per asset, not one portfolio per node. Every node's
189 subtree covers a contiguous slice of the leaf order, so after its children are
190 split that slice holds the node's own portfolio; scaling it in place by the
191 node's split turns it into the parent's share. Memory is ``O(n)`` for any tree
192 shape. Each node keeps only its split, from which its portfolio is rebuilt
193 when read (see :class:`~pyhrp.cluster.Cluster`).
195 Args:
196 root (Cluster): The root node of the cluster tree.
197 cov (pl.DataFrame | CovarianceOperator): Covariance matrix of asset returns.
198 node_variances (NodeVariances): Per-node rule mapping a node's left/right
199 child sub-portfolios (and the covariance operator) to the
200 ``(v_left, v_right)`` risk pair used to split that node.
201 assets (Sequence[str], optional): Asset names for an operator ``cov``.
203 Returns:
204 Cluster: The root node with portfolio weights assigned.
206 Raises:
207 ValueError: If the tree's leaves do not index the covariance columns one-to-one.
208 """
209 op, names = as_operator(cov, assets)
211 # Iterative layout rather than recursion: a chain-degenerate tree is as deep as
212 # the universe is wide. The layout validates every non-leaf node, so a malformed
213 # tree raises before any weight is written.
214 leaves, spans = root._layout()
216 # A leaf's value indexes into the assets, so a tree built for a different universe
217 # would either silently drop assets or fail with a bare IndexError deep in the walk.
218 leaf_values = sorted(int(leaf.value) for leaf in leaves)
219 if leaf_values != list(range(op.n)):
220 msg = (
221 f"Cluster tree does not match the covariance matrix: expected {op.n} leaves "
222 f"indexing columns 0..{op.n - 1}, got {len(leaf_values)} leaves with values {leaf_values}"
223 )
224 raise ValueError(msg)
226 position = np.array([int(leaf.value) for leaf in leaves], dtype=np.intp)
227 w = np.ones(len(leaves))
228 shares: list[float] = []
229 for _, lo, mid, hi in spans:
230 v_left, v_right = node_variances((position[lo:mid], w[lo:mid]), (position[mid:hi], w[mid:hi]), op)
231 share = _left_share(v_left, v_right)
232 w[lo:mid] *= share
233 w[mid:hi] *= 1.0 - share
234 shares.append(share)
236 # Every node's portfolio is replaced, never accumulated into, which keeps
237 # repeated allocations on the same tree idempotent.
238 for leaf in leaves:
239 leaf._set_allocation(1.0, names)
240 for (node, _, _, _), share in zip(spans, shares, strict=True):
241 node._set_allocation(share, names)
242 # The root's portfolio is the weight vector just computed; keep it rather than rebuild it.
243 root.portfolio = _portfolio(leaves, w, names)
244 return root
247def _block_variance(block: Block, op: CovarianceOperator) -> float:
248 """Compute the variance of a sub-portfolio from one block product with the covariance operator."""
249 idx, w = block
250 return float(w @ op.block_matvec(idx, idx, w))
253# A Schur-augmented variance within this fraction of the plain block variance is
254# round-off of an exact zero, and is treated as zero.
255SCHUR_RTOL = 1e-10
258def _schur_variance(variance: float, correction: float) -> float:
259 """The Schur-augmented variance ``variance - correction``, with cancellation set to zero.
261 The augmented block ``A - gamma B D^+ B^T`` is positive semidefinite, so the
262 true value is non-negative. It is exactly zero when ``gamma = 1`` and the
263 block adds no rank to its sibling (a Gram covariance with fewer observations
264 than assets), and then round-off of either sign is all that is left; a
265 negative value would push the split outside ``[0, 1]`` and the weights below
266 zero. Values within ``SCHUR_RTOL`` of the plain variance are therefore zero,
267 and two zero children split evenly (see :func:`_left_share`).
269 Args:
270 variance (float): The plain block variance ``w' A w``.
271 correction (float): The Schur term ``gamma w' B D^+ B^T w``.
273 Returns:
274 float: The augmented variance, ``0.0`` if it is within round-off of zero.
275 """
276 augmented = variance - correction
277 return augmented if augmented > SCHUR_RTOL * variance else 0.0
280def _left_share(v_left: float, v_right: float) -> float:
281 """The share of a node's weight that goes to its left child, inversely proportional to risk.
283 The split satisfies v_left * alpha_left == v_right * alpha_right with
284 alpha_left + alpha_right == 1. If both variances are zero (e.g. riskless
285 sub-portfolios), the weight is split equally.
287 Args:
288 v_left (float): Variance of the left sub-portfolio
289 v_right (float): Variance of the right sub-portfolio
291 Returns:
292 float: ``alpha_left``, in [0, 1]
293 """
294 total = v_left + v_right
295 return v_right / total if total > 0 else 0.5
298def one_over_n(root: Cluster, assets: list[str]) -> Generator[tuple[int, Portfolio]]:
299 """Generate 1/N (equal-weight) portfolios one tree level at a time.
301 This implements a hierarchical 1/N strategy where weights are distributed
302 equally among the leaves of each cluster, and the weight budget halves at
303 each successive level of the tree.
305 Unlike :func:`risk_parity` and :func:`schur_risk_parity` — which rebuild a
306 single final allocation and return the root ``Cluster`` — this allocator is
307 intentionally a **generator**: its purpose is to expose the equal-weight
308 allocation as the tree deepens, yielding one portfolio per level. It shares
309 the sibling input contract (a ``Cluster`` tree plus the asset names) and, like
310 them, does not mutate the tree: weights accumulate in a local buffer, so a
311 leaf that terminates at a shallow level keeps its weight in the deeper levels
312 (each yielded portfolio is therefore a complete allocation over all assets),
313 and re-running on the same tree yields an identical sequence.
315 Args:
316 root (Cluster): The root node of the cluster tree.
317 assets (list[str]): Asset names; a leaf's value indexes into this list.
319 Yields:
320 tuple[int, Portfolio]: The level number and the (cumulative) equal-weight
321 portfolio at that level.
323 Examples:
324 >>> import polars as pl
325 >>> from pyhrp.hrp import build_tree
326 >>> from pyhrp.algos import one_over_n
327 >>> cor = pl.DataFrame({"A": [1.0, 0.3], "B": [0.3, 1.0]})
328 >>> dg = build_tree(cor, method="ward")
329 >>> levels = list(one_over_n(dg.root, dg.assets))
330 >>> len(levels) > 0
331 True
332 """
333 # Accumulate into a local buffer so the input tree is never mutated.
334 portfolio = Portfolio()
336 # Initial weight to distribute
337 w: float = 1.0
339 # Process each level of the tree
340 for n, level in enumerate(root.levels):
341 for node in level:
342 # Distribute weight equally among all leaves in this node
343 for leaf in node.leaves:
344 portfolio[assets[leaf.value]] = w / node.leaf_count
346 # Reduce weight for the next level
347 w *= 0.5
349 # Yield the current level number and a deep copy of the portfolio
350 yield n, deepcopy(portfolio)