scipy.sparse.linalg.svds#
- scipy.sparse.linalg.svds(A, k=6, ncv=None, tol=0, which='LM', v0=None, maxiter=None, return_singular_vectors=True, solver='arpack', random_state=None, options=None)[source]#
Partial singular value decomposition of a sparse matrix.
Compute the largest or smallest k singular values and corresponding singular vectors of a sparse matrix A. The order in which the singular values are returned is not guaranteed.
In the descriptions below, let
M, N = A.shape
.- Parameters
- Andarray, sparse matrix, or LinearOperator
Matrix to decompose of a floating point numeric dtype.
- kint, default: 6
Number of singular values and singular vectors to compute. Must satisfy
1 <= k <= kmax
, wherekmax=min(M, N)
forsolver='propack'
andkmax=min(M, N) - 1
otherwise.- ncvint, optional
When
solver='arpack'
, this is the number of Lanczos vectors generated. See ‘arpack’ for details. Whensolver='lobpcg'
orsolver='propack'
, this parameter is ignored.- tolfloat, optional
Tolerance for singular values. Zero (default) means machine precision.
- which{‘LM’, ‘SM’}
Which k singular values to find: either the largest magnitude (‘LM’) or smallest magnitude (‘SM’) singular values.
- v0ndarray, optional
The starting vector for iteration; see method-specific documentation (‘arpack’, ‘lobpcg’), or ‘propack’ for details.
- maxiterint, optional
Maximum number of iterations; see method-specific documentation (‘arpack’, ‘lobpcg’), or ‘propack’ for details.
- return_singular_vectors{True, False, “u”, “vh”}
Singular values are always computed and returned; this parameter controls the computation and return of singular vectors.
True
: return singular vectors.False
: do not return singular vectors."u"
: ifM <= N
, compute only the left singular vectors and returnNone
for the right singular vectors. Otherwise, compute all singular vectors."vh"
: ifM > N
, compute only the right singular vectors and returnNone
for the left singular vectors. Otherwise, compute all singular vectors.
If
solver='propack'
, the option is respected regardless of the matrix shape.- solver{‘arpack’, ‘propack’, ‘lobpcg’}, optional
The solver used. ‘arpack’, ‘lobpcg’, and ‘propack’ are supported. Default: ‘arpack’.
- random_state{None, int,
numpy.random.Generator
, numpy.random.RandomState
}, optionalPseudorandom number generator state used to generate resamples.
If random_state is
None
(or np.random), thenumpy.random.RandomState
singleton is used. If random_state is an int, a newRandomState
instance is used, seeded with random_state. If random_state is already aGenerator
orRandomState
instance then that instance is used.- optionsdict, optional
A dictionary of solver-specific options. No solver-specific options are currently supported; this parameter is reserved for future use.
- Returns
- undarray, shape=(M, k)
Unitary matrix having left singular vectors as columns.
- sndarray, shape=(k,)
The singular values.
- vhndarray, shape=(k, N)
Unitary matrix having right singular vectors as rows.
Notes
This is a naive implementation using ARPACK or LOBPCG as an eigensolver on
A.conj().T @ A
orA @ A.conj().T
, depending on which one is more efficient, followed by the Rayleigh-Ritz method as postprocessing; see https://w.wiki/4zmsAlternatively, the PROPACK solver can be called.
form="array"
Choices of the input matrix
A
numeric dtype may be limited. Onlysolver="lobpcg"
supports all floating point dtypes real: ‘np.single’, ‘np.double’, ‘np.longdouble’ and complex: ‘np.csingle’, ‘np.cdouble’, ‘np.clongdouble’. Thesolver="arpack"
supports only ‘np.single’, ‘np.double’, and ‘np.cdouble’.Examples
Construct a matrix
A
from singular values and vectors.>>> from scipy.stats import ortho_group >>> from scipy.sparse import csc_matrix, diags >>> from scipy.sparse.linalg import svds >>> rng = np.random.default_rng() >>> orthogonal = csc_matrix(ortho_group.rvs(10, random_state=rng)) >>> s = [0.0001, 0.001, 3, 4, 5] # singular values >>> u = orthogonal[:, :5] # left singular vectors >>> vT = orthogonal[:, 5:].T # right singular vectors >>> A = u @ diags(s) @ vT
With only three singular values/vectors, the SVD approximates the original matrix.
>>> u2, s2, vT2 = svds(A, k=3) >>> A2 = u2 @ np.diag(s2) @ vT2 >>> np.allclose(A2, A.toarray(), atol=1e-3) True
With all five singular values/vectors, we can reproduce the original matrix.
>>> u3, s3, vT3 = svds(A, k=5) >>> A3 = u3 @ np.diag(s3) @ vT3 >>> np.allclose(A3, A.toarray()) True
The singular values match the expected singular values, and the singular vectors are as expected up to a difference in sign.
>>> (np.allclose(s3, s) and ... np.allclose(np.abs(u3), np.abs(u.toarray())) and ... np.allclose(np.abs(vT3), np.abs(vT.toarray()))) True
The singular vectors are also orthogonal. >>> (np.allclose(u3.T @ u3, np.eye(5)) and … np.allclose(vT3 @ vT3.T, np.eye(5))) True