psd-covariance 1.0.5

Computes several covariance matrix estimators that ensure positive semi-definiteness (PSD).

A Package for Postive Definite Covariance Matrix Estimation (psd-covariance)

Introduction

psd-covariance provides fast and reliable tools for estimating positive semidefinite (PSD) covariance matrices, even in challenging high-dimensional or noisy settings. It is designed for practitioners and researchers who frequently encounter unstable, ill-conditioned, or non-PSD covariance matrices. At the core of the package are the Posterior Mean (PM) and Fixed-Trace Posterior Mean (PM-FT) estimators based on Bayesian eigenvalue regularization, as introduced in Boudt (2025).

Whether you are building portfolio optimization models, performing principal component analysis, or simply need a numerically robust PSD covariance estimate, psd-covariance provides a practical and user-friendly implementation of these state-of-the-art estimators.

Quick Start

Installation

pip install psd-covariance

Imports

import pandas as pd
import numpy as np
from psd_covariance import PosteriorMeanEstimator

Example 1: Transforming non-PSD matrices

We construct a non-PSD estimated covariance matrix with d=10, such that the smallest two eigenvalues are negative. We update these eigenvalues using the (fixed-trace) posterior mean estimators.

eigvals = [-0.50014538, -0.47905291,  0.59185356,  0.69919373,  0.8262933,   0.89063562, 1.005641,    1.39291291,  1.66503103,  1.95878406]
# Create matrix
d = len(eigvals)
A = np.random.randn(d, d)
Q, _ = np.linalg.qr(A)
Sigma_tilde = Q @ np.diag(eigvals) @ Q.T

# 1. Posterior Mean 
pm = PosteriorMeanEstimator(fixed_trace=False)
pm_result = pm.fit(Sigma_tilde, sigma=0.5)  
eigvals_pm = np.linalg.eigvalsh(pm_result.cov)

# 2. Posterior Mean Fixed Trace
ft = PosteriorMeanEstimator(fixed_trace=True)
ft_result = ft.fit(Sigma_tilde, sigma=0.5)
eigvals_ft = np.linalg.eigvalsh(ft_result.cov)

print("\nEigenvalues of posterior mean matrices:")
print("PM Estimator            :", eigvals_pm)
print("PM Estimator (Fixed Tr.):", eigvals_ft)
# Eigenvalues of cleaned matrices:
# PM Estimator            : [0.2625387  0.26678995 0.70412859 0.78083976 0.87984287 0.93304454 1.03263025 1.39704147 1.66581096 1.9588768 ]
# PM Estimator (Fixed Tr.): [0.21390763 0.21737141 0.5737001  0.63620176 0.71686614 0.76021306 0.84135212 1.13826203 1.35724629 1.5960264 ]

Example 2: PM and PM-FT Estimation using Cross-Validation

We generate a covariance matrix with a Toeplitz structure with d=10 and we draw n=20 observations from a Normal distribution with mean 0. We compute PM and PM-FT estimators using 10-fold cross validation.

# Generate data
np.random.seed(0)
d = 10
n = 20
rho = 0.8
cov_matrix = np.fromfunction(lambda i, j: rho ** np.abs(i - j), (d, d))
X = np.random.multivariate_normal(np.zeros(d), cov_matrix, size=n)
S = sample_cov(X)

reg_range = np.linspace(0.01, 2.0, 150)

# PM cross validation
pm = PosteriorMeanEstimator(fixed_trace=False)
reg_pm = pm.cross_validate_sigma(X, reg_range)
print(reg_pm)
# 0.2504026845637584
result_pm = pm.fit(S, reg_pm)

# FT cross validation
ft = PosteriorMeanEstimator(fixed_trace=True)
reg_ft = ft.cross_validate_sigma(X, reg_range)
print(reg_ft)
# 0.2771140939597316
result_ft = ft.fit(S, reg_ft)

Available Classes with Methods

PosteriorMeanEstimator

Input for posterior mean estimation: a covariance matrix Sigma_tilde : ndarray (dxd) together with the regularization parameter sigma : float, and the optional argument fixed_trace : bool indicating whether the fixed-trace variant (PM-FT) should be applied.

• PosteriorMeanEstimator(fixed_trace=False): constructor that selects PM (False) or PM-FT (True).
• fit(Sigma_tilde, sigma): computes the posterior mean covariance estimator using the chosen regularization parameter. Returns the named tuple PosteriorMeanResult(cov, inv, sigma), where cov is the estimated covariance matrix, inv its precision matrix, and sigma the regularization parameter used in the estimator.
• cross_validate_sigma(X, regularization_range, n_splits=10, n_jobs=-1): selects the optimal regularization parameter by evaluating each value in regularization_range using K-fold predictive log-likelihood cross-validation. Here X is a data matrix, regularization_range is an array of candidate regularization values, n_splits sets the number of folds (default 10), and n_jobs controls parallel computing (-1 uses all CPU cores). Method returns optimal regularization value.

EigenvalueCleaning

Input for eigenvalue cleaning methods: a covariance matrix cov : ndarray (dxd) together with the optional parameters epsilon : float (replacement value for negative eigenvalues), fixed_trace : bool (whether to preserve the trace of the original matrix), and positive_definite : bool (whether to enforce strict positive definiteness).

• threshold_negative(cov, fixed_trace=False). Sets all negative eigenvalues to zero.
• replace_negative(cov, epsilon=1e-4, fixed_trace=False, positive_definite=False) Replaces negative eigenvalues with epsilon.
• absolute_negative(cov, fixed_trace=False, positive_definite=False). Takes the absolute value of all eigenvalues.

All eigenvalue cleaning methods return return the named tuple CleanedCovariance(cov, inv), with cov the estimated covariance matrix and inv its precision matrix.

ShrinkageEstimator

Input for shrinkage estimators: data matrix of observations X : pandas.DataFrame (nxd).

• linear_shrinkage(X): Ledoit & Wolf (2002) linear shrinkage. Returns covariance, precision, and shrinkage intensity alpha.
• quadratic_inverse_shrinkage(X): Ledoit & Wolf (2022) quadratic inverse shrinkage. Returns covariance and precision matrices.

All shrinkage estimators return ShrinkageResult(cov, inv, shrinkage), with cov the estimated covariance matrix, inv its precision matrix, and shrinkage the shrinkage intensity (only for linear shrinkage; None for QIS).

Authors

This package is based on the paper 'Well-Conditioned Covariance Estimation via Bayesian Eigenvalue Regularization', by Kris Boudt, Jesper Cremers, Kirill Dragun & Steven Vanduffel. The psd-covariance package is developed and maintained by Jesper Cremers.

References

• Boudt, K., J. Cremers, K. Dragun, and S. Vanduffel (2025). Well-conditioned covariance estimation via bayesian eigenvalue regularization. Working paper.
• Ledoit, O. and M. Wolf (2004). Honey, I shrunk the sample covariance matrix. The Journal of Portfolio Management 30 (4), 110-119.
• Ledoit, O. and M. Wolf (2022). Quadratic shrinkage for large covariance matrices. Bernoulli 28 (3), 1519-1547.
• Rousseeuw, P. J. and G. Molenberghs (1993). Transformation of non positive semidefinite correlation matrices. Communications in Statistics - Theory and Methods 22 (4), 965-984.
• Parts of the ShrinkageEstimator class contain adapted code from Michael Wolf's reference implementation: https://github.com/pald22/covShrinkage.