2017-07-28 76 views
5

我想遍历一个CSR矩阵的行和列的总和除以每个元素,类似这样的位置:大NumPy的SciPy的CSR矩阵,行明智的操作

numpy divide row by row sum

我的问题是我正在处理一个大矩阵:(96582,350138)

并且当从链接的帖子应用该操作时,由于返回的矩阵是密集的,所以它扩大了我的记忆。

所以这是我第一次尝试:

for row in counts: 
    row = row/row.sum() 

不幸的是,这并不影响基质可言,所以我想出了第二个想法,以创建一个新的CSR矩阵和CONCAT行使用vstack:

from scipy import sparse 
import time 

start_time = curr_time = time.time() 
mtx = sparse.csr_matrix((0, counts.shape[1])) 
for i, row in enumerate(counts): 
    prob_row = row/row.sum() 
    mtx = sparse.vstack([mtx, prob_row]) 
    if i % 1000 == 0: 
     delta_time = time.time() - curr_time 
     total_time = time.time() - start_time 
     curr_time = time.time() 
     print('step: %i, total time: %i, delta_time: %i' % (i, total_time, delta_time)) 

这种运作良好,但一些迭代之后它变得越来越慢:

step: 0, total time: 0, delta_time: 0 
step: 1000, total time: 1, delta_time: 1 
step: 2000, total time: 5, delta_time: 4 
step: 3000, total time: 12, delta_time: 6 
step: 4000, total time: 23, delta_time: 11 
step: 5000, total time: 38, delta_time: 14 
step: 6000, total time: 55, delta_time: 17 
step: 7000, total time: 88, delta_time: 32 
step: 8000, total time: 136, delta_time: 47 
step: 9000, total time: 190, delta_time: 53 
step: 10000, total time: 250, delta_time: 59 
step: 11000, total time: 315, delta_time: 65 
step: 12000, total time: 386, delta_time: 70 
step: 13000, total time: 462, delta_time: 76 
step: 14000, total time: 543, delta_time: 81 
step: 15000, total time: 630, delta_time: 86 
step: 16000, total time: 722, delta_time: 92 
step: 17000, total time: 820, delta_time: 97 

有何建议?任何想法为什么vstack变得越来越慢?

+1

见https://stackoverflow.com/a/45339754和https://stackoverflow.com/q/44080315 – hpaulj

+0

与密集数组一样,循环中的重复级联很慢。在列表中累积结果并执行'vstack'更快。 – hpaulj

回答

4

vstackO(n)操作,因为它需要为结果分配内存,然后将作为参数传递的所有数组的内容复制到结果数组中。

您可以简单地使用multiply的运做:

>>> res = counts.multiply(1/counts.sum(1)) # multiply with inverse 
>>> res.todense() 
matrix([[ 0.33333333, 0.  , 0.66666667], 
     [ 0.  , 0.  , 1.  ], 
     [ 0.26666667, 0.33333333, 0.4  ]]) 

但它也很容易使用np.lib.stride_tricks.as_strided做你想要的操作(相对高性能的)。这个as_strided函数还允许在数组上执行更复杂的操作(如果没有方法或函数)。使用scipy documentation的例子CSR

例如:

>>> from scipy.sparse import csr_matrix 
>>> import numpy as np 
>>> row = np.array([0,0,1,2,2,2]) 
>>> col = np.array([0,2,2,0,1,2]) 
>>> data = np.array([1.,2,3,4,5,6]) 
>>> counts = csr_matrix((data,(row,col)), shape=(3,3)) 
>>> counts.todense() 
matrix([[ 1., 0., 2.], 
     [ 0., 0., 3.], 
     [ 4., 5., 6.]]) 

您可以通过将其分割每行的总和是这样的:

>>> row_start_stop = np.lib.stride_tricks.as_strided(counts.indptr, 
                shape=(counts.shape[0], 2), 
                strides=2*counts.indptr.strides) 
>>> for start, stop in row_start_stop: 
... row = counts.data[start:stop] 
... row /= row.sum() 
>>> counts.todense() 
matrix([[ 0.33333333, 0.  , 0.66666667], 
     [ 0.  , 0.  , 1.  ], 
     [ 0.26666667, 0.33333333, 0.4  ]]) 
+1

很好的答案,您更新行的方式比我的效率高得多。像一百次!这应该是被接受的答案。 – Anis

2

编辑

@MSeifert答案更有效率,这应该是正确的做法。我认为写作counts[i, :]意味着一些列切片完成,我没有意识到。该文档明确指出这些对csr_matrix的操作确实是低效的。方式确实是一个很好的例子。

回答

的医生说排切片是有效的,我认为你应该做

for i in range(counts.shape[0]): 
    counts[i,:] /= counts[i,:].sum() 

这样,你在地方编辑矩阵,它保持稀疏,你不必使用vstack 。我不知道这是最高效的操作,但至少你不应该有内存问题,并没有因为被计算的行没有减速效果:

import time() 
s = time.time() 
for i in range(counts.shape[0]): 
    counts[i, :] /= (counts[i, :].sum() + 1) 
    if i % 1000 == 0: 
     e = time.time() 
     if i > 0: 
      print i, e-s 
     s = time.time() 

性能:

1000 6.00199794769 2000 6.02894091606 3000 7.44459486008 4000 7.10011601448 5000 6.16998195648 6000 7.79510307312 7000 7.00139117241 8000 7.37821507454 9000 7.28075814247 ...

性能为MSeifert的回答是:

row_start_stop = np.lib.stride_tricks.as_strided(counts.indptr, shape=(counts.shape[0], 2), 
               strides=2*counts.indptr.strides) 
for i, (start, stop) in enumerate(row_start_stop): 
    row = counts.data[start:stop] 
    row /= row.sum() 
    if i % 1000 == 0: 
     e = time.time() 
     if i > 0: 
      print i,e-s 
     s = time.time() 

1000 0.00735783576965 2000 0.0108380317688 3000 0.0102109909058 4000 0.0131571292877 5000 0.00670218467712 6000 0.00608897209167 7000 0.00663685798645 8000 0.0164499282837 9000 0.0061981678009 ... 至于为什么使用至很慢,@ MSeifert答案很棒。

+0

只是出于好奇,'为项目计数:item/= item.sum()'票价?我已经确信,接受的答案应该是Mseifert's,这是对隐式和显式切片之间性能差异的个人好奇心 – Uvar

+0

问题在于OP的问题。当在计数中为'item'做'''时,'''item'''实际上不是对实际计数行的引用。这就是为什么修改不会影响矩阵计数。所以我可以衡量表演,但我真的不明白这一点。 – Anis

+0

也许我对稀疏矩阵的理解只是不了解而已。 :/无论如何,正如你在我的答案中可以看到的那样,对于密集矩阵,你可以用这种方式对实际行进行引用。以为它会以相同的方式工作..:$ – Uvar

相关问题