low_rank_toolbox.randomized
Randomized algorithms for low-rank matrix approximation.
This submodule implements state-of-the-art randomized algorithms for computing low-rank approximations of large matrices, based on the seminal work by Halko, Martinsson, and Tropp (2010).
Functions
rangefinder : Approximate the range of a dense matrix with a target rank adaptive_rangefinder : Approximate the range of a dense matrix with a target tolerance randomized_svd : Implementation of the randomized SVD algorithm with a prescribed rank adaptive_randomized_svd : Adaptive randomized SVD with tolerance-based rank selection generalized_nystrom : Generalized Nyström method for low-rank approximation
References
Author: Benjamin Carrel, University of Geneva
- low_rank_toolbox.randomized.adaptive_randomized_svd(X, tol=1e-06, failure_prob=1e-06, max_rank=None, seed=1234, **extra_data)[source]
Adaptive randomized SVD algorithm.
- Reference:
“Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions”, Halko, Martinsson and Tropp 2010.
The adaptive randomized SVD computes an approximate SVD of a matrix X with the following steps: 1. Estimate the range with Q = adaptive_rangefinder(X, tol, failure_prob, seed) 2. Form the smaller matrix C = Q^H X 3. Compute the SVD of C: C = U_C S V^H 4. Form the approximate SVD of X: X_approx = Q U_C S V^H
- Parameters:
X (SVD | LowRankMatrix | ndarray | LinearOperator) – Matrix or operator to be decomposed
tol (float, optional) – Target approximation error (default: 1e-6)
failure_prob (float, optional) – Failure probability for the approximation error, must be in (0, 1) (default: 1e-6)
max_rank (int, optional) – Maximum rank for the approximation (default: None, no limit)
seed (int, optional) – Random seed for reproducibility (default: 1234)
**extra_data (dict) – Additional arguments passed to SVD constructor
- Returns:
Near-optimal approximation of X with ||X - X_approx||_F <= tol with probability >= 1 - failure_prob
- Return type:
- low_rank_toolbox.randomized.adaptive_rangefinder(A, tol=1e-06, failure_prob=1e-06, seed=1234)[source]
The adaptive (randomized) rangefinder method. The tolerance is the error made by the approximation space ||A - QQ^H A||_F <= tol The failure probability is the probability that the error is larger than tol
- Reference:
“Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions”, Halko, Martinsson and Tropp 2010.
NOTE: For efficiency, the method performs computations in blocks of size r (defined by the failure probability). Blocking allows to use BLAS3 operations.
- Parameters:
- Returns:
Q – The sketched matrix with orthonormal columns.
- Return type:
ndarray
- low_rank_toolbox.randomized.generalized_nystrom(X, r, oversampling_params=(10, 15), epsilon=None, seed=1234, **extra_data)[source]
Generalized Nyström method
- Reference:
“Fast and stable randomized low-rank matrix approximation” Nakatsukasa, 2019
- Approximation with formula
X ~= X J (K^T X J)^{dagger} K^T X
- Parameters:
X (LinearOperator) – Matrix to approximate (real-valued only, complex matrices not supported)
r (int) – Rank of approximation, must be positive
oversampling_params (tuple, optional) – Oversampling parameters (p1, p2) for the two sketch matrices (default: (10, 15))
epsilon (float, optional) – When given, perform stable GN with epsilon-truncation for SVD (default: None)
seed (int, optional) – Random seed for reproducibility (default: 1234)
**extra_data (dict) – Additional arguments passed to QuasiSVD constructor
- Returns:
Near optimal best rank r approximation of X in QuasiSVD format
- Return type:
Notes
This method returns a QuasiSVD (not SVD) because the middle matrix S is typically inverted, making it non-diagonal. Convert to SVD if needed:
result = generalized_nystrom(X, r).to_svd()
- low_rank_toolbox.randomized.randomized_svd(X, r, p=10, q=0, truncate=True, seed=1234, **extra_data)[source]
Randomized SVD algorithm.
- Reference:
“Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions”, Halko, Martinsson and Tropp 2010.
The randomized SVD computes an approximate SVD of a matrix X with the following steps: 1. Estimate the range with Q = rangefinder(X, r, p, q, seed) 2. Form the smaller matrix C = Q^H X 3. Compute the SVD of C: C = U_C S V^H 4. Form the approximate SVD of X: X_approx = Q U_C S V^H
- Parameters:
X (SVD | LowRankMatrix | ndarray | LinearOperator) – Matrix or operator to be decomposed
r (int) – Rank of the decomposition, must be positive
p (int, optional) – Oversampling parameter (default: 10)
q (int, optional) – Number of power iterations (default: 0)
truncate (bool, optional) – If True, truncate the SVD to rank r (default: True)
seed (int, optional) – Random seed for reproducibility (default: 1234)
**extra_data (dict) – Additional arguments passed to SVD constructor. If ‘Omega’ is provided, use it as the sketching matrix.
- Returns:
Near-optimal best rank r approximation of X
- Return type:
- low_rank_toolbox.randomized.rangefinder(A, r, p=5, q=0, seed=1234, **extra_args)[source]
The (randomized) rangefinder method, also called HMT method.
- Reference:
“Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions”, Halko, Martinsson and Tropp 2010.
- The mathematical formulation is:
Q, _ = qr((A A^H)^q A Omega, mode=’economic’)
where Omega is a Gaussian random matrix of shape (n, r + p).
When A is a low-rank matrix, the rangefinder provides a good approximation of the range of A with high probability.
- Parameters:
A (LinearOperator) – The matrix to sketch.
r (int) – The target rank.
p (int, optional) – The number of over-sampling (default: 5).
q (int, optional) – Number of power iterations (default: 0).
seed (int, optional) – The seed for the random number generator (default: 1234).
**extra_args (dict) – Additional arguments. If ‘Omega’ is provided, use it as the sketching matrix.
- Returns:
Q – Estimation of the range of A.
- Return type:
ndarray
Modules
|
Generalized Nyström method |
|
Randomized SVD algorithm. |
|
The (randomized) rangefinder method, also called HMT method. |