计算巨大的稀疏对角相关矩阵的最有效方法是什么？答案

【问题标题】：What is the most efficient way to calculate a huge sparse diagonal correlation matrix?计算巨大的稀疏对角相关矩阵的最有效方法是什么？
【发布时间】：2016-03-09 10:17:15
【问题描述】：

我需要计算一个巨大的相关矩阵，例如 200000x200000，它太大而无法存储在内存中。幸运的是，大多数值都是 0，除了接近矩阵对角线的值，我需要计算。因此，我很好奇 scipy/numpy 中的稀疏矩阵是否可以帮助我加快速度。

我目前构造数据的方式如下。

#Input variables are snps, and max_dist
num_snps, num_indivs = snps.shape    
corr_table = {}
for i in range(num_snps):
    corr_table[i] = {}

for i in range(0, num_snps - 1):
    start_i = i + 1
    end_i = min(start_i + max_dist, num_snps)
    corr_vec = sp.dot(snps[i], sp.transpose(snps[start_i:end_i])) / float(num_indivs)
    corr_vec = sp.array(corr_vec).flatten()
    for k in range(start_i, end_i):
        corr_vec_i = k - start_i
        corr_table[i][k] = corr_vec[corr_vec_i]
        corr_table[k][i] = corr_vec[corr_vec_i]
return corr_table

这里 snps 是一个 MxN 矩阵，具有标准化的行向量（均值 0 和方差 1），我想为此计算 MxM 相关矩阵。目前，相关矩阵存储为一个巨大的字典（corr_table）。 max_dist 表示我计算相关性的一对 SNP（snps 矩阵中的行）之间的最大距离。对于所有其他相关性（不在 corr_table 中），我假设它们为 0。

不幸的是，这在实践中仍然不是很有效，因此，我想知道是否可以将矩阵乘法与稀疏矩阵一起使用来更有效地计算相关矩阵，而无需使用更多内存。

任何建议将不胜感激。

【问题讨论】：

标签： python scipy sparse-matrix correlation bigdata

【解决方案1】：

我没有尝试理解或运行您的计算，但我可以添加一些关于稀疏矩阵的指针。

包中有六个不同的稀疏formats，每个都有不同的优点和缺点。并且很容易将一种格式转换为另一种格式；稀疏函数经常这样做。没有一个非常适合增量更新。

dok 格式实际上是一个字典子类。键是坐标元组，例如corr_table[(i,k)] = corr_vec[corr_vec_i]。我在其他 SO 问题中发现，使用这些键构建一个普通字典会更快，然后是 update 和 dok。 corr_matrix[i,k]=... 索引的开销更大。

lol 也比较适合增量更改。它将值存储在 2 个列表列表中，矩阵的每行一个子列表。

csr 适用于矩阵计算，但对于索引分配很慢。它最好使用coo 输入样式构建，它使用3 个一维数组data、i、j。传统上，稀疏矩阵是通过构造这 3 个数组来构建的，如果以增量方式进行，可能作为列表，然后将它们传递给 coo_matrix。

还有块格式和对角格式，这对于一些稀疏布局可能更好。

我怀疑这一步：

for k in range(start_i, end_i):
    corr_vec_i = k - start_i
    corr_table[i][k] = corr_vec[corr_vec_i]
    corr_table[k][i] = corr_vec[corr_vec_i]

可以执行numpy向量操作，类似

jj = np.arange(start_i, end_i)
ii = jj - start_i
data_list.append(corr_vec[ii])
row_list.append(ii)
col_list.append(jj)
data_list.append(corr_vec[ii])
row_list.append(jj)
col_list.append(ii)

data_list 等是我为coo 收集输入的列表。它们可能必须通过np.concatenate 来生成sparse.coo_matrix 可以使用的一维数组。

我没有测试过这段代码，所以有错误，但希望它能给你一些想法。

【讨论】：