样本协方差的特征值远离协方差的特征值

Question

[我编写了以下代码来创建具有指定特征值的随机矩阵Sigma，然后从零均值和协方差Sigma的多元正态分布中采样矢量。

def generate_data(N, d):
    eigenvalues = [0] * (d + 1)
    for k in range(1, d + 2):
        eigenvalues[k - 1] = k

    random_matrix = numpy.random.randn(d + 1, d + 1)
    random_orthogonal = numpy.linalg.qr(random_matrix)[0]
    sqrt_cov = random_orthogonal @ numpy.diag(numpy.sqrt(eigenvalues))

    X = numpy.zeros((N, d + 1))
    for i in range(N):
        vec = numpy.random.standard_normal(d + 1)
        X[i] = sqrt_cov @ vec

在此代码之后，X应该是从所需分布中采样的d + 1矩阵为N。

现在我想知道X的样本协方差矩阵的特征值是什么。如果我没记错的话，它应该类似于Sigma

def get_sample_covariance(data):
    data = data - data.mean(axis=0)
    sample_cov = data.T @ data / (data.shape[0] - 1)
    return sample_cov

然后我绘制了sample_cov的特征值

我期望一个大致线性的函数，从d（等于500）到1。我知道了

什么给了？错误在哪里？

Answer 1

样本的生成似乎是正确的。但是，仅当N >> d时，来自N个随机样本的协方差估计才是正确的。特别地，最高特征值倾向于系统地偏离。但是，如果您从所有特征值中减去平均特征值，那将是非常准确的。

def get_max_mean_eig_run(N, d):
    """Return highest and mean eigenvalue for random samples"""

    X = generate_data(N, d)
    cov = get_sample_covariance(X)
    eig = np.linalg.eigvals(cov)
    return eig.max(), eig.mean()

plt.close('all')
np.random.seed(1)
ds = np.arange(1, 501, 25)

# shape (N, 2) for max and mean eigenvalues
eigs = np.array([get_max_mean_eig_run(1000, d) for d in ds])

plt.close('all')
plt.xlabel('d')
plt.ylabel('Eigenvalue')
plt.plot(ds, eigs[:, 0], 'r-', label='Max')
plt.plot(ds, eigs[:, 1], 'r--', label='Mean')
plt.scatter(ds, ds+1, color='k', marker='o', label='max expected')
plt.scatter(ds, (ds+2)/2, color='k', marker='*', label='mean expected')
plt.legend()

对于此特征值的特定输入频谱，您需要N > 100*d才能使最大特征值接近最大输入特征值（+ 5％），但在更实际的情况下可能会有所不同。

这里是特征值的直方图（对于N = 1000，d = 500）

plt.close('all')
X = generate_data(1000, 500)
eigs = np.linalg.eigvals(get_sample_covariance(X))
hist, bin_edges = np.histogram(eigs, bins=20)
bin_centers = (bin_edges[1:] + bin_edges[:-1])/2
plt.step(bin_centers, hist, where='mid')
plt.xlabel('Eigenvalue')
plt.ylabel('Count (per bin)')