Ctrl+K

从 spmatrix 迁移到 sparray#

本文档提供了关于如何将 scipy.sparse 中的稀疏矩阵代码转换为稀疏数组代码的指南。

从稀疏矩阵到稀疏数组的更改，与从 np.matrix 到 np.ndarray 的转换类似。本质上，我们必须从一个以全二维矩阵乘法为中心的 matrix 对象，转换为支持矩阵乘法运算符和逐元素计算的一维或二维“数组”对象。

符号：在本指南中，我们通常将稀疏数组类表示为 sparray，将稀疏矩阵类表示为 spmatrix。稠密 NumPy 数组表示为 np.ndarray，稠密矩阵类表示为 np.matrix。支持的稀疏格式表示为 BSR、COO、CSC、CSR、DIA、DOK、LIL，并且所有格式都由 sparray 和 spmatrix 支持。术语 sparse 指的是 sparray 或 spmatrix，而 dense 指的是 np.ndarray 或 np.matrix。

概述和总体情况#

  • 构造函数名称 *_matrix（例如，csr_matrix）更改为 *_array。

  • spmatrix M 始终是二维的（行 x 列），即使是 M.min(axis=0) 也是如此。sparray A 可以是一维或二维的。NumPy 标量将返回完整（0 维）缩减，例如 M.min()。

  • 迭代 sparray 会产生一维 sparray。迭代 spmatrix 会产生二维行 spmatrix

  • 更改行为的运算符有：*、@、*=、@=、**

    • 标量乘法，例如 5 * A，使用 *，而 5 @ A 未实现。

    • sparray 使用 * 进行逐元素乘法，使用 @ 进行矩阵乘法，而 spmatrix 使用运算符 * 或 @ 进行矩阵乘法。两者都可以使用 A.multiply(B) 进行逐元素乘法。

    • 标量指数，例如 A**2，对于 sparray 使用逐元素幂，对于 spmatrix 使用矩阵幂。sparray 的矩阵幂使用 scipy.sparse.linalg.matrix_power(A, n)。

  • 当向构造函数提供索引数组时，spmatrix 会根据传入数组的 dtype 和值选择 dtype，而 sparray 仅根据传入数组的 dtype 选择。例如，M=csr_matrix((data, indices, indptr)) 会导致 M.indices 的 dtype 为 int32，只要 indices 和 indptr 中的值较小，即使传入数组的 dtype 为 int64 也是如此。相反，当输入数组为 int64 时，A=csr_array((data, indices, indptr)) 会导致 A.indices 的 dtype 为 int64。这为 sparray 提供了更可预测的、通常更大的索引 dtype，并减少了为匹配 dtype 而进行的强制转换。

  • 检查稀疏类型和格式

    • issparse(A) 对于任何稀疏数组/矩阵返回 True。

    • isspmatrix(M) 对于任何稀疏矩阵返回 True。

    • isspmatrix_csr(M) 检查具有特定格式的稀疏矩阵。应将其替换为与数组兼容的版本，例如

    • issparse(A) and A.format == 'csr'，它检查 CSR 稀疏数组/矩阵。

  • 使用稀疏输入/输出处理您的软件包 API

    • 输入很容易与 spmatrix 或 sparray 一起使用。只要您使用 A.multiply(B) 进行逐元素计算，使用 A @ B 进行矩阵乘法，并使用 sparse.linalg.matrix_power 进行矩阵幂运算，在完成下一节中描述的“第一遍”迁移步骤后，您应该没有问题。您的代码将可以互换地处理这两种类型的输入。

    • 迁移来自您函数的稀疏输出需要更多思考。列出返回 spmatrix 对象的所有公共函数。检查您是否可以接受返回 sparray。这取决于您的库及其用户。如果您希望这些函数继续返回 spmatrix 或 sparray 对象，您通常可以使用稀疏输入来做到这一点，该输入也可以作为应返回哪种类型输出的信号。设计您的函数以返回输入的类型。该方法可以扩展到稠密输入。如果输入是 np.matrix 或以 np.matrix 作为其 ._baseclass 属性的掩码数组，则返回 spmatrix。否则，返回 sparray。如果没有这些输入，则可以使用其他两种方法：创建一个关键字参数来指示要返回哪个，或者创建一个新的函数（就像我们使用 eye_array 所做的那样），该函数具有相同的基本语法，但返回 sparray。您选择哪种方法应取决于您的库、您的用户以及您的偏好。

详细信息：构造函数#

以下四个函数是新的，仅处理 sparray：block_array、diags_array、eye_array 和 random_array。它们的签名如下：

def block_array(blocks, format=None, dtype=None):
def diags_array(diagonals, /, *, offsets=0, shape=None, format=None, dtype=None):
def eye_array(m, n=None, *, k=0, dtype=float, format=None):
def random_array(shape, density=0.01, format='coo', dtype=None, rng=None, data_sampler=None):

random_array 函数有一个 shape（2 元组）参数，而不是两个整数。并且 random_state 参数默认为 NumPy 的新 default_rng()。这与 spmatrix 的 rand 和 random 不同，后者默认为全局 RandomState 实例。如果您不太关心这些，则保持默认设置应该可以正常工作。如果您关心随机数的种子，则在切换函数时，您可能应该在此调用中添加 random_state=... 关键字参数。总之，要迁移到 random_array，请更改函数名称，将 shape 参数切换为单个元组参数，保留任何其他参数不变，并考虑应使用哪种 random_state= 参数（如果有）。

diags_array 函数对参数使用仅关键字规则。因此，您必须在偏移量参数前面键入 offsets=。从使用 diags 迁移时，这似乎很麻烦，但它有助于避免混淆并方便阅读。单个 shape 参数也取代了此迁移的两个整数。

需要仔细迁移的现有函数#

以下函数根据它们接收的输入类型返回 sparray 或 spmatrix：kron、kronsum、hstack、vstack、block_diag、tril 和 triu。它们的签名如下：

def kron(A, B, format=None):
def kronsum(A, B, format=None):
def hstack(blocks, format=None, dtype=None):
def vstack(blocks, format=None, dtype=None):
def block_diag(mats, format=None, dtype=None):
def tril(A, k=0, format=None):
def triu(A, k=0, format=None):

应该检查这些函数的使用情况并调整输入，以确保返回值是 sparray。反过来，输出也应该被视为 sparray。要返回 sparray，至少一个输入必须是 sparray。如果您使用列表的列表或 numpy 数组作为输入，则应将其中一个转换为稀疏数组以获得稀疏数组输出。

迁移时更改名称的函数#

函数

新函数

注释

eye

eye_array

identity

eye_array

diags

diags_array

仅关键字输入

spdiags

dia_array

形状为 2 元组

bmat

block

rand

random_array

形状为 2 元组，且为 default_rng

random

random_array

形状为 2 元组，且为 default_rng

详细信息：形状更改和缩减#

  • 使用值的一维列表进行构造

    • csr_array([1, 2, 3]).shape == (3,) 一维输入创建一个一维数组。

    • csr_matrix([1, 2, 3]).shape == (1, 3) 一维输入创建一个二维矩阵。

  • 索引和迭代

    • sparray 的索引允许使用一维对象，这些对象可以使用 np.newaxis 或 None 制成二维。例如，A[3, None, :] 给出一个二维行。使用隐式（未给定）列索引对二维 sparray 进行索引会产生一维结果，例如 A[3]（注意：最好不要这样做 - 请将其写为 A[3, :]）。如果需要二维结果，请在索引中使用 np.newaxis 或 None，或者将整数索引包装为列表，这样使用花式索引会给出二维结果，例如 A[[3], :]。

    • 对稀疏对象进行迭代：next(M) 产生一个稀疏二维行矩阵，next(A) 产生一个稀疏一维数组。

  • 沿轴的缩减操作会减小形状

    • M.sum(axis=1) 通过沿轴 1 求和返回一个二维行矩阵。

    • A.sum(axis=1) 返回一个沿轴 1 求和的一维 coo_array。一些缩减返回密集数组/矩阵而不是稀疏数组/矩阵

      方法

      结果

      sum(axis)

      密集

      mean(axis)

      密集

      argmin(axis)

      密集

      argmax(axis)

      密集

      min(axis)

      稀疏

      max(axis)

      稀疏

      nanmin(axis)

      稀疏

      nanmax(axis)

      稀疏

    通常，二维 sparray 输入会导致一维结果。二维 spmatrix 输入会导致二维结果。

  • 一些缩减返回标量。这些应该像以前一样工作，并且在迁移期间不需要考虑。例如，A.sum()

已删除的方法和属性#

  • 方法 get_shape、getrow、getcol、asfptype、getnnz、getH 和属性 .A 和 .H 仅存在于 spmatrix 上，而不存在于 sparray 上。建议您在开始转换为 sparray 之前，用替代方法替换它们的使用。

    函数

    替代方法

    M.get_shape()

    A.shape

    M.getformat()

    A.format

    M.asfptype(…)

    A.astype(…)

    M.getmaxprint()

    A.maxprint

    M.getnnz()

    A.nnz

    M.getnnz(axis)

    A.count_nonzero(axis)

    M.getH()

    A.conj().T

    M.getrow(i)

    A[i, :]

    M.getcol(j)

    A[:, j]

    M.A

    A.toarray()

    M.H

    A.conj().T

  • 不允许对 sparray 进行形状赋值 (M.shape = (2, 6))。相反，您应该使用 A.reshape。

  • M.getnnz() 返回存储值的数量，而不是非零值的数量。A.nnz 执行相同的操作。要获取非零值的数量，请使用 A.count_nonzero()。这不是迁移的新内容，但可能会令人困惑。

    要从 M.getnnz(axis=...) 的 axis 参数迁移，可以使用 A.count_nonzero(axis=...)，但它不是完全的替代品，因为它计算的是非零值而不是存储的值。区别在于显式存储的零值的数量。如果确实想要按轴计算存储值的数量，则需要使用一些 NumPy 工具。

    NumPy 工具方法适用于 COO、CSR、CSC 格式，因此请转换为其中一种格式。对于 CSR 和 CSC，主轴被压缩，np.diff(A.indptr) 返回一个密集的 1D 数组，其中包含每个主轴值（CSR 为行，CSC 为列）的存储值数量。可以使用 np.bincount(A.indices, minlength=N) 计算次轴，其中 N 是次轴的长度（例如，CSR 为 A.shape[1]）。bincount 函数适用于 COO 格式的任何轴，可以使用 A.coords[axis] 代替 A.indices。

使用测试查找 * 和 ** 的位置#

  • 在迁移代码时，区分标量乘法 * 和矩阵乘法 * 可能会很棘手。在理论上，Python 通过引入矩阵乘法运算符 @ 解决了这个问题。* 用于标量乘法，而 @ 用于矩阵乘法。但是，转换同时使用 * 的表达式可能会很棘手，并且会导致眼睛疲劳。幸运的是，如果你的代码有一个测试套件，涵盖了你需要转换的表达式，你可以使用它来查找 * 被用于涉及稀疏矩阵的矩阵乘法的地方。将这些地方更改为 @。

    该方法通过 monkey-patch spmatrix 类中的 dunder 方法，以便在 * 用于矩阵乘法时引发异常（而对于标量乘法则不引发异常）。测试套件将在这些位置标记失败。这里测试失败是成功的，因为它显示了需要进行更改的位置。将出错的 * 更改为 @，在附近查找其他类似的更改，然后再次运行测试。同样，此方法有助于查找 ** 用于矩阵幂的位置。当 @ 用于标量乘法时，SciPy 会引发异常，这将捕获您不应该更改的位置。因此，带有此 monkey-patch 的测试套件也会检查更正是否正确。

    将以下代码添加到你的 conftest.py 文件的顶部附近。然后本地运行你的测试。如果有许多矩阵表达式，你可能希望一次测试代码库的一个部分。快速阅读代码会发现，只要 * 在两个类似矩阵的对象（稀疏或密集）之间使用，以及只要 ** 用于矩阵幂，它就会引发 ValueError。

    import scipy


class _strict_mul_mixin:
    def __mul__(self, other):
        if not scipy.sparse._sputils.isscalarlike(other):
            raise ValueError('Operator * used here! Change to @?')
        return super().__mul__(other)

    def __rmul__(self, other):
        if not scipy.sparse._sputils.isscalarlike(other):
            raise ValueError('Operator * used here! Change to @?')
        return super().__rmul__(other)

    def __imul__(self, other):
        if not scipy.sparse._sputils.isscalarlike(other):
            raise ValueError('Operator * used here! Change to @?')
        return super().__imul__(other)

    def __pow__(self, *args, **kwargs):
        raise ValueError('spmatrix ** found! Use linalg.matrix_power?')

class _strict_coo_matrix(_strict_mul_mixin, scipy.sparse.coo_matrix):
    pass

class _strict_bsr_matrix(_strict_mul_mixin, scipy.sparse.bsr_matrix):
    pass

class _strict_csr_matrix(_strict_mul_mixin, scipy.sparse.csr_matrix):
    pass

class _strict_csc_matrix(_strict_mul_mixin, scipy.sparse.csc_matrix):
    pass

class _strict_dok_matrix(_strict_mul_mixin, scipy.sparse.dok_matrix):
    pass

class _strict_lil_matrix(_strict_mul_mixin, scipy.sparse.lil_matrix):
    pass

class _strict_dia_matrix(_strict_mul_mixin, scipy.sparse.dia_matrix):
    pass

scipy.sparse.coo_matrix = scipy.sparse._coo.coo_matrix = _strict_coo_matrix
scipy.sparse.bsr_matrix = scipy.sparse._bsr.bsr_matrix = _strict_bsr_matrix
scipy.sparse.csr_matrix = scipy.sparse._csr.csr_matrix = _strict_csr_matrix
scipy.sparse.csc_matrix = scipy.sparse._csc.csc_matrix = _strict_csc_matrix
scipy.sparse.dok_matrix = scipy.sparse._dok.dok_matrix = _strict_dok_matrix
scipy.sparse.lil_matrix = scipy.sparse._lil.lil_matrix = _strict_lil_matrix
scipy.sparse.dia_matrix = scipy.sparse._dia.dia_matrix = _strict_dia_matrix

scipy.sparse._compressed.csr_matrix = _strict_csr_matrix

scipy.sparse._construct.bsr_matrix = _strict_bsr_matrix
scipy.sparse._construct.coo_matrix = _strict_coo_matrix
scipy.sparse._construct.csc_matrix = _strict_csc_matrix
scipy.sparse._construct.csr_matrix = _strict_csr_matrix
scipy.sparse._construct.dia_matrix = _strict_dia_matrix

scipy.sparse._extract.coo_matrix = _strict_coo_matrix

scipy.sparse._matrix.bsr_matrix = _strict_bsr_matrix
scipy.sparse._matrix.coo_matrix = _strict_coo_matrix
scipy.sparse._matrix.csc_matrix = _strict_csc_matrix
scipy.sparse._matrix.csr_matrix = _strict_csr_matrix
scipy.sparse._matrix.dia_matrix = _strict_dia_matrix
scipy.sparse._matrix.dok_matrix = _strict_dok_matrix
scipy.sparse._matrix.lil_matrix = _strict_lil_matrix

索引数组 DType#

如果你向构造函数提供压缩索引，例如 csr_array((data, indices, indptr))，稀疏数组仅通过检查索引数组的 dtype 来设置索引 dtype，而稀疏矩阵也会检查索引值，并可能向下转换为 int32（有关更多详细信息，请参见gh-18509）。这意味着你可能会得到 int64 索引，而以前你得到的是 int32。你可以通过在实例化之前设置 dtype，或者在构造后重新转换来控制这一点。

两个稀疏实用函数可以帮助处理索引 dtype。在创建索引时，使用 get_index_dtype(arrays, maxval, check_contents) 来查找适合你的压缩索引的 dtype（int32 或 int64）。

使用 safely_cast_index_arrays(A, idx_dtype) 在构造后重新转换，同时确保在向下转换期间不会创建溢出。此函数实际上不会更改输入数组。将返回转换后的数组。仅在需要时才创建副本。因此，你可以使用 if indices is not A.indices: 检查是否进行了转换。

函数签名如下

def get_index_dtype(arrays=(), maxval=None, check_contents=False):
def safely_cast_index_arrays(A, idx_dtype=np.int32, msg=""):

示例用法包括以下 get_index_dtype 的用法

.. code-block:: python

    # select index dtype before construction based on shape
    shape = (3, 3)
    idx_dtype = scipy.sparse.get_index_dtype(maxval=max(shape))
    indices = np.array([0, 1, 0], dtype=idx_dtype)
    indptr = np.arange(3, dtype=idx_dtype)
    A = csr_array((data, indices, indptr), shape=shape)

和 safely_cast_index_arrays 的用法

.. code-block:: python

    # rescast after construction, raising exception if shape too big
    indices, indptr = scipy.sparse.safely_cast_index_arrays(B, np.int32)
    B.indices, B.indptr = indices, indptr

其他#

  • 二元运算符 +, -, *, /, @, !=, > 对稀疏和/或密集操作数起作用

    • 如果所有输入都是稀疏的，则输出通常也是稀疏的。例外是 /，它返回密集值（除以默认值 0 是 nan）。

    • 如果输入是稀疏和密集的混合，则结果通常是密集的（即 np.ndarray）。例外是 *，它是稀疏的；以及 /，它没有为 dense / sparse 实现，并且对于 sparse / dense 返回稀疏值。

  • 带有数组和/或矩阵操作数的二元运算符 +, -, *, /, @, !=, >

    • 如果所有输入都是数组，则输出是数组，对于矩阵也是如此。

    • 当将稀疏数组与稀疏矩阵混合时，前导操作数提供输出的类型，例如，sparray + spmatrix 给出稀疏数组，而反转顺序则给出稀疏矩阵。

    • 当将密集矩阵与稀疏数组混合时，结果通常是数组，但比较除外，例如，> 返回密集矩阵。

    • 当将密集数组与稀疏矩阵混合时，结果通常是矩阵，但 array @ sparse matrix 返回密集数组除外。

在此页面上