了解csr格式

Question

我正在尝试理解scipy CSR是如何工作的。

https://docs.scipy.org/doc/scipy/reference/sparse.html

例如，https://en.wikipedia.org/wiki/Sparse_matrix上的以下矩阵

( 0 0 0 0 )
( 5 8 0 0 )
( 0 0 3 0 )
( 0 6 0 0 )

它说CSR表示如下所示。

V是否必须在行列表中从左到右列出具有非零元素的一行？

我可以理解COL_INDEX是对应于V中的元素的列索引(列1被索引为0)。

我不明白ROW_INDEX。谁能告诉我ROW_INDEX是如何从原始矩阵中创建的？谢谢。

   V         = [ 5 8 3 6 ]
   COL_INDEX = [ 0 1 2 1 ]
   ROW_INDEX = [ 0 0 2 3 4 ]

Answer 1

coo格式

我认为最好从coo定义开始。它更容易理解，并得到了广泛的应用：

In [90]: A = np.array([[0,0,0,0],[5,8,0,0],[0,0,3,0],[0,6,0,0]])                                 
In [91]: M = sparse.coo_matrix(A)                                                                

这些值存储在3个属性中：

In [92]: M.row                                                                                   
Out[92]: array([1, 1, 2, 3], dtype=int32)
In [93]: M.col                                                                                   
Out[93]: array([0, 1, 2, 1], dtype=int32)
In [94]: M.data                                                                                  
Out[94]: array([5, 8, 3, 6])

我们可以从这3个数组中生成一个新的矩阵：

In [95]: sparse.coo_matrix((_94, (_92, _93))).A                                                  
Out[95]: 
array([[0, 0, 0],
       [5, 8, 0],
       [0, 0, 3],
       [0, 6, 0]])

哦，我需要添加一个形状，因为有一列都是0：

In [96]: sparse.coo_matrix((_94, (_92, _93)), shape=(4,4)).A                                     
Out[96]: 
array([[0, 0, 0, 0],
       [5, 8, 0, 0],
       [0, 0, 3, 0],
       [0, 6, 0, 0]])

显示此矩阵的另一种方式：

In [97]: print(M)                                                                                
  (1, 0)    5
  (1, 1)    8
  (2, 2)    3
  (3, 1)    6

np.where(A)给出了相同的非零坐标。

In [108]: np.where(A)                                                                            
Out[108]: (array([1, 1, 2, 3]), array([0, 1, 2, 1]))

转换为csr

一旦我们有了coo，我们就可以很容易地将它转换成csr。事实上，sparse经常为我们做到这一点：

In [98]: Mr = M.tocsr()                                                                          
In [99]: Mr.data                                                                                 
Out[99]: array([5, 8, 3, 6], dtype=int64)
In [100]: Mr.indices                                                                             
Out[100]: array([0, 1, 2, 1], dtype=int32)
In [101]: Mr.indptr                                                                              
Out[101]: array([0, 0, 2, 3, 4], dtype=int32)

Sparse做了几件事--它对索引进行排序，对重复项求和，并用indptr数组替换row。在这里，它实际上比原始的更长，但通常它会更短，因为它每行只有一个值(加1)。也许更重要的是，大多数快速计算例程，特别是矩阵乘法，都是使用 csr 格式编写的。

我已经用过这个包很多次了。它的默认定义是coo样式，但内部存储是csc (但不像scipy那样向用户公开)。但是我从来没有尝试过从头开始派生indptr。我可以，但我不需要。

csr_matrix接受coo格式的输入，但也接受indptr etc格式的输入。我不推荐这样做，除非你已经计算了这些输入(比如从另一个矩阵)。它更容易出错，而且可能也不会快很多。

使用indptr进行迭代

但是，有时在intptr上迭代并直接在data上执行计算是很有用的。通常，这比使用提供的方法更快。

例如，我们可以按行列出非零值：

In [104]: for i in range(Mr.shape[0]): 
     ...:     pt = slice(Mr.indptr[i], Mr.indptr[i+1]) 
     ...:     print(i, Mr.indices[pt], Mr.data[pt]) 
     ...:                                                                                        
0 [] []
1 [0 1] [5 8]
2 [2] [3]
3 [1] [6]

保留初始0会使迭代变得更容易。当矩阵是(10000,90000)时，没有太多的动机将indptr的大小减少1。

lil格式

lil格式以类似的方式存储矩阵：

In [105]: Ml = M.tolil()                                                                         
In [106]: Ml.data                                                                                
Out[106]: array([list([]), list([5, 8]), list([3]), list([6])], dtype=object)
In [107]: Ml.rows                                                                                
Out[107]: array([list([]), list([0, 1]), list([2]), list([1])], dtype=object)

In [110]: for i,(r,d) in enumerate(zip(Ml.rows, Ml.data)): 
     ...:     print(i, r, d) 
     ...:                                                                                        
0 [] []
1 [0, 1] [5, 8]
2 [2] [3]
3 [1] [6]

由于行的存储方式，lil实际上允许我们获取view

In [167]: Ml.getrowview(2)                                                                       
Out[167]: 
<1x4 sparse matrix of type '<class 'numpy.longlong'>'
    with 1 stored elements in List of Lists format>
In [168]: for i in range(Ml.shape[0]): 
     ...:     print(Ml.getrowview(i)) 
     ...:                                                                                        

  (0, 0)    5
  (0, 1)    8
  (0, 2)    3
  (0, 1)    6

Answer 2

从scipy手册中：

数据索引((

，csr_matrix，indptr)，shape=(M，N))是标准的CSR表示，其中行i的列索引存储在索引[indptri:indptri+1]中，其相应的值存储在数据[indptri:indptri+1]中。如果未提供形状参数，则从索引数组推断矩阵尺寸。

indptr与ROW_INDEX相同，indicies与COL_INDEX相同。

下面是一个创建索引和值数组的简单方法的示例。本质上，ROW_INDICESi +1是从第0行到第i行的非零条目的总数，最后一个条目是非零条目的总数。

ROW_INDICES = [0]
COL_INDICES = []
VALS = []
for i in range(num_rows):
    ROW_INDICES.append(ROW_INDICES[i])
    for j in range(num_cols):
        if m[i, j] > 0:
            ROW_INDICES[i + 1] += 1
            COL_INDICES.append(j)
        VALS.append(m[i, j])