"""Empirical Prior estimator."""
# Copyright (c) 2023-2026
# Author: Hugo Delatte <hugo.delatte@skfoliolabs.com>
# SPDX-License-Identifier: BSD-3-Clause
from __future__ import annotations
import numbers
import warnings
import numpy as np
import sklearn.utils.metadata_routing as skm
import sklearn.utils.validation as skv
from skfolio.moments import BaseCovariance, BaseMu, EmpiricalCovariance, EmpiricalMu
from skfolio.prior._base import BasePrior
from skfolio.prior._model import ReturnDistribution
from skfolio.typing import ArrayLike, BoolArray, StrArray
from skfolio.utils._array_buffer import _ArrayBuffer
from skfolio.utils.tools import _call_estimator, check_estimator
_FITTED_ATTR = "return_distribution_"
# Fraction of zero-filled observations above which an investable asset triggers a
# UserWarning about understated risk in scenario-based measures.
_ZERO_FILL_WARNING_THRESHOLD = 0.05
# Maximum number of assets named in the zero-fill warning message.
_ZERO_FILL_WARNING_MAX_ASSETS = 10
[docs]
class EmpiricalPrior(BasePrior):
"""Empirical Prior estimator.
The Empirical Prior estimates the :class:`~skfolio.prior.ReturnDistribution` by
fitting a `mu_estimator` and a `covariance_estimator` separately.
**NaN handling:**
Missing data (NaN returns) caused by late listings, delistings and holidays is
accepted when both `mu_estimator` and `covariance_estimator` support it (for
example :class:`~skfolio.moments.EWMu` and
:class:`~skfolio.moments.EWCovariance`). The moment estimators receive the data
unchanged and apply their own NaN treatment.
In `return_distribution_.returns`, the scenario columns of non-investable
assets (NaN in the estimated `mu` and/or covariance diagonal) are left
unchanged and are removed downstream by
:meth:`~skfolio.prior.ReturnDistribution.investable_subset`. Missing
observations of investable assets are replaced by zero.
Zero-filling long gaps, such as the pre-listing history of a late-listed
asset, understates its risk in scenario-based measures (CVaR, EVaR, CDaR,
worst realization, ...). A `UserWarning` is emitted when more than 5% of an
investable asset's scenario history is zero-filled. The moments estimation is
not affected. To reduce the zero-filled share, set `max_history` or use a
factor model prior such as :class:`~skfolio.prior.CharacteristicsFactorModel`.
Parameters
----------
mu_estimator : BaseMu, optional
The assets :ref:`expected returns estimator <mu_estimator>`.
The default (`None`) is to use :class:`~skfolio.moments.EmpiricalMu`.
covariance_estimator : BaseCovariance , optional
The assets :ref:`covariance matrix estimator <covariance_estimator>`.
The default (`None`) is to use :class:`~skfolio.moments.EmpiricalCovariance`.
is_log_normal : bool, default=False
If this is set to True, the moments are estimated on the logarithmic returns
as opposed to the linear returns. Then the moments estimations of the
logarithmic returns are projected to the investment horizon and transformed
to obtain the moments estimation of the linear returns at the investment
horizon. If True, `investment_horizon` must be provided. The input `X` must be
**linear returns**. They will be converted into logarithmic returns only for the
moments estimation.
.. seealso::
:ref:`data preparation <data_preparation>`
investment_horizon : float, optional
The investment horizon used for the moments estimation of the linear returns
when `is_log_normal` is `True`.
max_history : int, optional
Maximum number of observations to keep in `return_distribution_.returns`.
This is useful for controlling memory usage during incremental learning
with :meth:`partial_fit`.
* If `None` (default), all returns are accumulated.
* If an integer, only the last `max_history` observations are kept
(rolling window).
Attributes
----------
return_distribution_ : ReturnDistribution
Fitted :class:`~skfolio.prior.ReturnDistribution` to be used by the optimization
estimators, containing the asset returns distribution and moments estimation.
mu_estimator_ : BaseMu
Fitted `mu_estimator`.
covariance_estimator_ : BaseCovariance
Fitted `covariance_estimator`.
n_features_in_ : int
Number of assets seen during `fit`.
feature_names_in_ : ndarray of shape (`n_features_in_`,)
Names of features seen during `fit`. Defined only when `X`
has feature names that are all strings.
References
----------
.. [1] "Linear vs. Compounded Returns - Common Pitfalls in Portfolio Management".
GARP Risk Professional.
Attilio Meucci (2010).
"""
mu_estimator_: BaseMu
covariance_estimator_: BaseCovariance
n_features_in_: int
feature_names_in_: StrArray
def __init__(
self,
mu_estimator: BaseMu | None = None,
covariance_estimator: BaseCovariance | None = None,
is_log_normal: bool = False,
investment_horizon: float | None = None,
max_history: int | None = None,
):
self.mu_estimator = mu_estimator
self.covariance_estimator = covariance_estimator
self.is_log_normal = is_log_normal
self.investment_horizon = investment_horizon
self.max_history = max_history
[docs]
def fit(self, X: ArrayLike, y=None, **fit_params) -> EmpiricalPrior:
"""Fit the Empirical Prior estimator.
Parameters
----------
X : array-like of shape (n_observations, n_assets)
Price returns of the assets. May contain NaN (holidays, late
listings, delistings) when both `mu_estimator` and
`covariance_estimator` handle missing data.
y : Ignored
Not used, present for API consistency by convention.
**fit_params : dict
Parameters to pass to the underlying estimators.
Only available if `enable_metadata_routing=True`, which can be
set by using `sklearn.set_config(enable_metadata_routing=True)`.
See :ref:`Metadata Routing User Guide <metadata_routing>` for
more details.
Returns
-------
self : EmpiricalPrior
Fitted estimator.
"""
self._reset()
return self._fit(X, y, method="fit", **fit_params)
[docs]
def partial_fit(self, X: ArrayLike, y=None, **fit_params) -> EmpiricalPrior:
"""Incrementally fit the Empirical Prior estimator.
This method allows for streaming/online updates to the prior estimate.
Each call updates the internal state with new observations.
Both `mu_estimator` and `covariance_estimator` must implement
`partial_fit` for this method to work.
Parameters
----------
X : array-like of shape (n_observations, n_assets)
Price returns of the assets. May contain NaN (holidays, late
listings, delistings) when both `mu_estimator` and
`covariance_estimator` handle missing data.
y : Ignored
Not used, present for API consistency by convention.
**fit_params : dict
Parameters to pass to the underlying estimators.
Only available if `enable_metadata_routing=True`, which can be
set by using `sklearn.set_config(enable_metadata_routing=True)`.
See :ref:`Metadata Routing User Guide <metadata_routing>` for
more details.
Returns
-------
self : EmpiricalPrior
Fitted estimator.
"""
return self._fit(X, y, method="partial_fit", **fit_params)
def _fit(self, X: ArrayLike, y, method: str, **fit_params) -> EmpiricalPrior:
"""Core fitting logic shared by fit and partial_fit.
Parameters
----------
X : array-like of shape (n_observations, n_assets)
Price returns of the assets.
y : Ignored
Not used, present for API consistency by convention.
method : str
Either "fit" or "partial_fit". Determines which method to call
on sub-estimators and how to handle returns accumulation.
**fit_params : dict
Parameters to pass to the underlying estimators.
Returns
-------
self : EmpiricalPrior
Fitted estimator.
"""
routed_params = skm.process_routing(self, method, **fit_params)
first_call = not hasattr(self, _FITTED_ATTR)
if first_call:
self._validate_params()
self._initialize()
if self.is_log_normal:
X_fit = np.log(1 + X)
y_fit = np.log(1 + y) if y is not None else None
else:
X_fit = X
y_fit = y
# Fit or partial_fit the mu estimator
_call_estimator(
self.mu_estimator_,
method,
X_fit,
y_fit,
routed_params=routed_params.mu_estimator,
)
# Fit or partial_fit the cov estimator
_call_estimator(
self.covariance_estimator_,
method,
X_fit,
y_fit,
routed_params=routed_params.covariance_estimator,
)
mu = self.mu_estimator_.mu_
covariance = self.covariance_estimator_.covariance_
# Transform log moments to linear if needed
if self.is_log_normal:
mu *= self.investment_horizon
covariance *= self.investment_horizon
# Convert to linear returns distribution
mu = np.exp(mu + 0.5 * np.diag(covariance))
covariance = np.outer(mu, mu) * (np.exp(covariance) - 1)
mu -= 1
# we validate and convert to numpy after all models have been fitted to keep
# features names information.
X = skv.validate_data(self, X, ensure_all_finite="allow-nan", reset=first_call)
# Accumulate returns with amortized O(1) appends
if first_call:
self._returns_buffer = _ArrayBuffer()
self._returns_buffer.append(X)
if self.max_history is not None:
self._returns_buffer.truncate_to_last(self.max_history)
# Zero-fill missing observations of investable assets so that the return
# scenarios stay usable by scenario-based risk measures. Columns of
# non-investable assets are left unchanged and are removed downstream by
# `ReturnDistribution.investable_subset`.
returns = self._returns_buffer.array
investable = np.isfinite(mu) & np.isfinite(np.diag(covariance))
missing = np.isnan(returns) & investable
if missing.any():
self._warn_zero_fill(missing)
returns = returns.copy()
returns[missing] = 0.0
self.return_distribution_ = ReturnDistribution(
mu=mu,
covariance=covariance,
returns=returns,
)
return self
def _warn_zero_fill(self, missing: BoolArray) -> None:
"""Warn when the zero-filled share of an investable asset is material.
A `UserWarning` is emitted for assets whose zero-filled fraction of the
current scenario history exceeds 5%. Each asset is reported once per
fitted state so that streaming `partial_fit` calls do not repeat the
warning.
Parameters
----------
missing : ndarray of shape (n_observations, n_assets)
Boolean mask of the zero-filled entries.
"""
fill_ratios = missing.mean(axis=0)
new = [
i
for i in np.flatnonzero(fill_ratios > _ZERO_FILL_WARNING_THRESHOLD)
if i not in self._zero_fill_warned_assets
]
if not new:
return
self._zero_fill_warned_assets.update(new)
names = getattr(self, "feature_names_in_", np.arange(missing.shape[1]))
details = ", ".join(
f"{names[i]} ({fill_ratios[i]:.1%})"
for i in new[:_ZERO_FILL_WARNING_MAX_ASSETS]
)
if len(new) > _ZERO_FILL_WARNING_MAX_ASSETS:
details += f", ... ({len(new)} assets in total)"
warnings.warn(
f"More than {_ZERO_FILL_WARNING_THRESHOLD:.0%} of the return scenarios "
f"of the following assets are zero-filled missing observations: "
f"{details}. This understates their risk in scenario-based measures "
"(CVaR, EVaR, CDaR, ...) but does not affect `mu` and `covariance`. "
"Set `max_history` or use a factor model prior to reduce the "
"zero-filled share.",
stacklevel=2,
)
def _validate_params(self) -> None:
"""Validate parameters."""
if self.is_log_normal:
if self.investment_horizon is None:
raise ValueError(
"`investment_horizon` must be provided when "
"`is_log_normal` is `True`"
)
else:
if self.investment_horizon is not None:
raise ValueError(
"`investment_horizon` must be `None` when "
"`is_log_normal` is `False`"
)
if self.max_history is not None:
if isinstance(self.max_history, bool) or not isinstance(
self.max_history, numbers.Integral
):
raise ValueError(
f"`max_history` must be a positive integer or None, "
f"got {self.max_history}"
)
if self.max_history < 1:
raise ValueError(
f"`max_history` must be a positive integer or None, "
f"got {self.max_history}"
)
def _initialize(self):
self._zero_fill_warned_assets = set()
self.mu_estimator_ = check_estimator(
self.mu_estimator,
default=EmpiricalMu(),
check_type=BaseMu,
)
self.covariance_estimator_ = check_estimator(
self.covariance_estimator,
default=EmpiricalCovariance(),
check_type=BaseCovariance,
)
def _reset(self) -> None:
"""Reset fitted state."""
if hasattr(self, _FITTED_ATTR):
delattr(self, _FITTED_ATTR)