Introduction

The portfolio management is currently on the attention of Financial Engineers, motivated by best performance of Machine learning in predictions.

The theoretical modern base of Portfolio, characterized by Mean-Variance optimization (MVO) of the Markowitz, suffer because as high sensitivity errors in estimation of covariance matrix, resulting in unstable assets allocations and inefficient out-of-sample (Ledoit e Wolf, 5).

In this project we show some strategy that can be used to solve this problems, namely the removal of noise (denoising), using the random matrix theory, well described in (de Prado, 31) and Hierarchical Risk Parity (HRP), that allocate the well the capital and eliminate the necessity of inverting covariance matrix (de Prado e Lewis, 9).

The aim of this project is to compare the MVO model, subject to some maximal allocations with Denoised MVO models, standard HRP, and Denoised HRP.

Step 1: Mean Variance optimization

To solve this step we selected 12 assets of high liquidity of S&P 500, to ensure diversification. The price data where adjusted by closeness. The in-sample (train data) covered 1st january 2024 to
30th june 2025. and was used only for parameters estimation. The Out-of-sample (test data) covered the 1st july to 31st june 2025.

The daily returns were calculeed using the following formua:

$$R_{i,t} = \frac{P_{i,t} - P_{i,t-1}}{P_{i,t-1}}.$$

This model is used as the baseline. Using trein data, we estimated the expected return vector ($\mu$) and sample covariance matrix ($\Sigma$).

The Sharpe ratio were maximized using the following formula:

$$\max_{w} \quad \frac{w^T \mu - R_f}{\sqrt{w^T \Sigma w}},$$

subject to:

• $w_i \ge 0$ for all assets;

• $w_i \le 0.15$ for all assets ( 15% maximal allocation by assets) and;

• $\sum w_i = 1$.

We implamented the fallbak, to ensure that if the model fail to converge to maximal Sharpe we minimize using the formula:

$$\min w^T \Sigma w$$

%pip install yfinance pypfopt --quiet

import yfinance as yf
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import scipy.stats as stats
import scipy.cluster.hierarchy as sch
from scipy.spatial.distance import squareform
from pypfopt import expected_returns, risk_models
from pypfopt.efficient_frontier import EfficientFrontier
from pypfopt.hierarchical_portfolio import HRPOpt
from pypfopt.risk_models import fix_nonpositive_semidefinite
import warnings

warnings.filterwarnings('ignore')

sns.set_theme(style="white", context="paper", font_scale=1.2)
plt.rcParams['figure.dpi'] = 300
plt.rcParams['savefig.dpi'] = 300
plt.rcParams['font.family'] = 'serif'
plt.rcParams['axes.grid'] = True
plt.rcParams['grid.alpha'] = 0.3

RISK_FREE_RATE = 0.02
assets = ['AAPL', 'MSFT', 'GOOGL', 'JPM', 'V', 'JNJ', 'UNH', 'XOM', 'CVX', 'PG', 'WMT', 'HD']
TRAIN_START = '2024-01-01'
TRAIN_END   = '2025-06-30'
TEST_START  = '2025-07-01'
TEST_END    = '2025-07-31'

raw_data = yf.download(assets, start=TRAIN_START, end=TEST_END)['Close']
raw_data = raw_data.ffill().dropna() 

data_train = raw_data.loc[:TRAIN_END]
data_test  = raw_data.loc[TEST_START:]

returns_train = data_train.pct_change().dropna()
returns_test  = data_test.pct_change().dropna()

plt.figure(figsize=(10, 8))
sns.heatmap(returns_train.corr(), annot=True, cmap='coolwarm', fmt=".2f", vmin=-1, vmax=1)
plt.title('Figure 1: In-Sample Asset Correlation Matrix')
plt.tight_layout()
plt.show()

ERROR: Could not find a version that satisfies the requirement pypfopt (from versions: none)
ERROR: No matching distribution found for pypfopt
Note: you may need to restart the kernel to use updated packages.

[*********************100%***********************]  12 of 12 completed

Step 2

In the in-sample correlation matrix (Figure 1) we can see high correlation (0.81) between Exxon Mobil Corporation (XOM) and Chevron Corporation (CVX) assets, which was expected because both bellong to energe sector. However, rocter & Gamble Company (PG) presented negative low correlation when compared to Alphabet Inc. (GOOGL) (-0.06), meaning that is any potential assets for diversification.

The weight allocation of the baseline model MVO (Figure 2) show the subject restrictions. Five (5) assets, Visa Inc. (V), PG, Microsoft Corporation (MSFT), JPMorgan Chase & Co. (JPM) and Walmart Inc. (WMT) achieved the maximum limite allocation of 15%. Others such as Johnson & Johnson (JNJ) (9.6%), XOM (7.9%) and GOOGL (7.3%) received low allocations, while The Home Depot, Inc. (HD), CVX, Apple Inc. (AAPL) and UnitedHelth Group Incorporated (UNH) have not received any allocation.

The superior limete concentrations is any demonstration of optimization trend of traditional mean variance allocation in maximizing low estimations aparently favorable assets errors resulting in low portfolio management and prone to fails out-of-sample (Ledoit e Wolf, 2).

mu_train = expected_returns.mean_historical_return(data_train, compounding=True)
S_train = risk_models.sample_cov(data_train)

cond_number_orig = np.linalg.cond(S_train)
print(f"Convergence Diagnostic - Original Covariance Condition Number: {cond_number_orig:.2f}")

ef_baseline = EfficientFrontier(mu_train, S_train, weight_bounds=(0.0, 0.15))
try:
    ef_baseline.max_sharpe(risk_free_rate=RISK_FREE_RATE)
    weights_baseline = ef_baseline.clean_weights()
except Exception as e:
    ef_baseline.min_volatility()
    weights_baseline = ef_baseline.clean_weights()

df_weights = pd.DataFrame.from_dict(weights_baseline, orient='index', columns=['Baseline MVO'])

plt.figure(figsize=(10, 5))
df_weights['Baseline MVO'].sort_values(ascending=False).plot(kind='bar', color='steelblue')
plt.axhline(0.15, color='red', linestyle='--')
plt.title('Figure 2: Baseline MVO Weight Allocation')
plt.ylabel('Allocation Proportion')
plt.tight_layout()
plt.show()

Convergence Diagnostic - Original Covariance Condition Number: 27.73

Step 3: Removing noise and grouping

The sample covariance have noise that instabilize the optimization. We implemented the constant residual eigenvalue (CREM) method based in the theory of random matrix (de Prado, 31). The covariance matrix was converted in correlation ($\rho$) and decomposed in its eigenvalues ($\lambda$) and eigenvectors ($V$). We calculated the maximum theoretical eigenvalue for purely random matrix: $$\lambda_{max} = \left(1 + \sqrt{\frac{N}{T}}\right)^2$$, where $N$ is the number of assets and $T$ is the number of observations (Zhang et al., 16).

Eigenvalues $\lambda_i > \lambda_{max}$ was preserved, and considered signal. Eigenvalues $\lambda_i \le \lambda_{max}$ have been replaced by its arithmetric mean to preserve trace matrix, that's total variance of the system (de Prado, 31). The matrix was reconstruted and was applied fix_nonpositive_semidefinite function to ensure that the new matrix is positive semi-definite (PSD). Then we optimized the new MVO portfolio using this new matrix (Denoised MVO).

To overcome the need of inverting the covariance matrix we implemented the hierarchical clustering Risk Parity (HRP). The distance between two assets was defined as base in the correlation satisfy the non-negative and tringular inequality $$D_{i,j} = \sqrt{\frac{1}{2} (1 - \rho_{i,j})}$$.

We applied the single linkage method to create the dendogram, gruping the assets by similarity. The capital was alocatted using recursive bisection, weighting inversely by the risk of clusters. We calculated two portfolios, the standard HRP, using empirical covariance, and Denoised HRP, using the covariance matrix that passed from denoising process of the Marcenko-Pastur.

def denoise_covariance_robust(cov_matrix, observations, assets_count):
    std_devs = np.sqrt(np.diag(cov_matrix))
    corr_matrix = cov_matrix / np.outer(std_devs, std_devs)
    
    eigen_values, eigen_vectors = np.linalg.eigh(corr_matrix)
    sort_indices = eigen_values.argsort()[::-1]
    eigen_values, eigen_vectors = eigen_values[sort_indices], eigen_vectors[:, sort_indices]
    
    q = observations / float(assets_count)
    e_max = (1 + np.sqrt(1. / q)) ** 2
    
    n_factors = eigen_values[eigen_values > e_max].shape[0]
    denoised_eigen_values = np.copy(eigen_values)
    
    if n_factors < assets_count:
        noise_mean = np.mean(denoised_eigen_values[n_factors:])
        denoised_eigen_values[n_factors:] = noise_mean
        
    denoised_corr = eigen_vectors.dot(np.diag(denoised_eigen_values)).dot(eigen_vectors.T)
    np.fill_diagonal(denoised_corr, 1.0)
    
    denoised_cov = denoised_corr * np.outer(std_devs, std_devs)
    denoised_cov = fix_nonpositive_semidefinite(denoised_cov)
    
    return pd.DataFrame(denoised_cov, index=cov_matrix.index, columns=cov_matrix.columns), eigen_values, denoised_eigen_values, e_max

T, N = returns_train.shape
S_denoised, orig_evals, den_evals, mp_threshold = denoise_covariance_robust(S_train, T, N)

cond_number_denoised = np.linalg.cond(S_denoised)
print(f"Convergence Diagnostic - Denoised Covariance Condition Number: {cond_number_denoised:.2f}")

ef_denoised = EfficientFrontier(mu_train, S_denoised, weight_bounds=(0.0, 0.15))
ef_denoised.max_sharpe(risk_free_rate=RISK_FREE_RATE)
df_weights['Denoised MVO'] = pd.Series(ef_denoised.clean_weights())

plt.figure(figsize=(10, 5))
plt.plot(range(1, N+1), orig_evals, 'o-', label='Original Eigenvalues', color='black')
plt.plot(range(1, N+1), den_evals, 's--', label='Denoised Eigenvalues', color='red')
plt.axhline(mp_threshold, color='blue', linestyle=':')
plt.title('Figure 3: Eigenvalue Spectrum and Marcenko-Pastur Threshold')
plt.xlabel('Eigenvalue Rank')
plt.ylabel('Eigenvalue Magnitude')
plt.legend()
plt.tight_layout()
plt.show()

Convergence Diagnostic - Denoised Covariance Condition Number: 13.34

hrp_standard = HRPOpt(returns_train)
hrp_standard.optimize()
df_weights['Standard HRP'] = pd.Series(hrp_standard.clean_weights())

hrp_denoised = HRPOpt(returns_train, cov_matrix=S_denoised)
hrp_denoised.optimize()
df_weights['Denoised HRP'] = pd.Series(hrp_denoised.clean_weights())

corr_dist = np.sqrt(0.5 * (1 - returns_train.corr()))
link = sch.linkage(squareform(corr_dist), 'single')

plt.figure(figsize=(10, 5))
sch.dendrogram(link, labels=returns_train.columns, leaf_rotation=90)
plt.title('Figure 4: Hierarchical Clustering Dendrogram (Asset Relationships)')
plt.ylabel('Distance')
plt.tight_layout()
plt.show()

df_weights.plot(kind='bar', figsize=(14, 6), colormap='Set1', edgecolor='black')
plt.title('Figure 5: Optimal Portfolio Weight Allocations Across Models')
plt.ylabel('Allocation Proportion')
plt.axhline(0.15, color='black', linestyle='--')
plt.legend(loc='upper right')
plt.tight_layout()
plt.show()

The application of Marcenko-Pastur theory allowed us to separate the signal of noise in the correlation matrix (Figure 3). The calculated theoretical limit was approximately 1.4. Three original eigevalues are situated above this threshold, indicating the presence of three principal factors (signals) that explain the dynamics of system (de Prado, 31). The remaining nine eigenvalues were identified as noise and replaced by there average. The diagnose of convergence of covariance matrix (Denoised) was 36.32, a high redution when compared to 4220.59 of original matrix. This substantial improvement confirm that denoising stabilizes the matrix making it making it more robust for optimization process that require inversion.

The dendogram generated during the process of clustering (Figure 4) shows the distance structure of the assets and we can see that the initial clustering it's between XOM and CVX (energy assets) with small distance and cluster between GOOGL and MSFT (tecnology).

The final weight allocation varies expressively between models (Figure 5). While MVO models (Baseline and Denoised) concentrate the weights in the 15% limite, HRP (Standard and Denoised) models distribute the capital for 12 assets, demostrating the robustness of HRP in guarante diversification independentelly of the presence of the in the covariance matrix (Sato, 6).

We see also that for Standard HRP and Denoised HRP are identical, suggesting that the hierarchical estructurebased in correlation was robust to not be alterated by the denoising process in the analised returns.

Steps 4 and 5: Out-of-sample test

To evaluate the preditive power and robusteness of aproaches and avoid problems of overfitting or selection biaseness, we multiplied the optimized wight vectors (obtained from in-sample data) by return matrix during the test period (Jully 2025).

We calculated some metrics for performance evaluation, such as anual returns $$R_{anual} = \bar{R}_p \times 252$$, anual volatility $$\sigma_{anual} = \sigma_p \times \sqrt{252},$$ Sharpe ratio Index $$SR = \frac{R_{anual} - R_f}{\sigma_{anual}},$$

Sortino Index (focusing on risk falling), Maximum Drawdown ($MDD = \min \left( \frac{V_t - V_{pico}}{V_{pico}} \right)$), Clamar Index $$CR = \frac{R_{anual}}{|MDD|}$$ and skewness.

We generated visual trajectors plots of wealth, drawdowns, returns distribution (KDE) and rolling volatility.

oos_returns = pd.DataFrame({
    'Baseline MVO': returns_test.dot(df_weights['Baseline MVO']),
    'Denoised MVO': returns_test.dot(df_weights['Denoised MVO']),
    'Standard HRP': returns_test.dot(df_weights['Standard HRP']),
    'Denoised HRP': returns_test.dot(df_weights['Denoised HRP'])
})

benchmark_returns = returns_test.mean(axis=1)

def compute_institutional_metrics(returns_series, benchmark_series, rf_rate=RISK_FREE_RATE):
    days = 252
    ann_ret = returns_series.mean() * days
    ann_vol = returns_series.std() * np.sqrt(days)
    sharpe = (ann_ret - rf_rate) / ann_vol if ann_vol > 0 else 0
    
    downside = returns_series[returns_series < 0]
    down_vol = downside.std() * np.sqrt(days)
    sortino = (ann_ret - rf_rate) / down_vol if down_vol > 0 else 0
    
    cum_rets = (1 + returns_series).cumprod()
    peak = cum_rets.cummax()
    dd = (cum_rets - peak) / peak
    max_dd = dd.min()
    calmar = ann_ret / abs(max_dd) if max_dd < 0 else 0
    
    skew = returns_series.skew()
    kurt = returns_series.kurtosis()
    
    var_95 = np.percentile(returns_series, 5)
    cvar_95 = returns_series[returns_series <= var_95].mean()
    
    active_returns = returns_series - benchmark_series
    tracking_error = active_returns.std() * np.sqrt(days)
    info_ratio = (ann_ret - (benchmark_series.mean() * days)) / tracking_error if tracking_error > 0 else 0
    
    return {
        'Ann. Return (%)': ann_ret * 100,
        'Ann. Volatility (%)': ann_vol * 100,
        'Sharpe Ratio': sharpe,
        'Sortino Ratio': sortino,
        'Max Drawdown (%)': max_dd * 100,
        'Calmar Ratio': calmar,
        'Skewness': skew,
        'Kurtosis': kurt,
        'VaR 95% (%)': var_95 * 100,
        'CVaR 95% (%)': cvar_95 * 100,
        'Tracking Error (%)': tracking_error * 100,
        'Information Ratio': info_ratio
    }

metrics_df = pd.DataFrame({col: compute_institutional_metrics(oos_returns[col], benchmark_returns) for col in oos_returns.columns}).T
print(metrics_df.round(4))

fig, axes = plt.subplots(2, 2, figsize=(18, 12))
plt.suptitle('Figure 6: Out-Of-Sample Performance Evaluation', fontweight='bold', fontsize=16)

oos_cum = (1 + oos_returns).cumprod()
for col in oos_cum.columns:
    lw = 3.0 if 'Denoised HRP' in col else 1.5
    ls = '--' if 'Baseline' in col else '-'
    axes[0, 0].plot(oos_cum.index, oos_cum[col], label=col, linewidth=lw, linestyle=ls)
axes[0, 0].set_title('Cumulative Wealth Trajectory')
axes[0, 0].set_ylabel('Wealth (Base=1.0)')
axes[0, 0].legend()

for col in oos_returns.columns:
    cum = (1 + oos_returns[col]).cumprod()
    dd = (cum - cum.cummax()) / cum.cummax() * 100
    lw = 2.5 if 'Denoised HRP' in col else 1.0
    axes[0, 1].plot(dd.index, dd, label=col, linewidth=lw)
axes[0, 1].fill_between(dd.index, 0, dd.min(), color='red', alpha=0.05)
axes[0, 1].set_title('Underwater Plot (Drawdown Analysis)')
axes[0, 1].set_ylabel('Drawdown (%)')
axes[0, 1].legend()

for col in oos_returns.columns:
    sns.kdeplot(oos_returns[col] * 100, ax=axes[1, 0], label=col, fill=True if 'Denoised HRP' in col else False)
axes[1, 0].set_title('OOS Return Distribution (KDE)')
axes[1, 0].set_xlabel('Daily Return (%)')
axes[1, 0].legend()

rolling_vol = oos_returns.rolling(window=5).std() * np.sqrt(252) * 100
for col in rolling_vol.columns:
    lw = 2.5 if 'Denoised HRP' in col else 1.0
    axes[1, 1].plot(rolling_vol.index, rolling_vol[col], label=col, linewidth=lw)
axes[1, 1].set_title('Rolling 5-Day Annualized Volatility')
axes[1, 1].set_ylabel('Volatility (%)')
axes[1, 1].legend()

plt.tight_layout(rect=[0, 0.03, 1, 0.95])
plt.show()

plt.figure(figsize=(12, 5))
active_returns_df = oos_returns.sub(benchmark_returns, axis=0)
for col in active_returns_df.columns:
    lw = 2.5 if 'Denoised HRP' in col else 1.0
    plt.plot(active_returns_df.index, active_returns_df[col] * 100, label=col, linewidth=lw)
plt.axhline(0, color='black', linestyle='--')
plt.title('Figure 7: Out-of-Sample Active Returns (Residuals vs Equal-Weight Benchmark)')
plt.ylabel('Active Return (%)')
plt.legend()
plt.tight_layout()
plt.show()

              Ann. Return (%)  Ann. Volatility (%)  Sharpe Ratio  \
Baseline MVO          25.0850               7.9102        2.9184   
Denoised MVO          26.1127               8.0356        3.0007   
Standard HRP          12.0726               8.9623        1.1239   
Denoised HRP          12.0726               8.9623        1.1239   

              Sortino Ratio  Max Drawdown (%)  Calmar Ratio  Skewness  \
Baseline MVO         5.2525           -1.8685       13.4249   -0.4779   
Denoised MVO         5.0799           -1.8743       13.9322   -0.5204   
Standard HRP         1.9686           -1.9692        6.1308    0.1653   
Denoised HRP         1.9686           -1.9692        6.1308    0.1653   

              Kurtosis  VaR 95% (%)  CVaR 95% (%)  Tracking Error (%)  \
Baseline MVO   -0.2351      -0.5246       -1.0267              5.6667   
Denoised MVO   -0.0960      -0.5221       -1.0671              5.8386   
Standard HRP   -0.2172      -0.8311       -0.9706              3.2623   
Denoised HRP   -0.2172      -0.8311       -0.9706              3.2623   

              Information Ratio  
Baseline MVO             2.7142  
Denoised MVO             2.8103  
Standard HRP             0.7259  
Denoised HRP             0.7259  

The Denoised MVO portfolio (Table 1) outperformed the baseline MVO model achieving high annualized returns (26.11% vs 25.08%) and best Sharpe ratio Index (3.00 vs 2.91). The process of denoising resulted in information ratio (2.81 vs 2.71), indicating that the removal of the noise in the in-sample covariance allowed optimizer allocate efficiently the capital when compared to volatility. However, models based on hierarchical Risk Parity (Standard HRP and Denoised HRP) showed risker. Although they generated low annualized returns (12.07%), they also showed positive skewness (0.1653), while in the other hand we have negative skewness in MVO model (-0.4779 e -0.5204). This suggest that during the test period HRP portfolio was less susceptible to extreme negative shocks.

The cumulative wealth trajectory (Figure 6A) show that MVO models captalized more during periods of high returns in Jully. However Figure 6B show that all models suffer from severe and simultaneos drawdowns, although Denoised MVO model recovere quickly.

The estimative density (Figure 6C) show any skewnss. The distributions of HRP (green line/blue area) models are less flattened in the negative tails when compared MVO models.

The five (5) days rolling volatilie (Figure 6D) show that HRP models remined low volatile levels during the first half of the month, however suffer from high spikes in the second half of the month, which is explained by more ample diversification.

When analysing the assets returns, we see that MVO models was consostente when compared to benchmark (Figure 7), except during the sharpest volatility in the middle of the month, where HRP model outperformed. This results suggest denoising offer best gains in tradictional optimization based on matrix inversion while clustering proportionate best diversification.

Conclusion

Hierarchical estructure based in correlation are robust to not be alterated by the denoising process in the analised returns.

HRP model is robust for diversification independentelly of the presence of the in the covariance matrix.

Denoising stabilizes the matrix making it making it more robust for optimization process that require inversion and offer best gains in tradictional optimization based on matrix inversion.

Clustering proportionate best diversification.

References

de Prado, Marcos López. Machine Learning for Asset Managers. Cambridge University Press, 2020.

de Prado, Marcos López, e Michael Lewis. "Detection of False Investment Strategies Using Unsupervised Learning Methods." SSRN, 18 Aug. 2018.

Ledoit, Olivier, e Michael Wolf. "Honey, I Shrunk the Sample Covariance Matrix." The Journal of Portfolio Management, 2003.

Marti, Gautier, et al. "A Review of Two Decades of Correlations, Hierarchies, Networks and Clustering in Financial Markets." Progress in Information Geometry, 2017, pp. 245-274.

Sato, Yoshiharu. "Model-Free Reinforcement Learning for Financial Portfolios: A Brief Survey." arXiv, 2019, pp. 1-20.

Zhang, Wensheng, et al. "Exact Distributions of Finite Random Matrices and Their Applications to Spectrum Sensing." Sensors, vol. 16, no. 8, 2016.