SciPy

Sparse matrices (scipy.sparse)

SciPy 2-D sparse matrix package for numeric data.

Contents

Sparse matrix classes

bsr_matrix(arg1[, shape, dtype, copy, blocksize]) Block Sparse Row matrix This can be instantiated in several ways: bsr_matrix(D, [blocksize=(R,C)]) where D is a dense matrix or 2-D ndarray.
coo_matrix(arg1[, shape, dtype, copy]) A sparse matrix in COOrdinate format.
csc_matrix(arg1[, shape, dtype, copy]) Compressed Sparse Column matrix This can be instantiated in several ways: csc_matrix(D) with a dense matrix or rank-2 ndarray D csc_matrix(S) with another sparse matrix S (equivalent to S.tocsc()) csc_matrix((M, N), [dtype]) to construct an empty matrix with shape (M, N) dtype is optional, defaulting to dtype=’d’.
csr_matrix(arg1[, shape, dtype, copy]) Compressed Sparse Row matrix This can be instantiated in several ways: csr_matrix(D) with a dense matrix or rank-2 ndarray D csr_matrix(S) with another sparse matrix S (equivalent to S.tocsr()) csr_matrix((M, N), [dtype]) to construct an empty matrix with shape (M, N) dtype is optional, defaulting to dtype=’d’.
dia_matrix(arg1[, shape, dtype, copy]) Sparse matrix with DIAgonal storage This can be instantiated in several ways: dia_matrix(D) with a dense matrix dia_matrix(S) with another sparse matrix S (equivalent to S.todia()) dia_matrix((M, N), [dtype]) to construct an empty matrix with shape (M, N), dtype is optional, defaulting to dtype=’d’.
dok_matrix(arg1[, shape, dtype, copy]) Dictionary Of Keys based sparse matrix.
lil_matrix(arg1[, shape, dtype, copy]) Row-based linked list sparse matrix This is a structure for constructing sparse matrices incrementally.
spmatrix([maxprint]) This class provides a base class for all sparse matrices.

Functions

Building sparse matrices:

eye(m[, n, k, dtype, format]) Sparse matrix with ones on diagonal Returns a sparse (m x n) matrix where the k-th diagonal is all ones and everything else is zeros.
identity(n[, dtype, format]) Identity matrix in sparse format Returns an identity matrix with shape (n,n) using a given sparse format and dtype.
kron(A, B[, format]) kronecker product of sparse matrices A and B
kronsum(A, B[, format]) kronecker sum of sparse matrices A and B
diags(diagonals[, offsets, shape, format, dtype]) Construct a sparse matrix from diagonals.
spdiags(data, diags, m, n[, format]) Return a sparse matrix from diagonals.
block_diag(mats[, format, dtype]) Build a block diagonal sparse matrix from provided matrices.
tril(A[, k, format]) Return the lower triangular portion of a matrix in sparse format Returns the elements on or below the k-th diagonal of the matrix A.
triu(A[, k, format]) Return the upper triangular portion of a matrix in sparse format Returns the elements on or above the k-th diagonal of the matrix A.
bmat(blocks[, format, dtype]) Build a sparse matrix from sparse sub-blocks :Parameters: blocks : array_like Grid of sparse matrices with compatible shapes.
hstack(blocks[, format, dtype]) Stack sparse matrices horizontally (column wise) :Parameters: blocks sequence of sparse matrices with compatible shapes format : str sparse format of the result (e.g.
vstack(blocks[, format, dtype]) Stack sparse matrices vertically (row wise) :Parameters: blocks sequence of sparse matrices with compatible shapes format : str, optional sparse format of the result (e.g.
rand(m, n[, density, format, dtype, ...]) Generate a sparse matrix of the given shape and density with uniformly distributed values.
random(m, n[, density, format, dtype, ...]) Generate a sparse matrix of the given shape and density with randomly distributed values.

Sparse matrix tools:

find(A) Return the indices and values of the nonzero elements of a matrix :Parameters: A : dense or sparse matrix Matrix whose nonzero elements are desired.

Identifying sparse matrices:

issparse(x)
isspmatrix(x)
isspmatrix_csc(x)
isspmatrix_csr(x)
isspmatrix_bsr(x)
isspmatrix_lil(x)
isspmatrix_dok(x)
isspmatrix_coo(x)
isspmatrix_dia(x)

Submodules

csgraph
linalg

Usage information

There are seven available sparse matrix types:

  1. csc_matrix: Compressed Sparse Column format
  2. csr_matrix: Compressed Sparse Row format
  3. bsr_matrix: Block Sparse Row format
  4. lil_matrix: List of Lists format
  5. dok_matrix: Dictionary of Keys format
  6. coo_matrix: COOrdinate format (aka IJV, triplet format)
  7. dia_matrix: DIAgonal format

To construct a matrix efficiently, use either dok_matrix or lil_matrix. The lil_matrix class supports basic slicing and fancy indexing with a similar syntax to NumPy arrays. As illustrated below, the COO format may also be used to efficiently construct matrices.

To perform manipulations such as multiplication or inversion, first convert the matrix to either CSC or CSR format. The lil_matrix format is row-based, so conversion to CSR is efficient, whereas conversion to CSC is less so.

All conversions among the CSR, CSC, and COO formats are efficient, linear-time operations.

Matrix vector product

To do a vector product between a sparse matrix and a vector simply use the matrix dot method, as described in its docstring:

>>> import numpy as np
>>> from scipy.sparse import csr_matrix
>>> A = csr_matrix([[1, 2, 0], [0, 0, 3], [4, 0, 5]])
>>> v = np.array([1, 0, -1])
>>> A.dot(v)
array([ 1, -3, -1], dtype=int64)

Warning

As of NumPy 1.7, np.dot is not aware of sparse matrices, therefore using it will result on unexpected results or errors. The corresponding dense array should be obtained first instead:

>>> np.dot(A.toarray(), v)
array([ 1, -3, -1], dtype=int64)

but then all the performance advantages would be lost.

The CSR format is specially suitable for fast matrix vector products.

Example 1

Construct a 1000x1000 lil_matrix and add some values to it:

>>> from scipy.sparse import lil_matrix
>>> from scipy.sparse.linalg import spsolve
>>> from numpy.linalg import solve, norm
>>> from numpy.random import rand
>>> A = lil_matrix((1000, 1000))
>>> A[0, :100] = rand(100)
>>> A[1, 100:200] = A[0, :100]
>>> A.setdiag(rand(1000))

Now convert it to CSR format and solve A x = b for x:

>>> A = A.tocsr()
>>> b = rand(1000)
>>> x = spsolve(A, b)

Convert it to a dense matrix and solve, and check that the result is the same:

>>> x_ = solve(A.toarray(), b)

Now we can compute norm of the error with:

>>> err = norm(x-x_)
>>> err < 1e-10
True

It should be small :)

Example 2

Construct a matrix in COO format:

>>> from scipy import sparse
>>> from numpy import array
>>> I = array([0,3,1,0])
>>> J = array([0,3,1,2])
>>> V = array([4,5,7,9])
>>> A = sparse.coo_matrix((V,(I,J)),shape=(4,4))

Notice that the indices do not need to be sorted.

Duplicate (i,j) entries are summed when converting to CSR or CSC.

>>> I = array([0,0,1,3,1,0,0])
>>> J = array([0,2,1,3,1,0,0])
>>> V = array([1,1,1,1,1,1,1])
>>> B = sparse.coo_matrix((V,(I,J)),shape=(4,4)).tocsr()

This is useful for constructing finite-element stiffness and mass matrices.

Further Details

CSR column indices are not necessarily sorted. Likewise for CSC row indices. Use the .sorted_indices() and .sort_indices() methods when sorted indices are required (e.g. when passing data to other libraries).