首页 文章

乘以稀疏矩阵的列元素

提问于
浏览
2

我有一个稀疏的csc矩阵,其中包含许多零元素,我想为每行计算所有列元素的乘积 .

即:

A = [[1,2,0,0],
      [2,0,3,0]]

应转换为:

V = [[2,
      6]]

使用numpy密集矩阵,可以通过用一个值替换所有零值并使用 A.prod(1) 来实现 . 然而,这不是一种选择,因为密集矩阵太大 .

有没有办法在不将稀疏矩阵转换为密集矩阵的情况下实现此目的?

3 回答

  • 1

    Approach #1: 我们可以使用稀疏元素的行索引作为ID,并使用np.multiply.reduceat执行这些元素的相应值的乘法,以获得所需的输出 .

    因此,实施将是 -

    from scipy import sparse
    from scipy.sparse import csc_matrix
    
    r,c,v = sparse.find(a) # a is input sparse matrix
    out = np.zeros(a.shape[0],dtype=a.dtype)
    unqr, shift_idx = np.unique(r,return_index=1)
    out[unqr] = np.multiply.reduceat(v, shift_idx)
    

    样品运行 -

    In [89]: # Let's create a sample csc_matrix
        ...: A = np.array([[-1,2,0,0],[0,0,0,0],[2,0,3,0],[4,5,6,0],[1,9,0,2]])
        ...: a = csc_matrix(A)
        ...: 
    
    In [90]: a
    Out[90]: 
    <5x4 sparse matrix of type '<type 'numpy.int64'>'
        with 10 stored elements in Compressed Sparse Column format>
    
    In [91]: a.toarray()
    Out[91]: 
    array([[-1,  2,  0,  0],
           [ 0,  0,  0,  0],
           [ 2,  0,  3,  0],
           [ 4,  5,  6,  0],
           [ 1,  9,  0,  2]])
    
    In [92]: out
    Out[92]: array([ -2,   0,   6, 120,   0,  18])
    

    Approach #2: 我们正在执行基于bin的乘法 . 我们有np.bincount的基于bin的求和解决方案 . 因此,这里可以使用的技巧是将数字转换为对数,执行基于bin的求和,然后使用 exponential (日志的反向)转换回原始格式,并且's it! For negative numbers, we might to add a step or more, but let'看到非实现的情况负数 -

    r,c,v = sparse.find(a)
    out = np.exp(np.bincount(r,np.log(v),minlength = a.shape[0]))
    out[np.setdiff1d(np.arange(a.shape[0]),r)] = 0
    

    运行非负数的样本 -

    In [118]: a.toarray()
    Out[118]: 
    array([[1, 2, 0, 0],
           [0, 0, 0, 0],
           [2, 0, 3, 0],
           [4, 5, 6, 0],
           [1, 9, 0, 2]])
    
    In [120]: out  # Using listed code
    Out[120]: array([   2.,    0.,    6.,  120.,   18.])
    
  • 2

    做一个样本:

    In [51]: A=np.array([[1,2,0,0],[0,0,0,0],[2,0,3,0]])
    In [52]: M=sparse.csr_matrix(A)
    

    lil 格式中,每行的值存储在列表中 .

    In [56]: Ml=M.tolil()
    In [57]: Ml.data
    Out[57]: array([[1, 2], [], [2, 3]], dtype=object)
    

    拿这些产品:

    In [58]: np.array([np.prod(i) for i in Ml.data])
    Out[58]: array([ 2.,  1.,  6.])
    

    csr 格式中,值存储为:

    In [53]: M.data
    Out[53]: array([1, 2, 2, 3], dtype=int32)
    In [54]: M.indices
    Out[54]: array([0, 1, 0, 2], dtype=int32)
    In [55]: M.indptr
    Out[55]: array([0, 2, 2, 4], dtype=int32)
    

    indptr 给出行值的开头 . csr (和 csc )矩阵上的计算代码通常执行这样的计算(但已编译):

    In [94]: lst=[]; i=M.indptr[0]
    In [95]: for j in M.indptr[1:]:
        ...:     lst.append(np.product(M.data[i:j]))
        ...:     i = j    
    In [96]: lst
    Out[96]: [2, 1, 6]
    

    使用Diavaker的测试矩阵:

    In [137]: M.A
    Out[137]: 
    array([[-1,  2,  0,  0],
           [ 0,  0,  0,  0],
           [ 2,  0,  3,  0],
           [ 4,  5,  6,  0],
           [ 1,  9,  0,  2]], dtype=int32)
    

    上面的循环产生:

    In [138]: foo(M)
    Out[138]: [-2, 1, 6, 120, 18]
    

    Divakar的代码 uniquereduceat

    In [139]: divk(M)
    Out[139]: array([ -2,   0,   6, 120,  18], dtype=int32)
    

    (空行的不同值) .

    使用 indptr 简化Reduceat只是:

    In [140]: np.multiply.reduceat(M.data,M.indptr[:-1])
    Out[140]: array([ -2,   2,   6, 120,  18], dtype=int32)
    

    需要修复空第2行的值( indptr 值为[2,2,...], reduceat 使用 M.data[2] ) .

    def wptr(M, empty_val=1):
        res = np.multiply.reduceat(M.data, M.indptr[:-1])
        mask = np.diff(M.indptr)==0
        res[mask] = empty_val
        return res
    

    有一个更大的矩阵

    Mb=sparse.random(1000,1000,.1,format='csr')
    

    这个 wptr 比Divaker的版本快30倍 .

    关于跨稀疏矩阵的行计算值的更多讨论:Scipy.sparse.csr_matrix: How to get top ten values and indices?

  • 0

    您可以使用numpy模块中的prod()方法计算A的每个子列表中所有元素的乘积,同时排除值0的元素 .

    import numpy as np
    print [[np.prod([x for x in A[i] if x!=0 ]) for i in range(len(A))]]
    

相关问题