按列相关性的函数？

Question

我有两个numpy数组

X.shape = (100, 10)
Y.shape = (100, 10)

我想找出X和Y的列之间的皮尔逊相关性

也就是说

from scipy.stats.stats import pearsonr

def corr( X, Y ):
    return np.array([ pearsonr( x, y )[0] for x,y in zip( X.T, Y.T ) ] )    

corr( X, Y ).shape = (10, )

有没有什么函数可以做到这一点？到目前为止，我找到的所有函数都是计算相关矩阵的。在Matlab里有一个成对相关的函数，所以我相信一定有人为Python写了类似的函数。

我不喜欢上面那个示例函数的原因是它似乎很慢。

Answer 1

我对 scipy.stats.pearsonr 进行了修改：

from scipy.stats import pearsonr

x = np.random.rand(100, 10)
y = np.random.rand(100, 10)

def corr( X, Y ):
    return np.array([ pearsonr( x, y )[0] for x,y in zip( X.T, Y.T) ] )

def pair_pearsonr(x, y, axis=0):
    mx = np.mean(x, axis=axis, keepdims=True)
    my = np.mean(y, axis=axis, keepdims=True)
    xm, ym = x-mx, y-my
    r_num = np.add.reduce(xm * ym, axis=axis)
    r_den = np.sqrt((xm*xm).sum(axis=axis) * (ym*ym).sum(axis=axis))
    r = r_num / r_den
    return r

np.allclose(pair_pearsonr(x, y, axis=0), corr(x, y))

Answer 2

如果把列看作变量，把行看作在 X 和 Y 中的观测值（也就是说你想要找出 X 和 Y 之间每一列的相关性）：

X = (X - X.mean(axis=0)) / X.std(axis=0)
Y = (Y - Y.mean(axis=0)) / Y.std(axis=0)
pearson_r = np.dot(X.T, Y) / X.shape[0]

要找到 p 值，可以把 pearson_r 转换成 t 统计量：

t = pearson_r * np.sqrt(X.shape[0] - 2) / np.sqrt(1 - pearson_r ** 2)

然后 p 值就是 2 × P(T > t)。