gc5*_*gc5 3 python parallel-processing scipy dataframe pandas
我想将斯皮尔曼相关性应用于具有相同列数的两个熊猫数据帧(每对行的相关性)。
我的目标是计算每对行 (r, s) 之间的斯皮尔曼相关性分布,其中 r 是第一个数据帧中的一行,s 是第二个数据帧中的一行。
我知道以前已经回答过类似的问题(请参阅此)。但是,这个问题有所不同,因为我想将第一个数据帧的每一行与第二个数据帧中的所有行进行比较。此外,由于数据大小,这需要大量计算,并且需要数小时的时间。我想对其进行并行化,并可能重写它以加快速度。
我尝试使用 numba 但不幸的是它失败了(与此类似的问题),因为它似乎无法识别scipy Spearmanr。我的代码如下:
def corr(a, b):
dist = []
for i in range(a.shape[0]):
for j in range(b.shape[0]):
dist += [spearmanr(a.iloc[i, :], b.iloc[j, :])[0]]
return dist
Run Code Online (Sandbox Code Playgroud)
from numba import njit
import pandas as pd
import numpy as np
@njit
def mean1(a):
n = len(a)
b = np.empty(n)
for i in range(n):
b[i] = a[i].mean()
return b
@njit
def std1(a):
n = len(a)
b = np.empty(n)
for i in range(n):
b[i] = a[i].std()
return b
Run Code Online (Sandbox Code Playgroud)
@njit
def c(a, b):
''' Correlation '''
n, k = a.shape
m, k = b.shape
mu_a = mean1(a)
mu_b = mean1(b)
sig_a = std1(a)
sig_b = std1(b)
out = np.empty((n, m))
for i in range(n):
for j in range(m):
out[i, j] = (a[i] - mu_a[i]) @ (b[j] - mu_b[j]) / k / sig_a[i] / sig_b[j]
return out
Run Code Online (Sandbox Code Playgroud)
r = df_test.rank(1).values
df_test.T.corr('spearman') == c(r, r)
Run Code Online (Sandbox Code Playgroud)
进行斯皮尔曼排名相关只是进行排名相关。
我们可以利用杠杆argsort
来获得排名。尽管 的argsort
确实argsort
为我们提供了排名,但我们可以通过切片分配将自己限制为一种排序。
def rank(a):
i, j = np.meshgrid(*map(np.arange, a.shape), indexing='ij')
s = a.argsort(1)
out = np.empty_like(s)
out[i, s] = j
return out
Run Code Online (Sandbox Code Playgroud)
在关联排名的情况下,均值和标准差均由数组第二维的大小预先确定。
你可以在没有 numba 的情况下完成同样的事情,但我假设你想要它。
from numba import njit
@njit
def c(a, b):
n, k = a.shape
m, k = b.shape
mu = (k - 1) / 2
sig = ((k - 1) * (k + 1) / 12) ** .5
out = np.empty((n, m))
a = a - mu
b = b - mu
for i in range(n):
for j in range(m):
out[i, j] = a[i] @ b[j] / k / sig ** 2
return out
Run Code Online (Sandbox Code Playgroud)
对于后代来说,我们可以完全避免内部循环,但这可能会产生内存问题。
@njit
def c1(a, b):
n, k = a.shape
m, k = b.shape
mu = (k - 1) / 2
sig = ((k - 1) * (k + 1) / 12) ** .5
a = a - mu
b = b - mu
return a @ b.T / k / sig ** 2
Run Code Online (Sandbox Code Playgroud)
np.random.seed([3, 1415])
a = np.random.randn(2, 10)
b = np.random.randn(2, 10)
rank_a = rank(a)
rank_b = rank(b)
c(rank_a, rank_b)
array([[0.32121212, 0.01818182],
[0.13939394, 0.55151515]])
Run Code Online (Sandbox Code Playgroud)
如果您正在与DataFrame
da = pd.DataFrame(a)
db = pd.DataFrame(b)
pd.DataFrame(c(rank(da.values), rank(db.values)), da.index, db.index)
0 1
0 0.321212 0.018182
1 0.139394 0.551515
Run Code Online (Sandbox Code Playgroud)
我们可以使用进行快速验证pandas.DataFrame.corr
pd.DataFrame(a.T).corr('spearman') == c(rank_a, rank_a)
0 1
0 True True
1 True True
Run Code Online (Sandbox Code Playgroud)