如何高效计算两个二维Numpy数组对应列的皮尔逊相关系数？

针对你需要计算两个同维度数组对应列皮尔逊相关系数的需求，直接基于皮尔逊系数的数学公式实现向量化计算是最优方案——既避免了循环的低效，也不会像np.corrcoef那样计算大量多余的相关系数，同时大幅降低内存占用。

实现思路 #

皮尔逊相关系数的核心公式为：
$$r = \frac{\text{cov}(x,y)}{\sigma_x \sigma_y}$$
其中：

• $\text{cov}(x,y)$ 是两列的协方差（无偏估计，除以样本量减1）
• $\sigma_x$、$\sigma_y$ 分别是两列的标准差（无偏估计，除以样本量减1）

我们可以将整个计算过程拆解为向量化操作，直接对所有列批量处理：

代码实现 #

import numpy as np

# 假设A和B是形状为(18000, 18000)的numpy数组
n_samples = A.shape[0]

# 1. 计算每列的均值
mu_A = A.mean(axis=0)
mu_B = B.mean(axis=0)

# 2. 对每列进行中心化（减去列均值）
A_centered = A - mu_A[np.newaxis, :]
B_centered = B - mu_B[np.newaxis, :]

# 3. 计算对应列的协方差（无偏估计）
cov_xy = (A_centered * B_centered).sum(axis=0) / (n_samples - 1)

# 4. 计算每列的标准差（无偏估计）
std_A = A.std(axis=0, ddof=1)
std_B = B.std(axis=0, ddof=1)

# 5. 计算对应列的皮尔逊相关系数
corr_vec = cov_xy / (std_A * std_B)

为什么这个方法更高效 #

1. 无冗余计算：仅针对A和B的对应列计算所需的协方差和标准差，不会像np.corrcoef那样生成包含A内部、B内部所有列对的(36000, 36000)巨型矩阵（这个矩阵仅存储就需要数十GB内存）。
2. 向量化操作：所有步骤都是Numpy底层优化的向量/矩阵运算，比Python循环快几个数量级。
3. 内存友好：中心化后的数组虽然和原数组同尺寸，但无需存储超大的相关系数矩阵，内存压力远小于np.corrcoef的方案。

特殊情况处理 #

如果某列的标准差为0（即该列所有元素相同），std_A或std_B会出现0，此时对应位置的相关系数会变为NaN，这和np.corrcoef的行为一致，符合统计学定义。

内容的提问来源于stack exchange，提问作者Long-Hung Pham