NumPy linalg svd: Python SVD Guide & Code

Q: What is the difference between full_matrices=True and full_matrices=False in numpy.linalg.svd?

When full_matrices=True (the default), NumPy returns complete unitary matrices where U has shape (M, M) and Vh has shape (N, N) . When full_matrices=False , NumPy computes the economy SVD, truncating U to shape (M, K) and Vh to shape (K, N) , where K = min(M, N) . For rectangular research datasets where observations substantially exceed variables ( M ≫ N ), full_matrices=False drastically lowers memory usage and computational execution time without losing any nonzero singular components.

Q: Why does numpy.linalg.svd return Vh instead of V?

NumPy returns Vh representing the conjugate transpose Vᴴ (or the transpose Vᵀ for real matrices) because the matrix factorization is algebraically defined as A = U · Σ · Vᴴ . By returning Vh , the rows of the returned array correspond directly to the right singular basis vectors, and the matrix product U @ np.diag(s) @ Vh can be executed directly without requiring an explicit transpose operation.

Q: How do I reconstruct the original matrix from U, s, and Vh?

Because s is returned as a 1D array of singular values, you must convert it into a diagonal matrix using np.diag(s) before performing matrix multiplication. In Python 3.5+, use the matrix multiplication operator @ : A_reconstructed = U @ np.diag(s) @ Vh . If you used full_matrices=True and M ≠ N , s must be padded into an M × N rectangular matrix before multiplying.

Q: When should I use scipy.linalg.svd instead of numpy.linalg.svd?

You should use scipy.linalg.svd when encountering convergence errors with NumPy’s default divide-and-conquer driver, as SciPy allows setting lapack_driver='gesvd' for traditional, highly robust QR iteration. SciPy also provides overwrite_a=True to conserve memory by reusing the input array buffer and integrates with scipy.sparse.linalg.svds for computing partial SVDs on sparse matrices.

The numpy.linalg.svd function computes the Singular Value Decomposition (SVD) of a real or complex matrix (or a batch stack of matrices) in Python. SVD factorizes an input matrix A into two unitary matrices, U and Vᴴ, and a one-dimensional array of singular values s, satisfying the relation A = U · diag(s) · Vᴴ. As a cornerstone of computational linear algebra, SVD enables quantitative researchers to perform dimensionality reduction, solve rank-deficient least-squares problems, isolate dominant dynamical modes, and filter measurement noise from experimental datasets.

What it is used for in research

Singular value decomposition provides an exact, coordinate-free factorization applicable to any rectangular matrix, making it indispensable across laboratory data analysis, physical modelling, and bioinformatics. Key research applications include:

Dimensionality Reduction and Principal Component Analysis: When applied to mean-centered observational data, SVD extracts principal axes without forming the full covariance matrix, minimizing rounding error and improving numerical stability on ill-conditioned data.
Low-Rank Matrix Approximation: By truncating singular values and singular vectors at rank k, researchers apply the Eckart-Young-Mirsky theorem to produce optimal low-rank reconstructions that denoise multi-channel sensor measurements and compress scientific imagery.
Pseudoinverse and Ill-Conditioned Inversion: Computing the Moore-Penrose pseudoinverse (A⁺ = V · Σ⁺ · Uᴴ) allows researchers to solve ill-posed inverse problems and identify regularized parameter solutions in linear regression models.
Matrix Condition Diagnostics: The ratio of the largest to smallest singular value (κ(A) = σ_max / σ_min) provides the matrix condition number, revealing whether experimental systems or sensitivity matrices are vulnerable to extreme numerical error amplification during matrix inversion.
Multivariate Feature Extraction: SVD furnishes orthogonal subspaces utilized directly in supervised dimensionality reduction workflows such as linear discriminant analysis, preserving core variance across research cohorts while upholding rigorous reproducibility across computational pipelines.

Install or access

The numpy.linalg.svd function is part of the standard NumPy library. Researchers can install NumPy via the official Python Package Index (PyPI) or Conda package repositories:

# Install via pip
pip install numpy

# Install via Conda from conda-forge
conda install -c conda-forge numpy

NumPy links against high-performance compiled linear algebra backends such as OpenBLAS, Intel MKL, or Apple Accelerate via standardized BLAS and LAPACK interfaces. To inspect which compiled BLAS and LAPACK libraries power your local environment, run:

import numpy as np

# Print the active BLAS and LAPACK configuration
np.show_config()

The standard function call requires importing the package and accessing the linear algebra module via import numpy as np and calling np.linalg.svd(...).

Step-by-step task

The following five-step task demonstrates computing the SVD, managing reduced-rank shapes, verifying numerical reconstruction, denoising via rank-k truncation, and executing batch decomposition across multi-trial experimental tensors.

Define the Research Observation Matrix:

Construct a sample rectangular matrix A containing M observations across N features (for example, 4 experimental timepoints across 3 sensor channels):

import numpy as np

# Create an experimental observation matrix (M=4, N=3)
A = np.array([
    [2.5, 1.2, 0.8],
    [1.8, 2.9, 0.4],
    [0.9, 1.1, 3.2],
    [3.1, 2.4, 1.5]
], dtype=np.float64)

print("Original Matrix Shape:", A.shape)

Compute Economy (Reduced) SVD:

Call np.linalg.svd setting full_matrices=False to compute the thin/economy factorization, avoiding unnecessary memory allocation for large research arrays where M ≫ N:

# Compute reduced SVD: U, s, Vh
U, s, Vh = np.linalg.svd(A, full_matrices=False)

print("U shape (Left singular vectors):", U.shape)    # (4, 3)
print("s shape (Singular values):", s.shape)            # (3,)
print("Vh shape (Right singular vectors):", Vh.shape)  # (3, 3)

Verify Exact Matrix Reconstruction:

Reconstruct the original matrix using matrix multiplication (@) and np.diag(s), validating numerical consistency with np.allclose:

# Reconstruct original matrix: A = U @ diag(s) @ Vh
A_reconstructed = U @ np.diag(s) @ Vh

# Verify equivalence within machine precision tolerance
is_reconstructed = np.allclose(A, A_reconstructed)
print("Matrix reconstructed accurately:", is_reconstructed)

Perform Rank-k Truncation (Noise Filtering):

Filter low-amplitude noise by truncating the decomposition to the dominant k singular values (e.g., k = 2 components):

k = 2  # Retain top 2 singular components
U_k = U[:, :k]
s_k = s[:k]
Vh_k = Vh[:k, :]

# Compute low-rank rank-2 approximation
A_rank_k = U_k @ np.diag(s_k) @ Vh_k
print("Rank-2 Approximated Matrix:
", np.round(A_rank_k, 3))

Execute Vectorized Batch Decomposition:

Pass multidimensional stacks of matrices (shape (B, M, N)) to decompose multiple experiment trials simultaneously without Python loops:

# Stack of 5 experimental trials with dimension (4, 3)
batch_trials = np.random.default_rng(42).standard_normal(size=(5, 4, 3))

# Vectorized batch SVD
U_batch, s_batch, Vh_batch = np.linalg.svd(batch_trials, full_matrices=False)

print("Batch U shape:", U_batch.shape)    # (5, 4, 3)
print("Batch s shape:", s_batch.shape)    # (5, 3)
print("Batch Vh shape:", Vh_batch.shape)  # (5, 3, 3)

Output and how to read it

The call U, s, Vh = np.linalg.svd(a) returns three outputs with distinct mathematical and computational properties:

U (Left Singular Vectors): An array of shape (..., M, M) when full_matrices=True, or (..., M, K) when full_matrices=False (where K = min(M, N)). The columns of U form an orthonormal basis spanning the column space of A. In physical modelling, these represent spatial modes or system state coordinates.
s (Singular Values): A one-dimensional array of shape (..., K) containing real, non-negative singular values sorted in strictly descending order (σ₁ ≥ σ₂ ≥ … ≥ σ_K ≥ 0). Each singular value corresponds to the square root of an eigenvalue of AᴴA, quantifying the variance or spectral energy captured along each principal component.
Vh (Right Singular Vectors Transposed): An array of shape (..., N, N) when full_matrices=True, or (..., K, N) when full_matrices=False. Crucial reading note: Vh is returned as the conjugate transpose Vᴴ (or Vᵀ for real matrices), not V. The rows of Vh are the orthonormal basis vectors for the row space of A. To project a sample row vector x into the singular coordinate space, researchers multiply by Vh.T (or Vh.conj().T for complex data), rather than Vh.

If singular vectors are not needed, setting compute_uv=False returns only the 1D array s, reducing memory usage and computational execution time.

Troubleshooting

Common issues encountered when executing numpy.linalg.svd in research workflows include:

LinAlgError: SVD did not converge: This error originates from the underlying LAPACK divide-and-conquer driver (DGESDD) when the singular value subproblem fails to converge. Common causes include unhandled NaN or Inf values, severe matrix ill-conditioning, or floating-point denormals. Verify data cleanliness using np.isfinite(A).all(). If ill-conditioning persists, compute SVD using SciPy’s alternative driver scipy.linalg.svd(A, lapack_driver='gesvd').
Out-of-Memory Errors with Large Rectangular Matrices: When factorizing a tall matrix where M = 50,000 and N = 200, leaving full_matrices=True causes NumPy to construct an M × M matrix U (50,000 × 50,000 float64 elements), requiring approximately 20 GB of memory. Always specify full_matrices=False on tall or wide matrices, or switch to iterative sparse solvers via scipy-sparse routines.
Incorrect Reconstruction from Misaligned Transposition: Attempting to reconstruct A using U @ np.diag(s) @ Vh.T produces either a dimension mismatch or invalid results because Vh is already transposed. The correct algebraic reconstruction is U @ np.diag(s) @ Vh.
Dimension Error (ValueError: array must be at least two-dimensional): Passing a 1D array to np.linalg.svd triggers a dimensional error. Reshape 1D observation vectors into 2D arrays (e.g., x[:, np.newaxis] or x.reshape(-1, 1)) prior to factorization.

Free alternatives

Depending on matrix sparsity, hardware acceleration requirements, and memory constraints, researchers can leverage several open-source alternatives to numpy.linalg.svd:

Tool / Library	Interface	Algorithm / Driver	Sparse Support	Best Research Use Case
NumPy (`numpy.linalg.svd`)	Python	LAPACK `_gesdd`	No (dense only)	General in-memory matrix decomposition, batch array processing, and desktop data analysis.
SciPy (`scipy.linalg.svd`)	Python	LAPACK `gesdd` or `gesvd`	No (dense only)	Advanced scientific pipelines requiring driver selection (`gesvd` for convergence stability) and in-place memory overwriting.
SciPy Sparse (`scipy.sparse.linalg.svds`)	Python	ARPACK / PROPACK	Yes (CSR, CSC)	Large-scale sparse research problems (e.g., single-cell genomics, document term matrices, citation networks) seeking top k singular values.
PyTorch (`torch.linalg.svd`)	Python / C++	cuSOLVER / LAPACK	Limited	GPU-accelerated tensor decompositions, deep learning latent representations, and automatic differentiation pipelines.
GNU Octave (`svd`)	Octave / MATLAB syntax	LAPACK `dgesvd` / `dgesdd`	Yes	Algorithmic prototyping, legacy research script migration, and interactive desktop matrix investigations.

FAQ

What is the difference between full_matrices=True and full_matrices=False in numpy.linalg.svd?

When full_matrices=True (the default), NumPy returns complete unitary matrices where U has shape (M, M) and Vh has shape (N, N). When full_matrices=False, NumPy computes the economy SVD, truncating U to shape (M, K) and Vh to shape (K, N), where K = min(M, N). For rectangular research datasets where observations substantially exceed variables (M ≫ N), full_matrices=False drastically lowers memory usage and computational execution time without losing any nonzero singular components.

Why does numpy.linalg.svd return Vh instead of V?

NumPy returns Vh representing the conjugate transpose Vᴴ (or the transpose Vᵀ for real matrices) because the matrix factorization is algebraically defined as A = U · Σ · Vᴴ. By returning Vh, the rows of the returned array correspond directly to the right singular basis vectors, and the matrix product U @ np.diag(s) @ Vh can be executed directly without requiring an explicit transpose operation.

How do I reconstruct the original matrix from U, s, and Vh?

Because s is returned as a 1D array of singular values, you must convert it into a diagonal matrix using np.diag(s) before performing matrix multiplication. In Python 3.5+, use the matrix multiplication operator @: A_reconstructed = U @ np.diag(s) @ Vh. If you used full_matrices=True and M ≠ N, s must be padded into an M × N rectangular matrix before multiplying.

When should I use scipy.linalg.svd instead of numpy.linalg.svd?

You should use scipy.linalg.svd when encountering convergence errors with NumPy’s default divide-and-conquer driver, as SciPy allows setting lapack_driver='gesvd' for traditional, highly robust QR iteration. SciPy also provides overwrite_a=True to conserve memory by reusing the input array buffer and integrates with scipy.sparse.linalg.svds for computing partial SVDs on sparse matrices.

How does numpy.linalg.svd relate to Principal Component Analysis (PCA)?

SVD is the primary computational engine underlying PCA. When an observation matrix X (samples by variables) is centered by subtracting column means, computing U, s, Vh = np.linalg.svd(X, full_matrices=False) yields right singular vectors Vh.T that are mathematically identical to the principal component loading vectors (eigenvectors of the covariance matrix). The principal component scores are obtained directly by computing X @ Vh.T or U * s.