Coverage for book/marimo/notebooks/Experiment5.py: 100%
100 statements
« prev ^ index » next coverage.py v7.16.0, created at 2026-09-09 08:57 +0000
1# /// script
2# requires-python = ">=3.12"
3# dependencies = [
4# "marimo==0.24.0",
5# "numpy==2.4.6",
6# "plotly==6.9.0",
7# "polars==1.44.1",
8# "jquantstats==0.11.0",
9# "tinycta==0.14.0"
10# ]
11#
12# [tool.ty.environment]
13# # ``from preamble import ...`` resolves at runtime via the sys.path.insert in the
14# # setup cell below. ty analyses a PEP 723 script in isolation from the project, so
15# # pyproject.toml's [tool.ty.environment] never reaches this file and the path has
16# # to be declared here. Preserve this table if marimo rewrites the header.
17# extra-paths = ["."]
18# ///
20"""Experiment 5: Advanced CTA strategy with correlation-based optimization.
22This module implements a sophisticated trend-following strategy that
23incorporates dynamic conditional correlation (DCC) and matrix optimization
24techniques to improve portfolio construction and risk management.
25"""
27import marimo
29__generated_with = "0.23.9"
30app = marimo.App()
32with app.setup:
33 import sys
34 from pathlib import Path
36 import marimo as mo
37 import numpy as np
38 import polars as pl
39 from jquantstats import Portfolio
40 from tinycta.linalg import inv_a_norm, solve
41 from tinycta.osc import osc
42 from tinycta.signal import shrink2id
43 from tinycta.util import vol_adj
45 sys.path.insert(0, str(Path(__file__).parent))
47 from preamble import date_col, load_prices
49 prices = load_prices(__file__)
50 prices_only = prices.drop(date_col)
51 assets = prices_only.columns
54@app.cell(hide_code=True)
55def _():
56 mo.md(r"""
57 # CTA 5.0 - Optimization 2.0
58 """)
59 return
62@app.function
63def f(price: "pl.Expr", fast: int = 32, slow: int = 96, vola: int = 32, clip: float = 4.2) -> "pl.Expr":
64 """Return the tanh oscillator of vol-adjusted cumulative price."""
65 return osc(vol_adj(price, vola=vola, clip=clip, min_samples=300).cum_sum(), fast=fast, slow=slow).tanh()
68@app.cell
69def _():
70 fast = mo.ui.slider(4, 192, step=4, value=32, label="Fast Moving Average")
71 slow = mo.ui.slider(4, 192, step=4, value=96, label="Slow Moving Average")
72 vola = mo.ui.slider(4, 192, step=4, value=32, label="Volatility")
73 winsor = mo.ui.slider(1.0, 6.0, step=0.1, value=4.2, label="Winsorizing")
74 corr = mo.ui.slider(50, 500, step=10, value=200, label="Correlation")
75 shrinkage = mo.ui.slider(0.0, 1.0, step=0.05, value=0.5, label="Shrinkage")
77 mo.vstack([fast, slow, vola, winsor, corr, shrinkage])
78 return corr, shrinkage, vola, winsor
81@app.function
82def ewm_covariance(returns_adj: "pl.DataFrame", *, corr: int) -> "np.ndarray":
83 """EWM covariance tensor of the vol-adjusted returns (the Engle-DCC numerator).
85 ``cov_t(i, j) = ewm_t(r_i * r_j) - ewm_t(r_i) * ewm_t(r_j)``, evaluated for
86 every day; the result has shape ``(n_rows, n_assets, n_assets)`` and is
87 symmetric in ``(i, j)``.
88 """
89 columns = returns_adj.columns
90 n_assets = len(columns)
91 n_rows = len(returns_adj)
93 ewm_means_np = returns_adj.select(pl.all().ewm_mean(com=corr, min_samples=int(corr))).to_numpy()
95 pair_indices = [(i, j) for i in range(n_assets) for j in range(i, n_assets)]
96 ewm_prod_np = returns_adj.select(
97 [
98 (pl.col(columns[i]) * pl.col(columns[j]))
99 .fill_nan(None)
100 .ewm_mean(com=corr, min_samples=int(corr))
101 .alias(f"p{i}_{j}")
102 for i, j in pair_indices
103 ]
104 ).to_numpy()
106 cov_np = np.full((n_rows, n_assets, n_assets), np.nan)
107 for _k, (_i, _j) in enumerate(pair_indices):
108 _cov = ewm_prod_np[:, _k] - ewm_means_np[:, _i] * ewm_means_np[:, _j]
109 cov_np[:, _i, _j] = _cov
110 cov_np[:, _j, _i] = _cov
111 return cov_np
114@app.function
115def correlation_from_covariance(cov_np: "np.ndarray") -> "np.ndarray":
116 """Normalize a per-day covariance tensor to a correlation tensor.
118 Days with non-positive variance keep their NaN off-diagonals; the diagonal
119 is forced to 1 wherever the variance is positive.
120 """
121 n_assets = cov_np.shape[1]
122 _var = cov_np[:, np.arange(n_assets), np.arange(n_assets)]
123 with np.errstate(invalid="ignore", divide="ignore"):
124 _denom = np.sqrt(_var[:, :, None] * _var[:, None, :])
125 cor_3d: np.ndarray = cov_np / _denom
126 for _k in range(n_assets):
127 cor_3d[_var[:, _k] > 0, _k, _k] = 1.0
128 return cor_3d
131@app.function
132def dcc_correlation(prices_only: "pl.DataFrame", *, vola: int, clip: float, corr: int) -> "np.ndarray":
133 """Engle-DCC per-day correlation tensor, shape ``(n_rows, n_assets, n_assets)``."""
134 returns_adj = prices_only.select(vol_adj(pl.all(), vola=vola, clip=clip, min_samples=300))
135 cov_np = ewm_covariance(returns_adj, corr=corr)
136 return correlation_from_covariance(cov_np)
139@app.function
140def positions(
141 cor_3d: "np.ndarray", mu: "np.ndarray", vo: "np.ndarray", prices_np: "np.ndarray", *, shrinkage: float
142) -> "np.ndarray":
143 """Per-day risk-parity positions from the shrunk DCC correlation tensor.
145 Each day, the correlation matrix is shrunk towards the identity, restricted
146 to the assets with a finite price, and used to solve for the risk-scaled
147 position. Days with no live assets or a singular/zero norm are left flat.
148 """
149 n_rows, n_assets = prices_np.shape
150 pos_matrix = np.zeros((n_rows, n_assets))
151 for _n in range(n_rows):
152 _mask = np.isfinite(prices_np[_n])
153 if _mask.sum() == 0:
154 continue
155 _full_shrunk = shrink2id(cor_3d[_n], lamb=shrinkage)
156 _matrix = _full_shrunk[_mask, :][:, _mask]
157 _expected_mu = np.nan_to_num(mu[_n][_mask])
158 _expected_vo = np.nan_to_num(vo[_n][_mask])
159 _norm = inv_a_norm(_expected_mu, _matrix)
160 if _norm == 0 or np.isnan(_norm):
161 continue
162 _risk_pos = solve(_matrix, _expected_mu) / _norm
163 pos_matrix[_n, _mask] = np.nan_to_num(1e6 * _risk_pos / _expected_vo, nan=0.0)
164 return pos_matrix
167@app.cell
168def _(corr, shrinkage, vola, winsor):
169 # EWM correlation (DCC by Engle), then per-day risk-parity positions.
170 cor_3d = dcc_correlation(prices_only, vola=vola.value, clip=winsor.value, corr=corr.value)
172 mu = prices_only.select(f(pl.all(), fast=32, slow=96, vola=vola.value, clip=winsor.value)).to_numpy()
173 vo = prices_only.select(
174 pl.all().fill_nan(None).pct_change().ewm_std(com=vola.value, min_samples=int(vola.value))
175 ).to_numpy()
177 pos_matrix = positions(cor_3d, mu, vo, prices_only.to_numpy(), shrinkage=shrinkage.value)
179 portfolio = Portfolio.from_cash_position(
180 prices=prices,
181 cash_position=pl.concat(
182 [prices.select(date_col), pl.from_numpy(pos_matrix, schema=dict.fromkeys(assets, pl.Float64))],
183 how="horizontal_extend",
184 ),
185 aum=1e8,
186 )
187 return (portfolio,)
190@app.cell
191def _(portfolio):
192 print(portfolio.stats.sharpe())
193 return
196@app.cell(hide_code=True)
197def _():
198 mo.md(r"""
199 # Conclusions
200 * Dramatic relative improvements observable despite using the same signals as in previous Experiment.
201 * Main difference here is to take advantage of cross-correlations in the risk measurement.
202 * Possible to add constraints on individual assets or groups of them.
203 * Possible to reflect trading costs in objective with regularization terms (Ridge, Lars, Elastic Nets, ...)
204 """)
205 return
208@app.cell
209def _(portfolio):
210 portfolio.plots.snapshot()
211 return
214@app.cell
215def _(portfolio):
216 portfolio.stats.summary()
217 return
220@app.cell
221def _():
222 return
225if __name__ == "__main__":
226 app.run()