在统计学中，估计数据集的协方差矩阵是一个重要任务。最大似然估计（MLE）虽然常用，但在某些情况下可能会产生不稳定的结果。为了解决这个问题，可以采用正则化技术，其中收缩（shrinkage）是一种流行的正则化方法。Ledoit和Wolf提出了一种计算渐近最优收缩参数的方法，该参数通过最小化均方误差（MSE）准则来确定，从而得到Ledoit-Wolf协方差估计。Chen等人进一步提出了OAS系数，这是对Ledoit-Wolf收缩参数的改进，其在假设数据为高斯分布的情况下收敛性能更佳。
本示例受到Chen等人发表的论文[1]的启发，展示了使用高斯分布数据估计LW和OAS方法的均方误差（MSE）的比较。文中提到的两种方法，Ledoit-Wolf和OAS，都旨在通过收缩技术来改善协方差矩阵的估计。收缩技术的核心思想是将样本协方差矩阵向一个目标矩阵收缩，以减少估计的不确定性。Ledoit-Wolf方法通过一个封闭公式计算收缩参数，而OAS方法则在Ledoit-Wolf的基础上进行了改进，以期获得更好的收敛性能。
为了比较这两种方法，生成了一系列不同样本大小的高斯分布数据，并分别使用Ledoit-Wolf和OAS方法估计协方差矩阵。通过计算每种方法的均方误差，可以评估它们在不同样本大小下的性能。实验结果表明，OAS方法在样本数量较少时，相较于Ledoit-Wolf方法，能够提供更准确的协方差估计。
以下是实现这一比较的Python代码示例。代码中首先设置了随机数生成器的种子，以确保结果的可重复性。然后定义了特征数量，并构建了一个AR(1)过程的模拟协方差矩阵。接着，为不同的样本数量范围和重复次数准备了数据结构，用于存储每种方法的均方误差和收缩参数。在循环中，生成了正态分布的随机样本，并分别使用Ledoit-Wolf和OAS方法进行拟合，计算它们的均方误差和收缩参数。最后，使用matplotlib库绘制了均方误差和收缩系数的图表，以直观地展示两种方法的性能差异。
import matplotlib.pyplot as plt
import numpy as np
from scipy.linalg import cholesky, toeplitz
from sklearn.covariance import OAS, LedoitWolf
np.random.seed(0)
n_features = 100 # 特征数量
r = 0.1 # AR(1)过程的参数
real_cov = toeplitz(r ** np.arange(n_features)) # 模拟协方差矩阵
coloring_matrix = cholesky(real_cov) # Cholesky分解
n_samples_range = np.arange(6, 31) # 不同的样本数量
repeat = 100 # 重复次数
lw_mse = np.zeros((len(n_samples_range), repeat)) # Ledoit-Wolf的MSE
oa_mse = np.zeros((len(n_samples_range), repeat)) # OAS的MSE
lw_shrinkage = np.zeros((len(n_samples_range), repeat)) # Ledoit-Wolf的收缩参数
oa_shrinkage = np.zeros((len(n_samples_range), repeat)) # OAS的收缩参数
for i, n_samples in enumerate(n_samples_range):
for j in range(repeat):
X = np.dot(np.random.normal(size=(n_samples, n_features)), coloring_matrix.T)
lw = LedoitWolf(store_precision=False, assume_centered=True)
lw.fit(X)
lw_mse[i, j] = lw.error_norm(real_cov, scaling=False)
lw_shrinkage[i, j] = lw.shrinkage_
oa = OAS(store_precision=False, assume_centered=True)
oa.fit(X)
oa_mse[i, j] = oa.error_norm(real_cov, scaling=False)
oa_shrinkage[i, j] = oa.shrinkage_
# 绘制MSE图表
plt.subplot(2, 1, 1)
plt.errorbar(n_samples_range, lw_mse.mean(1), yerr=lw_mse.std(1), label="Ledoit-Wolf", color="navy", lw=2)
plt.errorbar(n_samples_range, oa_mse.mean(1), yerr=oa_mse.std(1), label="OAS", color="darkorange", lw=2)
plt.ylabel("Squared error")
plt.legend(loc="upper right")
plt.title("协方差估计器比较")
plt.xlim(5, 31)
# 绘制收缩系数图表
plt.subplot(2, 1, 2)
plt.errorbar(n_samples_range, lw_shrinkage.mean(1), yerr=lw_shrinkage.std(1), label="Ledoit-Wolf", color="navy", lw=2)
plt.errorbar(n_samples_range, oa_shrinkage.mean(1), yerr=oa_shrinkage.std(1), label="OAS", color="darkorange", lw=2)
plt.xlabel("样本数量")
plt.ylabel("收缩系数")
plt.legend(loc="lower right")
plt.ylim(plt.ylim()[0], 1.0 + (plt.ylim()[1] - plt.ylim()[0]) / 10.0)
plt.xlim(5, 31)
plt.show()
通过上述代码，可以清楚地看到，在不同样本数量下，Ledoit-Wolf和OAS方法的均方误差和收缩系数的变化情况。这些图表为提供了直观的比较，帮助理解不同协方差估计方法在实际应用中的表现。