Sparse Arrays

ジュリアは、SparseArrays 標準ライブラリモジュールでスパースベクトルと sparse matrices をサポートしています。スパース配列は、ゼロが十分に含まれている配列であり、特別なデータ構造に格納することで、密な配列と比較して、スペースと実行時間の節約が可能になります。

異なるスパースストレージタイプ、マルチディメンショナルスパース配列などを実装した外部パッケージは、Noteworthy External Sparse Packages にあります。

Compressed Sparse Column (CSC) Sparse Matrix Storage

ジュリアでは、スパース行列は Compressed Sparse Column (CSC) format に格納されます。ジュリアのスパース行列は SparseMatrixCSC{Tv,Ti} という型を持ち、ここで Tv は格納される値の型、Ti は列ポインタと行インデックスを格納するための整数型です。SparseMatrixCSC の内部表現は次のようになります：

struct SparseMatrixCSC{Tv,Ti<:Integer} <: AbstractSparseMatrixCSC{Tv,Ti}
    m::Int                  # Number of rows
    n::Int                  # Number of columns
    colptr::Vector{Ti}      # Column j is in colptr[j]:(colptr[j+1]-1)
    rowval::Vector{Ti}      # Row indices of stored values
    nzval::Vector{Tv}       # Stored values, typically nonzeros
end

圧縮スパース列ストレージは、スパース行列の列内の要素に簡単かつ迅速にアクセスできるようにしますが、行によるスパース行列へのアクセスはかなり遅くなります。CSC構造内で、以前に保存されていなかったエントリを1つずつ挿入する操作は遅くなる傾向があります。これは、挿入ポイントを超えるスパース行列のすべての要素を1つずつ移動させる必要があるためです。

すべてのスパース行列に対する操作は、パフォーマンスを向上させるためにCSCデータ構造を活用し、高コストの操作を回避するように慎重に実装されています。

CSC形式のデータが他のアプリケーションやライブラリからあり、Juliaにインポートしたい場合は、1ベースのインデックスを使用することを確認してください。各列の行インデックスはソートされている必要があり、そうでない場合、行列は正しく表示されません。SparseMatrixCSCオブジェクトにソートされていない行インデックスが含まれている場合、ソートする簡単な方法はダブルトランスポーズを行うことです。トランスポーズ操作は遅延されるため、各トランスポーズを具現化するためにコピーを作成してください。

いくつかのアプリケーションでは、SparseMatrixCSCに明示的なゼロ値を保存することが便利です。これらはBaseの関数によって受け入れられます（ただし、変更操作で保持される保証はありません）。このように明示的に保存されたゼロは、多くのルーチンによって構造的な非ゼロとして扱われます。nnz関数は、構造的でないゼロを含むスパースデータ構造に明示的に保存された要素の数を返します。正確な数の数値非ゼロをカウントするには、count(!iszero, x)を使用してください。これはスパース行列のすべての保存された要素を検査します。dropzerosおよびインプレースのdropzeros!は、スパース行列から保存されたゼロを削除するために使用できます。

julia> A = sparse([1, 1, 2, 3], [1, 3, 2, 3], [0, 1, 2, 0])
3×3 SparseMatrixCSC{Int64, Int64} with 4 stored entries:
 0  ⋅  1
 ⋅  2  ⋅
 ⋅  ⋅  0

julia> dropzeros(A)
3×3 SparseMatrixCSC{Int64, Int64} with 2 stored entries:
 ⋅  ⋅  1
 ⋅  2  ⋅
 ⋅  ⋅  ⋅

Sparse Vector Storage

スパースベクトルは、スパース行列の圧縮スパース列形式に近い形で保存されます。Juliaでは、スパースベクトルは型 SparseVector{Tv,Ti} を持ち、ここで Tv は保存される値の型、Ti はインデックスの整数型です。内部表現は次のようになります:

struct SparseVector{Tv,Ti<:Integer} <: AbstractSparseVector{Tv,Ti}
    n::Int              # Length of the sparse vector
    nzind::Vector{Ti}   # Indices of stored values
    nzval::Vector{Tv}   # Stored values, typically nonzeros
end

SparseMatrixCSCについて、SparseVector型は明示的に保存されたゼロも含むことができます。（Sparse Matrix Storageを参照してください。）

Sparse Vector and Matrix Constructors

スパース配列を作成する最も簡単な方法は、Juliaが密な配列を扱うために提供するzeros関数に相当する関数を使用することです。代わりにスパース配列を生成するには、同じ名前にspプレフィックスを付けて使用できます：

julia> spzeros(3)
3-element SparseVector{Float64, Int64} with 0 stored entries

sparse 関数は、スパース配列を構築する便利な方法です。たとえば、スパース行列を構築するために、行インデックスのベクトル I、列インデックスのベクトル J、および格納された値のベクトル V を入力できます（これは COO (coordinate) format とも呼ばれます）。sparse(I,J,V) は、S[I[k], J[k]] = V[k] となるスパース行列を構築します。等価なスパースベクトルコンストラクタは sparsevec で、（行）インデックスベクトル I と格納された値のベクトル V を取り、R[I[k]] = V[k] となるスパースベクトル R を構築します。

julia> I = [1, 4, 3, 5]; J = [4, 7, 18, 9]; V = [1, 2, -5, 3];

julia> S = sparse(I,J,V)
5×18 SparseMatrixCSC{Int64, Int64} with 4 stored entries:
⎡⠀⠈⠀⠀⠀⠀⠀⠀⢀⎤
⎣⠀⠀⠀⠂⡀⠀⠀⠀⠀⎦

julia> R = sparsevec(I,V)
5-element SparseVector{Int64, Int64} with 4 stored entries:
  [1]  =  1
  [3]  =  -5
  [4]  =  2
  [5]  =  3

sparse および sparsevec 関数の逆は findnz であり、これはスパース配列を作成するために使用された入力（ゼロに等しい保存されたエントリを含む）を取得します。 findall(!iszero, x) は、x の非ゼロエントリのカーテシアンインデックスを返します（ゼロに等しい保存されたエントリは含まれません）。

julia> findnz(S)
([1, 4, 5, 3], [4, 7, 9, 18], [1, 2, 3, -5])

julia> findall(!iszero, S)
4-element Vector{CartesianIndex{2}}:
 CartesianIndex(1, 4)
 CartesianIndex(4, 7)
 CartesianIndex(5, 9)
 CartesianIndex(3, 18)

julia> findnz(R)
([1, 3, 4, 5], [1, -5, 2, 3])

julia> findall(!iszero, R)
4-element Vector{Int64}:
 1
 3
 4
 5

別の方法として、sparse 関数を使用して、密な配列を疎な配列に変換することができます：

julia> sparse(Matrix(1.0I, 5, 5))
5×5 SparseMatrixCSC{Float64, Int64} with 5 stored entries:
 1.0   ⋅    ⋅    ⋅    ⋅
  ⋅   1.0   ⋅    ⋅    ⋅
  ⋅    ⋅   1.0   ⋅    ⋅
  ⋅    ⋅    ⋅   1.0   ⋅
  ⋅    ⋅    ⋅    ⋅   1.0

julia> sparse([1.0, 0.0, 1.0])
3-element SparseVector{Float64, Int64} with 2 stored entries:
  [1]  =  1.0
  [3]  =  1.0

他の方向に進むには、Array コンストラクタを使用できます。issparse 関数は、行列がスパースであるかどうかを照会するために使用できます。

julia> issparse(spzeros(5))
true

Sparse matrix operations

Arithmetic operations on sparse matrices also work as they do on dense matrices. Indexing of, assignment into, and concatenation of sparse matrices work in the same way as dense matrices. Indexing operations, especially assignment, are expensive, when carried out one element at a time. In many cases it may be better to convert the sparse matrix into (I,J,V) format using findnz, manipulate the values or the structure in the dense vectors (I,J,V), and then reconstruct the sparse matrix.

Correspondence of dense and sparse methods

次の表は、スパース行列の組み込みメソッドと、それに対応する密行列タイプのメソッドとの対応関係を示しています。一般に、スパース行列を生成するメソッドは、与えられたスパース行列 S と同じスパースパターンに従うか、結果として得られるスパース行列が密度 d を持つ点で、密行列の対応メソッドとは異なります。つまり、各行列要素は非ゼロである確率 d を持っています。

詳細は標準ライブラリリファレンスの Sparse Vectors and Matrices セクションに記載されています。

Sparse	Dense	Description
`spzeros(m,n)`	`zeros(m,n)`	Creates a m-by-n matrix of zeros. (`spzeros(m,n)` is empty.)
`sparse(I,n,n)`	`Matrix(I,n,n)`	Creates a n-by-n identity matrix.
`sparse(A)`	`Array(S)`	Interconverts between dense and sparse formats.
`sprand(m,n,d)`	`rand(m,n)`	Creates a m-by-n random matrix (of density d) with iid non-zero elements distributed uniformly on the half-open interval $[0, 1)$.
`sprandn(m,n,d)`	`randn(m,n)`	Creates a m-by-n random matrix (of density d) with iid non-zero elements distributed according to the standard normal (Gaussian) distribution.
`sprandn(rng,m,n,d)`	`randn(rng,m,n)`	Creates a m-by-n random matrix (of density d) with iid non-zero elements generated with the `rng` random number generator

Sparse Linear Algebra

スパース行列ソルバーは、SuiteSparse から関数を呼び出します。以下の因子分解が利用可能です：

Type	Description
`CHOLMOD.Factor`	Cholesky and LDLt factorizations
`UMFPACK.UmfpackLU`	LU factorization
`SPQR.QRSparse`	QR factorization

これらの因数分解は、Sparse Linear Algebra API sectionでより詳細に説明されています:

SparseArrays API

SparseArrays.AbstractSparseArray — Type

AbstractSparseArray{Tv,Ti,N}

Tv型の要素とTi型のインデックスを持つN次元のスパース配列（または配列のような型）のスーパタイプです。 SparseMatrixCSC、SparseVectorおよびSuiteSparse.CHOLMOD.Sparseはこのサブタイプです。

コンポーネント	説明
`L`	`LU` の `L`（下三角）部分
`U`	`LU` の `U`（上三角）部分
`p`	右置換 `Vector`
`q`	左置換 `Vector`
`Rs`	スケーリングファクターの `Vector`
`:`	`(L,U,p,q,Rs)` コンポーネント

コンポーネント	説明
`F.L`	`LU` の `L`（下三角）部分
`F.U`	`LU` の `U`（上三角）部分
`F.p`	（右）置換 `Vector`
`F.P`	（右）置換 `Matrix`

サポートされる関数	`LU`	`LU{T,Tridiagonal{T}}`
`/`	✓
`\`	✓	✓
`inv`	✓	✓
`det`	✓	✓
`logdet`	✓	✓
`logabsdet`	✓	✓
`size`	✓	✓