有谁知道如何从 python 中非常大的稀疏矩阵计算相关矩阵？基本上，我正在寻找可以numpy.corrcoef在 scipy 稀疏矩阵上工作的东西。
问问题
9982 次
4 回答
12
您可以从协方差矩阵中相当直接地计算相关系数，如下所示：
import numpy as np
from scipy import sparse
def sparse_corrcoef(A, B=None):
if B is not None:
A = sparse.vstack((A, B), format='csr')
A = A.astype(np.float64)
n = A.shape[1]
# Compute the covariance matrix
rowsum = A.sum(1)
centering = rowsum.dot(rowsum.T.conjugate()) / n
C = (A.dot(A.T.conjugate()) - centering) / (n - 1)
# The correlation coefficients are given by
# C_{i,j} / sqrt(C_{i} * C_{j})
d = np.diag(C)
coeffs = C / np.sqrt(np.outer(d, d))
return coeffs
检查它是否工作正常：
# some smallish sparse random matrices
a = sparse.rand(100, 100000, density=0.1, format='csr')
b = sparse.rand(100, 100000, density=0.1, format='csr')
coeffs1 = sparse_corrcoef(a, b)
coeffs2 = np.corrcoef(a.todense(), b.todense())
print(np.allclose(coeffs1, coeffs2))
# True
被警告：
计算协方差矩阵所需的内存量C将在很大程度上取决于A（和B，如果给定）的稀疏结构。例如，如果A是一个(m, n)仅包含一列非零值的矩阵，那么C将是一个(n, n)包含所有非零值的矩阵。如果n很大，那么就内存消耗而言，这可能是一个非常坏的消息。
于 2013-10-07T22:31:43.920 回答
10
你不需要引入一个大的密集矩阵。只需使用 Numpy 保持稀疏：
import numpy as np
def sparse_corr(A):
N = A.shape[0]
C=((A.T*A -(sum(A).T*sum(A)/N))/(N-1)).todense()
V=np.sqrt(np.mat(np.diag(C)).T*np.mat(np.diag(C)))
COR = np.divide(C,V+1e-119)
return COR
测试性能：
A = sparse.rand(1000000, 100, density=0.1, format='csr')
sparse_corr(A)
于 2016-08-02T18:12:42.170 回答
1
我提出了一个并行运行的 scipy 稀疏矩阵的答案。这不是返回一个巨大的相关矩阵，而是在检查所有字段的正和负 Pearson 相关后返回要保留的字段特征掩码。
我还尝试使用以下策略最小化计算：
- 处理每一列
- 从当前列 + 1 开始并计算向右移动的相关性。
- 对于任何 abs(correlation) >= 阈值，将当前列标记为要删除，并且不再计算相关性。
- 对数据集中除最后一列之外的每一列执行这些步骤。
这可以通过保留标记为要删除的列的全局列表并跳过对此类列的进一步相关性计算来进一步加快速度，因为列将无序执行。但是，我对 python 中的竞争条件知之甚少，无法在今晚实现这一点。
返回列掩码显然允许代码处理比返回整个相关矩阵更大的数据集。
使用此函数检查每一列：
def get_corr_row(idx_num, sp_mat, thresh):
# slice the column at idx_num
cols = sp_mat.shape[1]
x = sp_mat[:,idx_num].toarray().ravel()
start = idx_num + 1
# Now slice each column to the right of idx_num
for i in range(start, cols):
y = sp_mat[:,i].toarray().ravel()
# Check the pearson correlation
corr, pVal = pearsonr(x,y)
# Pearson ranges from -1 to 1.
# We check both positive and negative correlations >= thresh using abs(corr)
if abs(corr) >= thresh:
# stop checking after finding the 1st correlation > thresh
return False
# Mark column at idx_num for removal in the mask
return True
并行运行列级相关性检查：
from joblib import Parallel, delayed
import multiprocessing
def Get_Corr_Mask(sp_mat, thresh, n_jobs=-1):
# we must make sure the matrix is in csc format
# before we start doing all these column slices!
sp_mat = sp_mat.tocsc()
cols = sp_mat.shape[1]
if n_jobs == -1:
# Process the work on all available CPU cores
num_cores = multiprocessing.cpu_count()
else:
# Process the work on the specified number of CPU cores
num_cores = n_jobs
# Return a mask of all columns to keep by calling get_corr_row()
# once for each column in the matrix
return Parallel(n_jobs=num_cores, verbose=5)(delayed(get_corr_row)(i, sp_mat, thresh)for i in range(cols))
一般用法：
#Get the mask using your sparse matrix and threshold.
corr_mask = Get_Corr_Mask(X_t_fpr, 0.95)
# Remove features that are >= 95% correlated
X_t_fpr_corr = X_t_fpr[:,corr_mask]
于 2020-12-17T04:18:18.697 回答
0
不幸的是，Alt 的回答对我没有用。赋予np.sqrt函数的值大多为负数，因此得到的协方差值为 nan。
我也无法使用ali_m 的答案，因为我的矩阵太大，无法将centering = rowsum.dot(rowsum.T.conjugate()) / n矩阵放入记忆中（我的矩阵尺寸为：3.5*10^6 x 33）
相反，我使用scikit-learnStandardScaler来计算标准稀疏矩阵，然后使用乘法来获得相关矩阵。
from sklearn.preprocessing import StandardScaler
def compute_sparse_correlation_matrix(A):
scaler = StandardScaler(with_mean=False)
scaled_A = scaler.fit_transform(A) # Assuming A is a CSR or CSC matrix
corr_matrix = (1/scaled_A.shape[0]) * (scaled_A.T @ scaled_A)
return corr_matrix
我相信这种方法比其他提到的方法更快、更健壮。此外，它还保留了输入矩阵的稀疏模式。
于 2019-12-30T15:03:06.153 回答