From 8a61b0fa778e771b357472e1a8c63719d464249e Mon Sep 17 00:00:00 2001 From: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> Date: Thu, 22 Oct 2026 21:53:37 +0500 Subject: [PATCH 1/3] Handle missing return observations explicitly in shrinkage and semicovariance Signed-off-by: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> --- docs/RiskModels.rst | 14 +++++ pypfopt/risk_models.py | 36 +++++++++--- tests/test_black_litterman.py | 40 ++++++------- tests/test_cla.py | 2 +- tests/test_efficient_frontier.py | 8 +-- tests/test_missing_return_estimation.py | 75 +++++++++++++++++++++++++ 6 files changed, 143 insertions(+), 32 deletions(-) create mode 100644 tests/test_missing_return_estimation.py diff --git a/docs/RiskModels.rst b/docs/RiskModels.rst index 3fc8b152..5c56bf88 100644 --- a/docs/RiskModels.rst +++ b/docs/RiskModels.rst @@ -163,3 +163,17 @@ References .. [3] Ledoit, O., & Wolf, M. (2003). `Honey, I Shrunk the Sample Covariance Matrix `_ The Journal of Portfolio Management, 30(4), 110–119. https://doi.org/10.3905/jpm.2004.110 .. [4] Ledoit, O., & Wolf, M. (2001). `Improved estimation of the covariance matrix of stock returns with an application to portfolio selection `_, 10, 603–621. .. [5] Chen et al. (2010), `Shrinkage Algorithms for MMSE Covariance Estimation `_, IEEE Transactions on Signals Processing, 58(10), 5016-5029. + +Missing returns +--------------- + +Automatic Ledoit-Wolf and OAS estimates use the common window of fully observed +returns. They warn when incomplete rows are discarded and require at least two +complete observations. Missing values are not zero returns; explicitly choose +an imputation model before calling these estimators if dropping rows would lose +too much information. Manual shrinkage retains the pairwise sample covariance. + +``semicovariance`` averages each pair over jointly observed returns, including +observed zero returns. It raises an error if any pair has no joint observation. +As with pairwise sample covariance, the resulting matrix may need the existing +positive-semidefinite repair, which can alter individual entries. diff --git a/pypfopt/risk_models.py b/pypfopt/risk_models.py index 171ce0d0..6db87ae9 100644 --- a/pypfopt/risk_models.py +++ b/pypfopt/risk_models.py @@ -246,10 +246,13 @@ def semicovariance( returns = prices else: returns = returns_from_prices(prices, log_returns) + observed = returns.notna().astype(float) + counts = observed.T @ observed + if (counts == 0).any().any(): + raise ValueError("Every asset pair must have at least one joint observation") drops = np.fmin(returns - benchmark, 0) - T = drops.shape[0] return fix_nonpositive_semidefinite( - (drops.T @ drops) / T * frequency, kwargs.get("fix_method", "spectral") + (drops.T @ drops) / counts * frequency, kwargs.get("fix_method", "spectral") ) @@ -462,6 +465,24 @@ def __init__(self, prices, returns_data=False, frequency=252, log_returns=False) self.S = self.X.cov().values self.delta = None # shrinkage constant + def _complete_returns(self): + """Select common observations for estimators requiring a complete sample.""" + complete = self.X.dropna(how="any") + if len(complete) < 2: + raise ValueError( + "Shrinkage estimation requires at least two complete observations" + ) + dropped = len(self.X) - len(complete) + if dropped: + warnings.warn( + "Dropping {} incomplete return observations for shrinkage estimation".format( + dropped + ), + UserWarning, + stacklevel=3, + ) + return complete.values + def _format_and_annualize(self, raw_cov_array): """ Helper method which annualises the output of shrinkage calculations, @@ -529,7 +550,7 @@ def ledoit_wolf(self, shrinkage_target="constant_variance"): shrunk sample covariance matrix """ if shrinkage_target == "constant_variance": - X = np.nan_to_num(self.X.values) + X = self._complete_returns() shrunk_cov, self.delta = self.covariance.ledoit_wolf(X) elif shrinkage_target == "single_factor": shrunk_cov, self.delta = self._ledoit_wolf_single_factor() @@ -553,7 +574,7 @@ def _ledoit_wolf_single_factor(self): np.ndarray, float shrunk sample covariance matrix, shrinkage constant """ - X = np.nan_to_num(self.X.values) + X = self._complete_returns() # De-mean returns t, n = np.shape(X) @@ -610,10 +631,11 @@ def _ledoit_wolf_constant_correlation(self): np.ndarray, float shrunk sample covariance matrix, shrinkage constant """ - X = np.nan_to_num(self.X.values) + X = self._complete_returns() t, n = np.shape(X) - S = self.S # sample cov matrix + # Use the same complete observations for covariance and higher moments. + S = np.atleast_2d(np.cov(X, rowvar=False)) # Constant correlation target var = np.diag(S).reshape(-1, 1) @@ -663,6 +685,6 @@ def oracle_approximating(self): np.ndarray shrunk sample covariance matrix """ - X = np.nan_to_num(self.X.values) + X = self._complete_returns() shrunk_cov, self.delta = self.covariance.oas(X) return self._format_and_annualize(shrunk_cov) diff --git a/tests/test_black_litterman.py b/tests/test_black_litterman.py index 25bdff39..5e241048 100644 --- a/tests/test_black_litterman.py +++ b/tests/test_black_litterman.py @@ -180,26 +180,26 @@ def test_bl_returns_all_views(): posterior_rets, np.array( [ - 0.11168648, - 0.16782938, - 0.12516799, - 0.24067997, - 0.32848296, - -0.22789895, - 0.16311297, - 0.11928542, - 0.25414308, - 0.11007738, - 0.06282615, - -0.03140218, - -0.16977172, - 0.05254821, - -0.10463884, - 0.32173375, - 0.26399864, - 0.1118594, - 0.22999558, - 0.08977448, + 0.0659959179376, + 0.142395737648, + 0.0897272713478, + 0.203540995698, + 0.329077173837, + -0.194678572889, + 0.140774034267, + 0.106320649238, + 0.206079787474, + 0.110853248026, + 0.0604675037668, + -0.0334756812694, + -0.127710177297, + 0.0364737295266, + -0.0939289294804, + 0.294362055441, + 0.239905888459, + 0.0852597727498, + 0.210724156209, + 0.061193023013, ] ), ) diff --git a/tests/test_cla.py b/tests/test_cla.py index c5f69e23..a96cff78 100644 --- a/tests/test_cla.py +++ b/tests/test_cla.py @@ -106,7 +106,7 @@ def test_cla_max_sharpe_semicovariance(): np.testing.assert_almost_equal(cla.weights.sum(), 1) np.testing.assert_allclose( cla.portfolio_performance(risk_free_rate=0.02), - (0.2721798377099145, 0.07258537193305141, 3.474251505420551), + (0.328899626820958, 0.1707193786308789, 1.809399901160873), atol=1e-4, rtol=1e-4, ) diff --git a/tests/test_efficient_frontier.py b/tests/test_efficient_frontier.py index 038788e0..119eba3d 100644 --- a/tests/test_efficient_frontier.py +++ b/tests/test_efficient_frontier.py @@ -1257,7 +1257,7 @@ def test_max_sharpe_semicovariance(): assert all([i >= 0 for i in w.values()]) np.testing.assert_allclose( ef.portfolio_performance(risk_free_rate=0.02), - (0.2762965426962885, 0.07372667096108301, 3.476307004714425), + (0.3337156666549726, 0.1733001108784157, 1.810245042919043), ) @@ -1271,7 +1271,7 @@ def test_max_sharpe_short_semicovariance(): np.testing.assert_almost_equal(ef.weights.sum(), 1) np.testing.assert_allclose( ef.portfolio_performance(risk_free_rate=0.02), - (0.42444834528495234, 0.0898263632679403, 4.50255727350929), + (0.6050331853361551, 0.1756415312028965, 3.330836285299403), ) @@ -1291,7 +1291,7 @@ def test_min_volatility_shrunk_L2_reg(): assert all([i >= 0 for i in w.values()]) np.testing.assert_allclose( ef.portfolio_performance(risk_free_rate=0.02), - (0.17358178582309983, 0.19563960638632416, 0.7850239972361532), + (0.1811698577736233, 0.1511922398792176, 1.065992923329772), ) @@ -1308,7 +1308,7 @@ def test_efficient_return_shrunk(): assert all([i >= 0 for i in w.values()]) np.testing.assert_allclose( ef.portfolio_performance(risk_free_rate=0.02), - (0.22, 0.08892163937693882, 2.2491713085967766), + (0.22, 0.1263908584841322, 1.582392923022269), ) diff --git a/tests/test_missing_return_estimation.py b/tests/test_missing_return_estimation.py new file mode 100644 index 00000000..748b2d90 --- /dev/null +++ b/tests/test_missing_return_estimation.py @@ -0,0 +1,75 @@ +import numpy as np +import pandas as pd +import pytest + +from pypfopt import risk_models + + +def estimate(data, method): + model = risk_models.CovarianceShrinkage(data, returns_data=True, frequency=1) + return ( + model.oracle_approximating() if method == "oas" else model.ledoit_wolf(method) + ) + + +@pytest.mark.parametrize( + "method", ["constant_variance", "single_factor", "constant_correlation", "oas"] +) +def test_shrinkage_uses_complete_observations(method): + rng = np.random.default_rng(52) + returns = pd.DataFrame(rng.normal(size=(120, 3)) / 100, columns=list("ABC")) + returns.loc[:79, "B"] = np.nan + returns.loc[95, "C"] = np.nan + original = returns.copy() + with pytest.warns(UserWarning, match="81.*incomplete"): + actual = estimate(returns, method) + expected = estimate(returns.dropna(), method) + pd.testing.assert_frame_equal(actual, expected) + pd.testing.assert_frame_equal(returns, original) + + +@pytest.mark.parametrize( + "method", ["constant_variance", "single_factor", "constant_correlation", "oas"] +) +@pytest.mark.parametrize("rows", [0, 1]) +def test_shrinkage_requires_two_complete_observations(method, rows): + returns = pd.DataFrame({"A": [0.1, 0.2, 0.3], "B": [np.nan, np.nan, np.nan]}) + if rows: + returns.loc[2, "B"] = 0.1 + with pytest.raises(ValueError, match="two complete"): + estimate(returns, method) + + +def test_semicovariance_uses_pairwise_observation_counts(): + returns = pd.DataFrame( + {"A": [-0.1, 0.2, np.nan, -0.2], "B": [np.nan, -0.1, -0.1, 0.4]} + ) + actual = risk_models.semicovariance( + returns, returns_data=True, benchmark=0, frequency=1 + ) + expected = pd.DataFrame( + [[1 / 60, 0], [0, 1 / 150]], index=returns.columns, columns=returns.columns + ) + pd.testing.assert_frame_equal(actual, expected) + + +def test_semicovariance_keeps_observed_zero_returns(): + returns = pd.DataFrame({"A": [-0.1, 0.0, -0.2], "B": [0.0, -0.1, 0.0]}) + actual = risk_models.semicovariance( + returns, returns_data=True, benchmark=0, frequency=1 + ) + np.testing.assert_allclose(actual, [[1 / 60, 0], [0, 1 / 300]]) + + +def test_semicovariance_rejects_pairs_without_joint_observations(): + returns = pd.DataFrame({"A": [0.1, np.nan], "B": [np.nan, 0.2]}) + with pytest.raises(ValueError, match="joint observation"): + risk_models.semicovariance(returns, returns_data=True) + + +def test_semicovariance_uses_joint_counts_off_diagonal(): + returns = pd.DataFrame({"A": [-0.1, -0.1, np.nan], "B": [-0.1, np.nan, -0.1]}) + actual = risk_models.semicovariance( + returns, returns_data=True, benchmark=0, frequency=1 + ) + np.testing.assert_allclose(actual, np.full((2, 2), 0.01)) From 43c7b5b45d6f06abd0fb5c94b12dde0b6e8c75c2 Mon Sep 17 00:00:00 2001 From: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> Date: Sat, 24 Oct 2026 18:11:51 +0500 Subject: [PATCH 2/3] Test pairwise PSD repair and document missing-data compatibility Signed-off-by: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> --- docs/RiskModels.rst | 26 ++++++++++++++ tests/test_missing_return_estimation.py | 45 +++++++++++++++++++++++++ 2 files changed, 71 insertions(+) diff --git a/docs/RiskModels.rst b/docs/RiskModels.rst index 5c56bf88..a897cb19 100644 --- a/docs/RiskModels.rst +++ b/docs/RiskModels.rst @@ -177,3 +177,29 @@ too much information. Manual shrinkage retains the pairwise sample covariance. observed zero returns. It raises an error if any pair has no joint observation. As with pairwise sample covariance, the resulting matrix may need the existing positive-semidefinite repair, which can alter individual entries. + +.. warning:: + + This changes results for inputs containing missing returns. Previous + versions treated missing returns as zero downside observations and divided + all semicovariance entries by the total number of rows. Saved portfolio + weights, reported risk and Sharpe ratios can therefore change materially; + recompute them and review the chosen missing-data policy. Fully observed + inputs retain the same semicovariance formula. + + On the repository stock-price fixture with ``benchmark=0`` and a 2% risk-free + rate, the long-only ``EfficientFrontier.max_sharpe`` example's reported risk + changes from approximately 0.0737 to 0.1733 and its Sharpe ratio from 3.48 to + 1.81. These are in-sample outputs for different optimized portfolios, using + a semicovariance risk matrix; they are not realized volatility, a universal + multiplier, or evidence of better out-of-sample performance. + +Pairwise deletion does not guarantee a positive-semidefinite matrix. For example, +the three return rows ``(-.1, -.1)``, ``(0, NaN)``, ``(NaN, 0)`` at benchmark zero +and frequency one produce raw entries ``[[.005, .01], [.01, .005]]``. Its +eigenvalues are -.005 and .015. The default spectral repair warns and sets the +negative eigenvalue to zero, producing .0075 in every entry: each diagonal +increases by .0025 (50%). Thus a repaired diagonal is no longer necessarily the +observed marginal downside moment. ``fix_method="diag"`` instead adds .0055 to +each diagonal in this example. Neither repair replaces a considered missing-data +model or restores observations that were never recorded. diff --git a/tests/test_missing_return_estimation.py b/tests/test_missing_return_estimation.py index 748b2d90..2b719a66 100644 --- a/tests/test_missing_return_estimation.py +++ b/tests/test_missing_return_estimation.py @@ -73,3 +73,48 @@ def test_semicovariance_uses_joint_counts_off_diagonal(): returns, returns_data=True, benchmark=0, frequency=1 ) np.testing.assert_allclose(actual, np.full((2, 2), 0.01)) + + +@pytest.mark.parametrize("frequency", [1, 252]) +@pytest.mark.parametrize( + "fix_method, expected", + [ + ("spectral", [[0.0075, 0.0075], [0.0075, 0.0075]]), + ("diag", [[0.0105, 0.01], [0.01, 0.0105]]), + ], +) +def test_semicovariance_repairs_indefinite_pairwise_estimate( + fix_method, expected, frequency +): + returns = pd.DataFrame({"A": [-0.1, 0.0, np.nan], "B": [-0.1, np.nan, 0.0]}) + original = returns.copy() + # Each variance uses two observations; the cross-product uses one. + # Raw matrix [[.005, .01], [.01, .005]] has eigenvalues -.005 and .015. + # Spectral clipping gives .015/2 in every entry (diagonal increases 50%). + # Diagonal repair instead adds 1.1 * .005 to both diagonal entries. + with pytest.warns(UserWarning, match="non positive semidefinite") as caught: + actual = risk_models.semicovariance( + returns, + returns_data=True, + benchmark=0, + frequency=frequency, + fix_method=fix_method, + ) + assert len(caught) == 1 + expected = pd.DataFrame( + np.array(expected) * frequency, index=returns.columns, columns=returns.columns + ) + pd.testing.assert_frame_equal(actual, expected) + assert np.linalg.eigvalsh(actual).min() >= -1e-12 + pd.testing.assert_frame_equal(returns, original) + + +def test_semicovariance_complete_data_retains_full_sample_formula(): + returns = pd.DataFrame({"A": [-0.1, 0.02, 0.0], "B": [0.03, -0.05, -0.02]}) + benchmark = 0.01 + downside = np.minimum(returns.to_numpy() - benchmark, 0) + expected = downside.T @ downside / len(returns) * 252 + actual = risk_models.semicovariance( + returns, returns_data=True, benchmark=benchmark, frequency=252 + ) + np.testing.assert_allclose(actual, expected, rtol=1e-12, atol=1e-12) From 2aef7fad180b664660ce8cb66b5c7dd3420f7d92 Mon Sep 17 00:00:00 2001 From: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> Date: Mon, 26 Oct 2026 19:59:03 +0500 Subject: [PATCH 3/3] DOC: clarify joint observation counts and PSD repair limits Signed-off-by: kadyrbekovhamit-cyber <288885044+kadyrbekovhamit-cyber@users.noreply.github.com> --- docs/RiskModels.rst | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/docs/RiskModels.rst b/docs/RiskModels.rst index a897cb19..3f7178e7 100644 --- a/docs/RiskModels.rst +++ b/docs/RiskModels.rst @@ -203,3 +203,13 @@ increases by .0025 (50%). Thus a repaired diagonal is no longer necessarily the observed marginal downside moment. ``fix_method="diag"`` instead adds .0055 to each diagonal in this example. Neither repair replaces a considered missing-data model or restores observations that were never recorded. + +There is no configurable minimum joint-observation count: a pair with one +shared return is allowed, while a pair with no shared returns raises an error. +Different entries can therefore be estimated from very different sample sizes. +In the example above, the off-diagonal uses one observation and each diagonal +uses two; dividing the off-diagonal by the marginal scales gives 2, not a valid +correlation. Spectral repair restores positive semidefiniteness but does not +make that one-observation estimate reliable or add missing information. Check +joint-observation counts and choose a missing-data policy before using the +matrix for portfolio optimization.