在R语言中验证样本协方差矩阵计算方法的正确性

核心问题：定义差异 #

• 你看到的文献公式$$S = \sum_{j=1}^n(\mathbf{X}_j - \bar{\mathbf{X}})(\mathbf{X}_j-\bar{\mathbf{X}})'$$，是未做归一化的交叉乘积和矩阵，仅计算了离均差外积的总和，没有除以样本相关系数。
• R内置的cov()函数默认计算的是无偏样本协方差矩阵，即上述求和结果除以n-1（n是样本数量，这里n=4，所以除以3）。

你的手动代码错误 #

循环写的是for (i in 1:3)，但你的样本有4个（X是4行3列），少算了第4个样本的离均差外积项，导致求和结果偏小。

修正后的手动实现与验证 #

步骤1：正确计算文献公式的求和结果 #

x <- c(1, 3, 5, 2)
y <- c(2, 3, 8, 7)
z <- c(22, 1, 3, 3)
X <- cbind(x, y, z)

xbar <- apply(X, 2, mean)
d <- as.matrix(t(t(X) - xbar))

# 遍历全部4个样本，而非3个
s_cross <- matrix(0, 3, 3)
for (i in 1:nrow(X)) {
  s_cross <- s_cross + d[i, ] %*% t(d[i, ])
}
s_cross

运行结果：

x     y        z
[1,]   8.750  10.0 -33.750
[2,]  10.000  26.0 -53.000
[3,] -33.750 -53.0 292.750

步骤2：对比不同归一化结果 #

• 除以n-1（和cov()结果一致，无偏样本协方差）：

s_cross / (nrow(X)-1)

结果与cov(X)完全匹配：

x          y         z
x   2.916667   3.333333 -11.25000
y   3.333333   8.666667 -17.66667
z -11.250000 -17.666667  97.58333

• 除以n（有偏总体协方差估计）：

s_cross / nrow(X)

结果：

x      y        z
x   2.187500  2.500 -8.43750
y   2.500000  6.500 -13.25000
z -8.437500 -13.250 73.18750

总结 #

• 文献公式对应的是修正后手动代码的s_cross结果（未归一化的交叉乘积和）。
• R的cov()是无偏样本协方差，等于交叉乘积和除以n-1。
• 你原始手动代码因循环次数错误（少算一个样本），结果既不符合文献公式，也不符合协方差定义。

内容的提问来源于stack exchange，提问作者John Smith