SciPy 稀疏矩阵（2）：COO

上回说到，计算机存储稀疏矩阵的核心思想就是对矩阵中的非零元素的信息进行一个必要的管理。然而，我们都知道在稀疏矩阵中零元素的分布通常情况下没有什么规律，因此仅仅存储非零元素的值是不够的，我们还需要非零元素的其他信息，具体需要什么信息很容易想到：考虑到在矩阵中的每一个元素不仅有值，同时对应的信息还有矩阵的行和列。因此，将非零元素的值外加上其对应的行和列构成一个三元组（行索引，列索引，值）。然后再按照某种规律存储这些三元组。

三元组的存储策略

如果存储一个稀疏矩阵对应的多个三元组可以有非常多的实现方式，针对每一种都进行讲解是非常不现实的，而且完全没有这个必要，因为三元组的存储策略可以分为 2 大类：三元组容器法以及三个序列法。

01

三元组容器法

我们首先看到三元组容器法，三元组容器法指的是把多个三元组放到一个容器中，按照这种方法可以定义出如下类：

class Element:
    def __init__(self, row, col, value):
        self.row, self.col, self.value = row, col, value


class SparseMatrix:
    def __init__(self, elements):
        self.elements = list(elements)

三元组在这里直接被我定义成了一个类，非常简单。至于稀疏矩阵就是多个三元组类的实例构成的一个容器，因此其属性初始化函数的参数就是多个三元组类的实例。需要注意的是我在属性初始化的时候使用 list 把多个三元组的实例转换成了序列，当然也可以转换成集合或者其他数据结构，做法不唯一。与此同时，针对稀疏矩阵类我们还可以添加一些功能，比如获取矩阵的行和列、多个三元组的行索引外加上列索引均重复该如何处理等等。

02

三个序列法

三个序列法就没有太多可供自由发挥的实现方式，它的实现方式非常的单一，先把多个三元组按照一定的顺序排好，然后依次读取三元组的行索引构成第 1 个序列（记作行索引序列），依次读取三元组的列索引构成第 2 个序列（记作列索引序列），依次读取三元组的值构成第 3 个序列（记作元素值序列），按照这种方法可以定义出如下类：

class Element:
    def __init__(self, row, col, value):
        self.row, self.col, self.value = row, col, value


class SparseMatrix:
    def __init__(self, elements):
        self.rows, self.cols, self.values = [], [], []
        for element in elements:
            self.rows.append(element.row)
            self.cols.append(element.col)
            self.values.append(element.value)

三个序列法实现起来就是这么简单。需要注意的是我采用 3 个列表来表示对应的 3 个序列，当然换成其他前驱和后继均只有一个的数据结构（比如链表）也是可以的。与此同时，针对稀疏矩阵类我们还可以添加一些功能，比如获取矩阵的行和列、多个三元组的行索引外加上列索引均重复该如何处理等等。

SciPy COO 格式的稀疏矩阵

在开始 SciPy COO 格式的稀疏矩阵之前我花了一些篇幅讲解稀疏矩阵的三元组存储策略，这主要是因为 SciPy COO 格式的稀疏矩阵用的存储策略就是三元组存储策略的第 2 种方法：三个序列法。和上述定义除了属性名有一点不同（意思是一样的），其他几乎没有什么区别。在 SciPy COO 格式的稀疏矩阵中，行索引序列的属性名就是 row，列索引序列的属性名就是 col，元素值序列的属性名就是 data。还有就是这 3 个序列并不是使用 Python 列表，而是 NumPy 数组。之所以这样的格式叫做 COO，是因为 COO 是英文 coordinate 的前 3 个字母，很明显这种存储格式只存储矩阵中非零元素的坐标和非零元素值。

01

实例化

SciPy COO 格式的稀疏矩阵类的定义位于 scipy.sparse 包中的 coo_matrix 类，对其进行实例化就能获取一个 SciPy COO 格式的稀疏矩阵的实例。当然，构造实例的方法主要有 4 种：

1. coo_matrix(D)：D 是一个普通矩阵（二维数组）。
2. coo_matrix(S)：S 是一个稀疏矩阵。
3. coo_matrix((M, N), [dtype])：会实例化一个 M 行 N 列元素类型为 dtype 的全 0 矩阵。dtype 是一个可选参数，默认值为双精度浮点数。
4. coo_matrix((data, (i, j)), [shape=(M, N)])：data、i 和 j 分别对应着元素值序列、行索引序列以及列索引序列。shape 参数如果没有被指定，则会通过行索引序列以及列索引序列进行推断。

02

案例

实例化一个 3 行 4 列元素类型为有符号 8 位整数的全 0 矩阵：

>>> from scipy import sparse
>>> import numpy as np
>>> mtx = sparse.coo_matrix((3, 4), dtype=np.int8)
>>> mtx.todense()
matrix([[0, 0, 0, 0],
        [0, 0, 0, 0],
        [0, 0, 0, 0]], dtype=int8)

通过元素值序列、行索引序列以及列索引序列来实例化一个稀疏矩阵：

>>> row = np.array([0, 3, 1, 0])
>>> col = np.array([0, 3, 1, 2])
>>> data = np.array([4, 5, 7, 9])
>>> mtx = sparse.coo_matrix((data, (row, col)), shape=(4, 4))
>>> mtx
<4x4 sparse matrix of type '<class 'numpy.int32'>'
        with 4 stored elements in COOrdinate format>
>>> mtx.todense()
matrix([[4, 0, 9, 0],
        [0, 7, 0, 0],
        [0, 0, 0, 0],
        [0, 0, 0, 5]])

对于上述实例化方法，可能存在一种比较极端的情况：非零元素的行列索引可能会重复多次。我们来看一下遇到这种情况会不会有什么问题：

>>> row = np.array([0, 0, 1, 3, 1, 0, 0])
>>> col = np.array([0, 2, 1, 3, 1, 0, 0])
>>> data = np.array([1, 1, 1, 1, 1, 1, 1])
>>> mtx = sparse.coo_matrix((data, (row, col)), shape=(4, 4))
>>> mtx.todense()
matrix([[3, 0, 1, 0],
        [0, 2, 0, 0],
        [0, 0, 0, 0],
        [0, 0, 0, 1]])

显然，这种情况是没有问题的，还有就是在转普通矩阵之后根据普通矩阵的元素可以看出它会把重复的行列索引对应的元素值做一个求和得到普通矩阵对应位置的元素。

同样是上述实例化方法，它还有一个极端情况：data 数组中有零元素。我们来看一下遇到这种情况会不会有什么问题：

>>> row = np.array([0, 0, 1, 3, 1, 0, 0])
>>> col = np.array([0, 2, 1, 3, 1, 0, 0])
>>> data = np.array([1, 1, 1, 1, 1, 1, 0])
>>> mtx = sparse.coo_matrix((data, (row, col)), shape=(4, 4))
>>> mtx
<4x4 sparse matrix of type '<class 'numpy.int32'>'
        with 7 stored elements in COOrdinate format>

显然，这种情况是没有问题的。需要注意的是 mtx 的 repr 表示说存储了 7 个元素，这是因为存储的元素数量和 data 的元素个数是一样的，不管行列索引有没有重复，不管元素是什么。我们都知道，一直带着这样的零元素或者重复的行列索引并不合理，如何消除这两者很简单，消除零元素可以通过调用 eliminate_zeros() 方法得以实现，消除重复的行列索引可以通过 sum_duplicates() 方法得以实现。这 2 个方法都是原地操作，无返回值。现在方法有了，怎么消除零元素以及重复的行列索引无非就是两个方法的调用顺序的问题。显然我们应该先消除重复的行列索引，再消除零元素。反过来之所以不行是因为可能存在重复 2 次的行列索引，一个地方元素值为 1，另一个地方元素值为 -1，显然它们都不是 0，所以先消除零元素不能把它们消去，然后消除重复的行列索引把它们加在一起又出现了零元素。

不支持随机存取：

>>> mtx[2, 3]
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: 'coo_matrix' object is not subscriptable

03

优缺点

SciPy COO 格式的稀疏矩阵有着以下优点：

1. 有利于各种稀疏矩阵格式的快速转换。
2. 允许重复的行列索引。
3. 可以高效地构造稀疏矩阵。
4. 在借助稀疏工具的情况下，可以高效地进行矩阵左乘列向量的操作。
5. 考虑到 data 属性是一个数组，所以针对逐元素操作（这里的逐元素操作要求零元素必须映射成零元素）效率很高，毕竟直接修改 data 属性即可。

当然，SciPy COO 格式的稀疏矩阵也有缺点：

1. 不支持元素访问以及切片访问。
2. 无法直接支持算术运算（进行算术运算之前通常会隐式地构造其他格式的稀疏矩阵，一般来说基本上都是 CSR 格式或者 CSC 格式）。

下回预告

COO 格式的稀疏矩阵因为只存储非零元素的信息，因此空间复杂度就是 O(k)，其中 k 表示非零元素的个数。当然，COO 格式的稀疏矩阵不支持元素访问是其中的一个不足之处，当然也没必要自己实现一个它的元素访问操作，因为在不改动 COO 属性定义的情况下我们实现的这一操作时间复杂度是 O(k)，毕竟要考虑到重复的行列索引，不能找到就停止查找！至于如何优化元素访问这一操作，继续使用这样的格式可能不好办了，需要从格式上进行改进。至于怎么去改进也不需要我们去实现，SciPy 已经实现了这样的稀疏矩阵格式，它就是下一回要介绍的 DOK 格式的稀疏矩阵！

bilibili 账号：新时代的运筹帷幄，喜欢的可以关注一下，看完视频不要忘记一键三连啊！

针对 SciPy 稀疏矩阵有比我这个更容易、门槛更低的学习路线的可以后台回复“加群”，备注：Python 机器学习算法说书人，不备注可是会被拒绝的哦~