scipy.sparse行列の要素・行・列・部分行列の値を取得・更新

Posted: | Tags: Python, SciPy

SciPyのscipy.sparse行列(疎行列)はNumPy配列ndarrayと同様にインデックス[]やスライス[:]で要素や行、列、部分行列の値を取得したり更新(変更)したりすることができる。

ここでは以下の内容について説明する。

  • scipy.sparse行列の種類(クラス)による違い
  • 行・列を取得: getrow(), getcol()
  • インデックスで要素を取得・更新
  • スライスで行・列・部分行列を取得・更新
  • 処理速度の比較

scipy.sparseの基本については以下の記事を参照。

以下の説明内容およびサンプルコードはSciPy1.3.1のもの。バージョンが違うと異なる場合があるので注意。

scipy.sparse行列の種類(クラス)による違い

NumPy配列ndarrayと異なり、scipy.sparseにはデータの格納形式によってCSR(csr_matrix)、LIL(lil_matrix)などのいくつかのクラスがあり、クラスによって同じ処理でも処理速度が異なる。

また、COO(coo_matrix)はインデックスやスライスなど[]を使う操作はサポートされていないといった違いもある。

処理速度比較の実験結果は後述するが、ざっくりまとめると以下のような使い分け。

  • 列の取得: CSC
  • 要素・行・部分行列の取得: LIL
  • 値の更新: LIL

以下、CSR(csr_matrix)、CSC(csc_matrix)、 LIL(lil_matrix)、COO(coo_matrix)を中心に説明する。

行・列を取得: getrow(), getcol()

scipy.sparse行列のクラスには行・列を取得するためのgetrow(), getcol()メソッドが用意されている。基底クラスであるscipy.sparse.spmatrixのメソッドなので、すべてのクラスで使用可能。

getrow()の例。引数に0始まりの行番号を指定する。1行x列の行列が取得できる。

from scipy.sparse import csr_matrix, csc_matrix, coo_matrix, lil_matrix

l = [[1, 0, 0, 0],
     [0, 2, 0, 0],
     [0, 0, 3, 0],
     [0, 0, 0, 4]]

csr = csr_matrix(l)
csc = csc_matrix(l)
coo = coo_matrix(l)
lil = lil_matrix(l)

print(csr.getrow(0))
#   (0, 0)  1

print(type(csr.getrow(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

print(csr.getrow(0).shape)
# (1, 4)

print(csr.getrow(0).toarray())
# [[1 0 0 0]]

csc_matrix, coo_matrixは自分自身のクラスではなくcsr_matrixを返すので注意。

print(type(csc.getrow(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(coo.getrow(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(lil.getrow(0)))
# <class 'scipy.sparse.lil.lil_matrix'>

getcol()の例。引数に0始まりの列番号を指定する。x行1列の行列が取得できる。

print(csr.getcol(0))
#   (0, 0)  1

print(type(csr.getcol(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

print(csr.getcol(0).shape)
# (4, 1)

print(csr.getcol(0).toarray())
# [[1]
#  [0]
#  [0]
#  [0]]

coo_matrix, lil_matrixは自分自身のクラスではなくcsr_matrixを返すので注意。

print(type(csc.getcol(0)))
# <class 'scipy.sparse.csc.csc_matrix'>

print(type(coo.getcol(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(lil.getcol(0)))
# <class 'scipy.sparse.csr.csr_matrix'>

getrow(), getcol()が返すのは元の行列のコピー。いずれかの値を更新しても片方の値はそのまま。

lil_row = lil.getrow(0)

lil_row[0, 0] = 100

print(lil.toarray())
# [[1 0 0 0]
#  [0 2 0 0]
#  [0 0 3 0]
#  [0 0 0 4]]

print(lil_row.toarray())
# [[100   0   0   0]]

インデックスで要素を取得・更新

numpy.ndarrayのようにscipy.sparseでも[行のインデックス, 列のインデックス]で要素を指定して値を取得できる。

coo_matrixではこのような操作がサポートされていない。

l = [[1, 0, 0, 0],
     [0, 2, 0, 0],
     [0, 0, 3, 0],
     [0, 0, 0, 4]]

csr = csr_matrix(l)
csc = csc_matrix(l)
coo = coo_matrix(l)
lil = lil_matrix(l)

print(csr[1, 1])
# 2

print(csc[1, 1])
# 2

print(lil[1, 1])
# 2

# print(coo[1, 1])
# TypeError: 'coo_matrix' object is not subscriptable

指定した位置に新たな値を代入することも可能。新たな要素を追加したり、既存の要素の値を変更したりできる。

lil[0, 0] = 10
lil[0, 1] = 100

print(lil)
#   (0, 0)  10
#   (0, 1)  100
#   (1, 1)  2
#   (2, 2)  3
#   (3, 3)  4

print(lil.toarray())
# [[ 10 100   0   0]
#  [  0   2   0   0]
#  [  0   0   3   0]
#  [  0   0   0   4]]

lil_matrixの場合、0を代入するとデータが削除される。

lil[1, 1] = 0

print(lil)
#   (0, 0)  10
#   (0, 1)  100
#   (2, 2)  3
#   (3, 3)  4

print(lil.toarray())
# [[ 10 100   0   0]
#  [  0   0   0   0]
#  [  0   0   3   0]
#  [  0   0   0   4]]

csr_matrix, csc_matrixでは新たな要素を追加しようとするとSparseEfficiencyWarning(非効率だという警告)が出る。ここではコメントアウトしているが、エラーではなく警告なので実行することは可能。

csr[0, 0] = 10

# csr[0, 1] = 100
# SparseEfficiencyWarning: Changing the sparsity structure of a csr_matrix is expensive. lil_matrix is more efficient.

また、csr_matrix, csc_matrixでは、0を代入した場合、lil_matrixと異なり値が0の要素として残る。toarray()などで変換した結果は同じだが、保持しているデータ数が変わってくるので注意。

csr[1, 1] = 0

print(csr)
#   (0, 0)  10
#   (1, 1)  0
#   (2, 2)  3
#   (3, 3)  4

print(csr.toarray())
# [[10  0  0  0]
#  [ 0  0  0  0]
#  [ 0  0  3  0]
#  [ 0  0  0  4]]

このような既存の要素の値の変更の場合は特に警告は出ない。

スライスで行・列・部分行列を取得・更新

numpy.ndarrayのようにscipy.sparseでもスライスで範囲を指定して行や列、部分行列を取得できる。

上述のインデックスと同じく、coo_matrixではこのような操作がサポートされていない。

いくつか例を示す。

行の取得。numpy.ndarrayと同様に末尾の, :は省略できる。

print(csr[0, :])
#   (0, 0)  1

print(csr[0, :].toarray())
# [[1 0 0 0]]

print(csr[0].toarray())
# [[1 0 0 0]]

列の取得。

print(csr[:, 0])
#   (0, 0)  1

print(csr[:, 0].toarray())
# [[1]
#  [0]
#  [0]
#  [0]]

部分行列の取得。start:stop:stepstepを指定することも可能。

print(csr[1:3, 1:3])
#   (0, 0)  2
#   (1, 1)  3

print(csr[1:3, 1:3].toarray())
# [[2 0]
#  [0 3]]

print(csr[:, ::2])
#   (0, 0)  1
#   (2, 1)  3

print(csr[:, ::2].toarray())
# [[1 0]
#  [0 0]
#  [0 3]
#  [0 0]]

上述のgetrow(), getcol()メソッドと異なり、どのような場合もスライスは元のオブジェクトと同じクラスのオブジェクトを返す。

print(type(csr[0]))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(csc[0]))
# <class 'scipy.sparse.csc.csc_matrix'>

print(type(lil[0]))
# <class 'scipy.sparse.lil.lil_matrix'>

print(type(csr[:, 0]))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(csc[:, 0]))
# <class 'scipy.sparse.csc.csc_matrix'>

print(type(lil[:, 0]))
# <class 'scipy.sparse.lil.lil_matrix'>

print(type(csr[1:3, 1:3]))
# <class 'scipy.sparse.csr.csr_matrix'>

print(type(csc[1:3, 1:3]))
# <class 'scipy.sparse.csc.csc_matrix'>

print(type(lil[1:3, 1:3]))
# <class 'scipy.sparse.lil.lil_matrix'>

また、numpy.ndarrayと異なり、scipy.sparseのスライスはコピーを返す模様。スライスで取得したオブジェクトの要素を更新しても元のオブジェクトの要素は元のまま。

例はcsr_matrixだが、csc_matrixでもlil_matrixでも同様。

csr_slice = csr[1:3, 1:3]
csr_slice[0, 0] = 100

print(csr.toarray())
# [[1 0 0 0]
#  [0 2 0 0]
#  [0 0 3 0]
#  [0 0 0 4]]

print(csr_slice.toarray())
# [[100   0]
#  [  0   3]]

スライスで指定した部分に新たな値を代入することもできる。

リストを代入。

lil[0] = [10, 20, 30, 40]

print(lil)
#   (0, 0)  10
#   (0, 1)  20
#   (0, 2)  30
#   (0, 3)  40
#   (1, 1)  2
#   (2, 2)  3
#   (3, 3)  4

print(lil.toarray())
# [[10 20 30 40]
#  [ 0  2  0  0]
#  [ 0  0  3  0]
#  [ 0  0  0  4]]

numpy.ndarrayを代入。

lil[1:3, 1:3] = np.arange(4).reshape(2, 2) * 100

print(lil)
#   (0, 0)  10
#   (0, 1)  20
#   (0, 2)  30
#   (0, 3)  40
#   (1, 2)  100
#   (2, 1)  200
#   (2, 2)  300
#   (3, 3)  4

print(lil.toarray())
# [[ 10  20  30  40]
#  [  0   0 100   0]
#  [  0 200 300   0]
#  [  0   0   0   4]]

scipy.sparse行列を代入。異なるクラスでもOK。

lil[:, 0] = csr[:, 3]

print(lil)
#   (0, 1)  20
#   (0, 2)  30
#   (0, 3)  40
#   (1, 2)  100
#   (2, 1)  200
#   (2, 2)  300
#   (3, 0)  4
#   (3, 3)  4

print(lil.toarray())
# [[  0  20  30  40]
#  [  0   0 100   0]
#  [  0 200 300   0]
#  [  4   0   0   4]]

形状が一致していないとエラーとなる。

# lil[1:3, 1:3] = [10, 20, 30, 40]
# ValueError: shape mismatch: objects cannot be broadcast to a single shape

上述のインデックスと同様に、csr_matrix, csc_matrixでは新たな要素を追加しようとするとSparseEfficiencyWarningが出る。エラーではなく警告なので実行することは可能。

csr[0] = [0, 0, 0, 100]
# /usr/local/lib/python3.7/site-packages/scipy/sparse/_index.py:127: SparseEfficiencyWarning: Changing the sparsity structure of a csr_matrix is expensive. lil_matrix is more efficient.
#   self._set_arrayXarray(i, j, x)

print(csr)
#   (0, 0)  0
#   (0, 1)  0
#   (0, 2)  0
#   (0, 3)  100
#   (1, 1)  2
#   (2, 2)  3
#   (3, 3)  4

print(csr.toarray())
# [[  0   0   0 100]
#  [  0   2   0   0]
#  [  0   0   3   0]
#  [  0   0   0   4]]

例の結果から分かるように、0が代入された場合、lil_matrixでは要素が削除されるが、csr_matrix, csc_matrixでは値が0の要素として残る。

処理速度の比較

scipy.sparseの各クラスにおける要素・行・列・部分行列の値取得の処理速度比較の結果を示す。

csr_matrix, csc_matrix, lil_matrixおよび要素のアクセスが高速とされるdok_matrixで比較を行った。

結果を先に示しておくと最速は以下の通り。

  • 1行取得: lil_matrixgetrow()
  • 1列取得: csc_matrixgetcol()
  • 要素取得: lil_matrix
  • 複数行取得: lil_matrix
  • 複数列取得: csc_matrix
  • 部分行列取得: lil_matrix

行列のサイズ、非ゼロ要素の割合によって変わる可能性があるので参考程度ではあるが、列は csc_matrix、それ以外はlil_matrixという使い分け。dok_matrixはあまりメリットが感じられなかった。

lil_matrixは列取得が極端に遅いので要注意。

以下、実験結果。

1000行1000列の単位行列(対角成分が1、ほかは0)を用いた場合。ここではscipy.sparseeye()を使っている。

n = 1000

csr = eye(n, format='csr')
csc = eye(n, format='csc')
lil = eye(n, format='lil')
dok = eye(n, format='dok')

%%timeit
csr.getrow(0)
# 48.3 µs ± 2.71 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc.getrow(0)
# 139 µs ± 11.6 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil.getrow(0)
# 18.1 µs ± 522 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
dok.getrow(0)
# 700 µs ± 60.7 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
csr[0]
# 78.3 µs ± 5.64 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[0]
# 74.5 µs ± 715 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[0]
# 43.9 µs ± 4.06 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
dok[0]
# 462 µs ± 16.7 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
csr.getcol(0)
# 67.5 µs ± 7.8 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc.getcol(0)
# 49.4 µs ± 2.21 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil.getcol(0)
# 955 µs ± 15.5 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
dok.getcol(0)
# 12.3 ms ± 154 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

%%timeit
csr[:, 0]
# 86.7 µs ± 7.77 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[:, 0]
# 77 µs ± 6.75 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[:, 0]
# 442 µs ± 31.9 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
dok[:, 0]
# 929 µs ± 36.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
csr[0, 0]
# 31.8 µs ± 454 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[0, 0]
# 104 µs ± 5.93 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[0, 0]
# 3.37 µs ± 123 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

%%timeit
dok[0, 0]
# 18.1 µs ± 134 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

%%timeit
csr[:10]
# 71.4 µs ± 929 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[:10]
# 184 µs ± 8.37 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[:10]
# 47.2 µs ± 1.57 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
dok[:10]
# 632 µs ± 7.13 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
csr[:, :10]
# 77.5 µs ± 1.32 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[:, :10]
# 74.9 µs ± 5.32 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[:, :10]
# 420 µs ± 14.9 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
dok[:, :10]
# 905 µs ± 6.62 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%%timeit
csr[:10, :10]
# 73.9 µs ± 1.63 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
csc[:10, :10]
# 178 µs ± 1.79 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
lil[:10, :10]
# 47.5 µs ± 925 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%%timeit
dok[:10, :10]
# 210 µs ± 6.29 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

1000行1000列に10000個の要素をランダムに配置した行列でも最速のクラスは変わらず。ここでの掲載は省略する。結果の確認は以下のリンクから。

関連カテゴリー

関連記事