從 spmatrix 遷移到 sparray#

本文檔為將程式碼從 scipy.sparse 中的稀疏矩陣轉換為稀疏陣列提供指南。

從稀疏矩陣到稀疏陣列的變更，反映了從 np.matrix 到 np.ndarray 的轉換。本質上，我們必須從以全 2D 矩陣乘法為中心的 matrix 物件，轉移到支援矩陣乘法運算子和逐元素計算的 1D 或 2D「陣列」物件。

符號：在本指南中，我們將稀疏陣列類別統稱為 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 始終是 2D（行 x 列），即使例如 M.min(axis=0)。sparray A 可以是 1D 或 2D。Numpy 純量會針對完整 (0D) 約簡傳回，即 M.min()。
迭代 sparray 會產生 1D sparray。迭代 spmatrix 會產生 2D 列 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。您選擇哪種方法應取決於您的程式庫、您的使用者以及您的偏好。

建議的遷移步驟#

第一階段（在程式碼中保留 spmatrix）
- 在您的 spmatrix 程式碼中，將 * 變更為 @ 以進行矩陣乘法。請注意，稀疏純量乘法應使用 *。（請參閱下方的輔助程式碼使用測試來尋找 * 和 ** 的位置）
- 矩陣乘冪，例如 M**3，應轉換為 scipy.sparse.linalg.matrix_power(A, 3)
- 實作不受支援的函式/方法的替代方案，例如 A.getnnz() -> A.nnz（請參閱下方的已移除的方法和屬性）。
- 根據需要變更關於 issparse() 和 isspmatrix() 的任何邏輯。通常，這表示將 isspmatrix 替換為 issparse，將 isspmatrix_csr(G) 替換為 issparse(G) and G.format == "csr"。此外，isspmatrix_csr(G) or isspmatrix_csc(G) 變成 issparse(G) and G.format in ['csr', 'csc']。git 搜尋慣用語 git grep 'isspm[a-z_]*(' 可以協助尋找這些。
- 將所有 spdiags 呼叫轉換為 dia_matrix。請參閱 spdiags 中的文件。搜尋 spdiags 是您在這裡需要的全部。
- 在產生的程式碼上執行所有測試。您仍然在使用 spmatrix，而不是 sparray。但您的程式碼和測試已為變更做好準備，您應該能夠將 sparray 作為輸入傳遞到您的程式碼中，並讓它們大致「正常運作」。
第二階段（切換到 sparray）
- 將建構函式（如 diags 和 triu）轉換為陣列版本（請參閱下方的詳細資訊：建構函式）。
- 將所有 *_matrix 建構函式呼叫重新命名為 *_array。
- 檢查遷移導致 1D 傳回值的所有函式/方法。這些主要是索引和約簡函式（請參閱下方的詳細資訊：形狀變更和約簡）。
- 檢查您迭代 spmatrix 的所有位置，並變更它們以說明 sparray 產生 1D sparray 而不是 2D spmatrix 的情況。
- 尋找並變更程式碼中使用 np.matrix 功能的位置。將這些轉換為 np.ndarray 功能。
- 如果您的程式碼使用 mmread、hb_read 或 loadmat 從檔案讀取稀疏資料，請在這些函式中使用新的關鍵字引數 spmatrix=False 以讀取到 sparray。
- 如果您使用的稀疏程式庫僅接受稀疏表示的 int32 索引陣列，我們建議使用即時轉換。在您呼叫需要 int32 的程式碼之前，轉換為 int32。
- sparray 根據輸入陣列的 dtype 而不是陣列中的值來選取索引 dtype。因此，如果您希望索引陣列為 int32，您需要確保每個索引陣列（例如您傳遞給 csr_array 的 indptr）的 dtype 為 int32。使用 spmatrix 時，很容易使用 numpy 陣列的預設 int64 dtype，並依靠稀疏建構函式在值很小時進行向下轉換。但這種向下轉換在處理其他矩陣、切片或算術運算式時會導致額外的重新轉換。對於 sparray，您仍然可以依靠建構函式來選擇 dtype。但您也可以透過傳入索引陣列的 dtype 而不是它們的值來選擇索引 dtype。因此，如果您想要 int32，請在建立索引陣列時設定 dtype，例如 indices = np.array([1,3,6], dtype=np.int32) 或 indptr = np.arange(9, dtype=np.int32)。如需更多資訊，請參閱下方的索引陣列 DType。在許多設定中，索引陣列 dtype 並不重要，您可以讓建構函式為 sparray 和 spmatrix 選擇 dtype。
- 測試您的程式碼。並閱讀您的程式碼。您已遷移到 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 遷移似乎很麻煩，但它有助於避免混淆並簡化閱讀。單一形狀參數也取代了此遷移的兩個整數。

需要謹慎遷移的現有函式#

以下函式會傳回 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

函式	新函式	註解
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

詳細資訊：形狀變更和約簡#

使用 1 維值清單建構
- csr_array([1, 2, 3]).shape == (3,) 1D 輸入會建立 1D 陣列。
- csr_matrix([1, 2, 3]).shape == (1, 3) 1D 輸入會建立 2D 矩陣。
索引和迭代
- sparray 的索引允許使用 1D 物件，這些物件可以使用 np.newaxis 或 None 設為 2D。例如，A[3, None, :] 會產生 2D 列。使用隱含（未給定）列索引對 2D sparray 進行索引會產生 1D 結果，例如 A[3]（注意：最好不要這樣做 - 請將其寫為 A[3, :]）。如果您需要 2D 結果，請在索引中使用 np.newaxis 或 None，或將整數索引包裝為清單，花式索引會產生 2D，例如 A[[3], :]。
- 迭代稀疏物件：next(M) 會產生稀疏 2D 列矩陣，next(A) 會產生稀疏 1D 陣列。
沿軸的約簡運算會縮減形狀
- M.sum(axis=1) 會傳回 2D 列矩陣，方法是沿軸 1 進行求和。
- A.sum(axis=1) 會傳回 1D coo_array，方法是沿軸 1 進行求和。某些約簡會傳回稠密陣列/矩陣，而不是稀疏陣列/矩陣
  
  方法
  
  結果
  
  sum(axis)
  
  稠密
  
  mean(axis)
  
  稠密
  
  argmin(axis)
  
  稠密
  
  argmax(axis)
  
  稠密
  
  min(axis)
  
  稀疏
  
  max(axis)
  
  稀疏
  
  nanmin(axis)
  
  稀疏
  
  nanmax(axis)
  
  稀疏
一般而言，2D sparray 輸入會產生 1D 結果。2D spmatrix 輸入會產生 2D 結果。
某些約簡會傳回純量。這些應與之前一樣運作，在遷移期間不應考慮。例如 A.sum()

方法	結果
sum(axis)	稠密
mean(axis)	稠密
argmin(axis)	稠密
argmax(axis)	稠密
min(axis)	稀疏
max(axis)	稀疏
nanmin(axis)	稀疏
nanmax(axis)	稀疏

已移除的方法和屬性#

方法 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.indices 的位置使用 A.coords[axis]。

函式	替代方案
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

使用測試來尋找 * 和 ** 的位置#

在您遷移程式碼時，區分純量乘法 * 與矩陣乘法 * 可能很棘手。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）。這表示您可能會在過去取得 int32 索引時取得 int64 索引。您可以透過在具現化之前設定 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）。例外情況是 * 運算子，其結果為稀疏格式；而 / 運算子對於 密集 / 稀疏 格式未實作，但對於 稀疏 / 密集 格式則回傳稀疏格式。
二元運算子 +, -, *, /, @, !=, > 具有陣列和/或矩陣運算元
- 如果所有輸入都是陣列，則輸出也是陣列；對於矩陣來說情況亦同。
- 當混合稀疏陣列和稀疏矩陣時，前導運算元會決定輸出的類型，例如 sparray + spmatrix 會產生稀疏陣列，而反轉順序則會產生稀疏矩陣。
- 當混合密集矩陣和稀疏陣列時，結果通常是陣列，但比較運算子是例外，例如 > 會回傳密集矩陣。
- 當混合密集陣列和稀疏矩陣時，結果通常是矩陣，但 array @ sparse matrix 是例外，它會回傳密集陣列。