Chapter 11 Singular Value Decomposition and Its Applications

The history of the singular value decomposition (SVD) shows how theory, computation, and applications push one another forward.

In 1873 the Italian mathematician Eugenio Beltrami studied real bilinear forms; in 1874 the French mathematician Camille Jordan obtained related results. Erhard Schmidt studied integral operators and approximation problems in 1907, Carl Eckart and Gale Young established the result on the best low-rank approximation of a matrix in 1936, and Mirsky extended the approximation theory to unitarily invariant norms in 1960. These results dealt separately with the existence of the decomposition, with operators, and with best approximation, and together they form a theory that developed step by step; without historical evidence, one should not claim that these scholars were unaware of each other or that their results were repeatedly forgotten. History and original sources

Beyond the theory, computational methods are just as essential. The SVD can be built by finding the eigenvalues of $\mathbf{A}^{*}\mathbf{A}$, but forming this Gram matrix directly may lose relative accuracy in the small singular values. For a matrix of full column rank, $\kappa_2(\mathbf{A}^{*}\mathbf{A})=\kappa_2(\mathbf{A})^2$; this is a matter of numerical conditioning, and it does not mean that every input error is amplified quadratically.

Algorithms did not first appear in 1970. Jacobi-type methods already existed in the 1950s; Golub and Kahan developed the bidiagonalization approach in 1965, and Golub and Reinsch went on in 1970 to give algorithms for the SVD and for least squares. These methods avoid making the explicit formation of the Gram matrix a required step, and they laid an important foundation for reliable numerical computation. Hestenes 1958, Golub–Kahan 1965

The teaching power of the SVD comes from a clear geometric picture: once orthonormal bases are chosen for the input and the output, $\mathbf{A}\mathbf{v}_i=\sigma_i\mathbf{u}_i$ ($i<r$), and the remaining input basis vectors fall into the null space. The SVD provides compatible orthogonal coordinates for the four fundamental subspaces; the nonzero singular values give the stretching factors—the SVD does not turn every vector into another unit vector.

Applications have also carried these ideas beyond the classroom. The Netflix Prize spurred regularized low-rank factor models for ratings with missing entries; these models share the low-rank viewpoint of the SVD, but they do not fill in the missing values arbitrarily and then run an ordinary SVD, nor are they the same as exact nuclear-norm completion. PCA, the pseudoinverse, and data compression have made singular values a common language of data analysis.

This line of development reminds us that theoretical correctness, computational feasibility, and the needs of applications together determine how much power a mathematical tool can exert.

Chapter Structure and Learning Objectives¶

This chapter follows the path of theory, computation, and applications, breaking the SVD into material that can be understood section by section.

§11.1 first answers the most fundamental question: what exactly does this theorem say, and why does it hold for every matrix? We rigorously prove the existence of the decomposition $\mathbf{A} = \mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{*}$ and examine its meaning again and again from three angles: algebra, subspaces, and geometry. $\mathbf{A}^{*}\mathbf{A}$ is the bridge connecting the SVD with the spectral theorem of Chapter 9; under the SVD, each of the four fundamental subspaces receives an orthonormal basis; and the geometric picture of “rotation—stretch—rotation” makes the action of a matrix clear at a glance. By the end of this section you will hold a key that opens the door to every application that follows.

In §11.2 this key faces its first test: in what sense is the low-rank approximation given by the truncated SVD “the best”? This is precisely the question Eckart and Young studied in 1936. Three matrix norms—the spectral norm, the Frobenius norm, and the nuclear norm—are all completely determined by the singular values, and the theorem’s answer is clean and decisive: among all matrices of the same rank bound, the truncated SVD has the smallest error, and the discarded singular values tell exactly how much is lost. This is the theoretical peak of the chapter, and the three application sections that follow are all its echoes.

§11.3 brings this optimality into the setting of data science. The “directions that explain the most variation in the data” sought by principal component analysis (PCA) are, in the language of the SVD, the right singular vectors; numerically, performing the SVD directly on the data matrix rather than on the covariance matrix avoids the very problem Golub worked to solve: the condition number being quietly squared during the computation. A complete demonstration on the Iris dataset takes this procedure from formulas into practice.

§11.4 changes the question from “finding directions” to “solving equations”: rectangular or singular matrices have no inverse, but the Moore–Penrose pseudoinverse, through the SVD, gracefully extends the meaning of “inversion” and provides the geometrically most natural solution to least-squares problems. As you read this section, heed the warning that “the reciprocals of small singular values blow up”—it leads to ridge regression regularization, and it also gives you your first reliable tool for handling any ill-conditioned linear problem in the future.

§11.5 then carries the SVD to an unexpected shore. Quantum entanglement—what Einstein called “spooky action at a distance”—is, mathematically, exactly a singular value problem for the coefficient matrix of a bipartite pure state. The decomposition Schmidt built for integral operators in 1907 became, after his death, a core language of quantum information theory, which is probably one of the most surprising “late applications” in the history of mathematics. This section is a mathematical deepening of the quantum mechanics framework of Chapter 10, and also the book’s final statement of the theme of “the deep unity of linear algebra and the physical world.”

Having read the whole chapter, you will find that singular values are not isolated numbers but a ranking of the “signal strength” of a matrix in each direction—truncate them, and you are compressing; invert them, and you are solving; measure them, and you are quantifying entanglement. The same set of numbers ties together three worlds—mathematics, data, and physics—and also the long development of theory and applications from 1873 to the present day.

11.1 Definition and Proof of the Singular Value Decomposition¶

From Diagonalization to the SVD: Breaking the Square-Matrix Barrier¶

In Chapter 8 we witnessed the power of diagonalization: for a diagonalizable square matrix $\mathbf{A} \in \mathbb{K}^{n \times n}$, there exists an invertible matrix $\mathbf{P}$ such that

where $\boldsymbol{\Lambda}$ is a diagonal matrix whose diagonal entries are the eigenvalues of $\mathbf{A}$. This decomposition reveals the essential structure of the matrix—in a suitable basis, the linear transformation is nothing more than independent stretching along each direction.

Chapter 9 went further: for a symmetric matrix $\mathbf{A} = \mathbf{A}^{\top}$ (real field) or a Hermitian matrix $\mathbf{A} = \mathbf{A}^{*}$ (complex field), the situation is even better—we can choose an orthogonal matrix $\mathbf{Q}$ (or a unitary matrix $\mathbf{U}$) to diagonalize it:

This is the spectral theorem—the orthogonal diagonalization of symmetric/Hermitian matrices. The beauty of orthogonal diagonalization is that the transformations $\mathbf{Q}$ and $\mathbf{Q}^{\top}$ are both geometric rotations (they preserve lengths and angles), so the whole linear transformation can be understood as “rotate—stretch—rotate back.”

But the real world presents a harsh fact: the vast majority of important matrices are neither square nor symmetric.

Consider the following practical scenarios:

• Image matrices: for a $1920 \times 1080$ picture, the matrix of pixel values is rectangular

• Data matrices: recording the ratings of 1000 customers on 50 products gives a $1000 \times 50$ matrix, far from square

• Document-term matrices: 10000 documents containing 5000 terms give a $10000 \times 5000$ matrix

Clearly these matrices cannot use the eigendecomposition directly. We need a more general decomposition that can handle an arbitrary $m \times n$ matrix.

The Core Idea of the SVD: Starting from $\mathbf{A}^{*}\mathbf{A}$¶

The key observation is this: even if $\mathbf{A}$ is not square, $\mathbf{A}^{*}\mathbf{A}$ and $\mathbf{A}\mathbf{A}^{*}$ are symmetric/Hermitian square matrices!

Let $\mathbf{A} \in \mathbb{C}^{m \times n}$. Then:

• $\mathbf{A}^{*}\mathbf{A} \in \mathbb{C}^{n \times n}$ and $(\mathbf{A}^{*}\mathbf{A})^{*} = \mathbf{A}^{*}\mathbf{A}$ (a Hermitian matrix)

• $\mathbf{A}\mathbf{A}^{*} \in \mathbb{C}^{m \times m}$ and $(\mathbf{A}\mathbf{A}^{*})^{*} = \mathbf{A}\mathbf{A}^{*}$ (a Hermitian matrix)

More importantly, these two matrices are also positive semidefinite (recall Chapter 9): for any vector $\mathbf{x}$,

By the spectral theorem of Chapter 9, $\mathbf{A}^{*}\mathbf{A}$ and $\mathbf{A}\mathbf{A}^{*}$ can both be orthogonally diagonalized, and all their eigenvalues are nonnegative real numbers. This is the theoretical foundation of the SVD.

Historical note. The work on bilinear forms by the Italian mathematician Beltrami in 1873 and the French mathematician Camille Jordan in 1874 is the early source of the SVD. Schmidt studied integral operators and approximation theory in 1907; Eckart and Young gave the result on the best low-rank approximation of a matrix in 1936, and Mirsky extended it to unitarily invariant norms in 1960. These works developed, respectively, the theory of the decomposition, of operators, and of approximation; they should not be lumped together as a single theorem that was repeatedly forgotten.

11.1.1 Definition of the SVD: The Real and Complex Cases¶

Uniqueness theorem. For a given matrix $\mathbf{A}$, its singular values are uniquely determined (counting multiplicity). The singular vectors, however, are not unique.

Three Forms of the SVD¶

Depending on practical needs, the SVD comes in three different forms:

1. Full SVD

This is the most standard form, and the direct statement of the theorem above:

where $\mathbf{U} \in \mathbb{C}^{m \times m}$, $\boldsymbol{\Sigma} \in \mathbb{R}^{m \times n}$, and $\mathbf{V} \in \mathbb{C}^{n \times n}$ are all full matrices.

Matrix dimensions (for $m=4, n=3$):

2. Reduced (Thin) SVD

When $m > n$, the lower part of $\boldsymbol{\Sigma}$ has many all-zero rows, and the corresponding last $m-n$ columns of $\mathbf{U}$ do not actually take part in the matrix multiplication. We can omit these redundant parts:

where:

• $\mathbf{U}_r \in \mathbb{C}^{m \times n}$ keeps only the first $n$ columns

• $\boldsymbol{\Sigma}_r \in \mathbb{R}^{n \times n}$ is a square matrix

• $\mathbf{V} \in \mathbb{C}^{n \times n}$ stays the same

Matrix dimensions (for $m=4, n=3$):

3. Truncated SVD

When the rank $r$ of the matrix is much smaller than $\min(m,n)$, or when we care only about the $k$ largest singular values, we can truncate further:

where:

• $\mathbf{U}_k \in \mathbb{C}^{m \times k}$ contains only the first $k$ left singular vectors

• $\boldsymbol{\Sigma}_k \in \mathbb{R}^{k \times k}$ contains only the first $k$ singular values

• $\mathbf{V}_k \in \mathbb{C}^{n \times k}$ contains only the first $k$ right singular vectors

Matrix dimensions (for $m=4, n=3, k=2$):

This gives a rank-$k$ approximation of $\mathbf{A}$, which is the basis of the Eckart–Young theorem (§11.2).

Outer-Product Expansion (Dyadic Decomposition)¶

Rewriting the truncated SVD as a sum of outer products gives a very illuminating decomposition of the matrix:

where $\mathbf{u}_i \mathbf{v}_i^{*}$ is a rank-one matrix (also called a dyad or an outer product).

Geometric interpretation. Every matrix $\mathbf{A}$ can be decomposed into a weighted sum of $r$ rank-one matrices, with the singular values $\sigma_i$ as weights. These rank-one components are arranged in decreasing order of importance (size of the singular value), so:

• $\sigma_0 \mathbf{u}_0 \mathbf{v}_0^{*}$ is the “most important” component of $\mathbf{A}$

• $\sigma_1 \mathbf{u}_1 \mathbf{v}_1^{*}$ is the next most important component

• and so on

11.1.2 The SVD and the Four Fundamental Subspaces¶

Recall from Chapter 6 that every matrix $\mathbf{A} \in \mathbb{C}^{m \times n}$ defines four fundamental subspaces:

1. Row space: $\text{Row}(\mathbf{A}) = \text{span}\{\text{row}_0(\mathbf{A}), \ldots, \text{row}_{m-1}(\mathbf{A})\} \subseteq \mathbb{C}^n$

2. Null space: $\ker(\mathbf{A}) = \{\mathbf{x} \in \mathbb{C}^n : \mathbf{A}\mathbf{x} = \mathbf{0}\} \subseteq \mathbb{C}^n$

3. Image (column space): $\text{Col}(\mathbf{A}) = \text{span}\{\text{col}_0(\mathbf{A}), \ldots, \text{col}_{n-1}(\mathbf{A})\} \subseteq \mathbb{C}^m$

4. Left null space: $\ker(\mathbf{A}^{*}) = \{\mathbf{y} \in \mathbb{C}^m : \mathbf{A}^{*}\mathbf{y} = \mathbf{0}\} \subseteq \mathbb{C}^m$

Here $\text{Row}(\mathbf{A})$ is spanned by row vectors written horizontally; to place it in the input space of upright (column) vectors we must take the conjugate transpose, obtaining $\text{Col}(\mathbf{A}^{*})$—an ordinary transpose is not enough. The corresponding two orthogonal decompositions are

What makes the SVD so convenient is that it provides explicit orthonormal bases for these four fundamental subspaces!

The Orthogonal Bases Provided by the SVD¶

Let $\mathbf{A} = \mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{*}$ be the SVD of a matrix of rank $r$, and let $\mathbf{v}_i$ be the $i$-th column of $\mathbf{V}$. Then:

11.1.3 Geometric Intuition: From the Unit Sphere to an Ellipsoid¶

The most beautiful geometric interpretation of the SVD is this: every linear transformation maps the unit sphere to an ellipsoid, and the SVD describes precisely the directions and lengths of the principal axes of this ellipsoid.

The Geometric Effect of a Linear Transformation¶

Consider a linear transformation $\mathbf{A}: \mathbb{R}^n \to \mathbb{R}^m$. The unit sphere in the domain is

and its image under $\mathbf{A}$ is

Question. What is the geometric shape of this image set?

Answer. It is an ellipsoid in $\mathbb{R}^m$! More precisely, it is an $r$-dimensional ellipsoid (where $r = \text{rank}(\mathbf{A})$), embedded in an $r$-dimensional subspace of the $m$-dimensional space (namely the image $\text{Col}(\mathbf{A})$).

The Geometric Interpretation of the SVD¶

Let $\mathbf{A} = \mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{\top}$ be an SVD. Then for any unit vector $\mathbf{x} \in S^{n-1}$,

We understand this in three steps:

Step 1. Rotation ($\mathbf{V}^{\top}$)

$\mathbf{V}^{\top}$ is an orthogonal matrix, corresponding to a rigid rotation (or reflection). It rotates $\mathbf{x}$ from the original coordinate system into the coordinate system of the right singular vectors $\{\mathbf{v}_0, \mathbf{v}_1, \ldots, \mathbf{v}_{n-1}\}$:

Since $\mathbf{V}$ is an orthogonal matrix, $\|\mathbf{y}\| = \|\mathbf{x}\| = 1$, so $\mathbf{y}$ is still on the unit sphere.

Step 2. Anisotropic stretching ($\boldsymbol{\Sigma}$)

$\boldsymbol{\Sigma}$ is a diagonal matrix that performs independent stretching along each coordinate direction:

This operation stretches the unit sphere into an ellipsoid:

• stretching by a factor of $\sigma_0$ in the direction $\mathbf{e}_0$ (the largest singular value, the longest principal axis)

• stretching by a factor of $\sigma_1$ in the direction $\mathbf{e}_1$

• ...

• stretching by a factor of $\sigma_{r-1}$ in the direction $\mathbf{e}_{r-1}$ (the smallest nonzero singular value)

• compressing to 0 in the directions $\mathbf{e}_r, \ldots, \mathbf{e}_{n-1}$ (corresponding to the null space)

Step 3. Rotation ($\mathbf{U}$)

$\mathbf{U}$ rotates the ellipsoid from the standard coordinate system into the coordinate system of the left singular vectors $\{\mathbf{u}_0, \mathbf{u}_1, \ldots, \mathbf{u}_{m-1}\}$:

In the end, the directions of the principal axes of the ellipsoid are the left singular vectors $\mathbf{u}_0, \mathbf{u}_1, \ldots, \mathbf{u}_{r-1}$, and the lengths of the principal axes are the corresponding singular values $\sigma_0, \sigma_1, \ldots, \sigma_{r-1}$.

The Variational Viewpoint: Directions of Maximal Stretching¶

Another way to understand the SVD is through optimization: the singular vectors are the directions in which $\mathbf{A}$ produces the largest/smallest stretching.

Geometric meaning:

• $\mathbf{v}_0$ is the direction that $\mathbf{A}$ stretches the most

• when $r>0$, among the unit vectors in $\ker(\mathbf{A})^\perp$, $\mathbf{v}_{r-1}$ attains the smallest stretching $\sigma_{r-1}$

• $\mathbf{v}_r, \ldots, \mathbf{v}_{n-1}$ are compressed to zero by $\mathbf{A}$

This variational viewpoint plays a central role in principal component analysis (PCA, §11.3)—PCA is essentially a search for the directions of greatest variation in the data, and these directions are exactly the right singular vectors of the data matrix!

11.2 Singular Values and Matrix Norms¶

In §11.1 we established the basic theory of the SVD and saw how the singular values describe the geometric structure by which a matrix maps the unit sphere to an ellipsoid. This section explores the relationship between singular values and matrix norms in depth and proves one of the most important applications of the SVD—the Eckart–Young theorem: the truncated SVD gives the best low-rank approximation of a matrix.

From Vector Norms to Matrix Norms¶

In Chapter 2 we became familiar with the concept of the $l_2$ (Euclidean) norm of a vector, which abstracts the intuitive notion of length. In fact, $l_2$ is just one special case of a whole family of $l_p$ norms. We introduce them here together, since they will be used later in the analysis of matrix norms.

Different $l_p$ norms define different “balls” (the set of all vectors with $\|\mathbf{v}\|_p \leq 1$). In the plane, the shapes of these unit balls reveal the geometric essence of each norm:

• $p=1$: a diamond (Manhattan distance)

• $p=2$: the standard circle (Euclidean distance)

• $p \to \infty$: a square (maximum-coordinate distance)

The figure drawn by the code below clearly shows the geometric process by which, as $p$ increases, the unit ball gradually “inflates” from a diamond into a circle and finally approaches a square. This intuition is essential for understanding the properties of matrix norms later on.

Question. How should we define the “size” of a matrix?

There are two most natural lines of generalization:

1. Induced norm: regard the matrix as an operator and measure its “maximum amplification factor” on vectors

2. Entrywise norm: regard the matrix as a high-dimensional vector and define a norm directly on all of its entries

These two lines lead to the two great families of matrix norms. The spectral norm, the Frobenius norm, and the nuclear norm introduced below are all unitarily invariant norms and can be described uniformly by the singular values; not every induced norm or entrywise norm has this property.

11.2.1 Characterizing Matrix Norms by Singular Values¶

The Spectral Norm (Operator Norm)¶

Definition. The spectral norm (also called the operator norm or the 2-norm) of a matrix $\mathbf{A} \in \mathbb{C}^{m \times n}$ is defined as

Geometric meaning. The spectral norm measures the maximum stretching ratio of unit vectors under the matrix $\mathbf{A}$.

Geometric intuition. Recall from §11.1.3 that the matrix $\mathbf{A}$ maps the unit sphere to an ellipsoid, and the length of the longest principal axis of the ellipsoid is precisely $\sigma_0$. The spectral norm captures exactly this “maximum stretching.”

The Frobenius Norm¶

Geometric meaning. The Frobenius norm measures the “total energy” of a matrix—it combines the contributions of all singular values (the square root of the sum of squares). This will play a key role in the low-rank approximation that follows.

The Nuclear Norm (Trace Norm)¶

Applications. The nuclear norm plays a central role in matrix completion and low-rank optimization problems. It is a convex relaxation of the rank of a matrix—the rank function $\text{rank}(\mathbf{A})$ is nonconvex, whereas the nuclear norm $\|\mathbf{A}\|_*$ is convex and can serve as a surrogate for the rank in optimization.

The following is an idealized noise-free matrix completion model. The Netflix Prize competition spurred regularized low-rank factor models for data with missing entries; its methods should not be equated directly with the exact nuclear-norm completion constraint below:

Relations with other norms:

Summary of the Geometric and Physical Meaning of the Norms¶

11.2.2 The Eckart–Young Theorem: Best Low-Rank Approximation¶

One of the most important applications of the SVD is to provide the best low-rank approximation of a matrix. This result was proved by Carl Eckart and Gale Young in 1936 and is now called the Eckart–Young theorem (or the Eckart–Young–Mirsky theorem, which includes Leon Mirsky’s 1960 extension to unitarily invariant norms).

Given a matrix $\mathbf{A} \in \mathbb{C}^{m \times n}$ of rank $r$, we wish to find a matrix $\mathbf{B}$ of rank at most $k$ ($k < r$) such that $\mathbf{B}$ is “closest” to $\mathbf{A}$.

Mathematical formulation:

The crux of the problem is:

1. What does “closest” mean? That is, which norm should serve as the measure of distance?

2. What is the optimal solution? Can it be constructed explicitly?

The Eckart–Young theorem answers both questions perfectly.

The Frobenius-Norm Version of the Eckart–Young Theorem¶

A complete proof of the general case requires the von Neumann trace inequality or the Mirsky inequality and is beyond the scope of this book; still, it is worth marking clearly exactly where the gap in the proof lies and which tools this section has already prepared.

The Spectral-Norm Version of the Eckart–Young Theorem¶

Comparing the two norm versions:

An important observation. The Frobenius norm takes into account the contributions of all discarded singular values, while the spectral norm cares only about the largest discarded singular value. In practice:

• if the singular values decay rapidly (for example $\sigma_i \sim e^{-i}$), both norms give a very good approximation

• if the singular values decay slowly (for example $\sigma_i \sim i^{-1}$), more terms may need to be kept

The Practical Meaning of Low-Rank Structure¶

The Eckart–Young theorem reveals the central place of the SVD in data science:

1. Image compression

An $m \times n$ grayscale image can be regarded as a matrix $\mathbf{A} \in \mathbb{R}^{m \times n}$. Storing it in full requires $mn$ numbers. But if we keep the first $k$ singular values, we only need to store

(namely $k$ vectors $\mathbf{u}_i$ ($km$ numbers in all), $k$ vectors $\mathbf{v}_i$ ($kn$ numbers in all), and $k$ singular values)

Compression ratio:

If $k \ll \min(m,n)$, the compression is substantial. For example, with $m=n=1000$ and $k=50$:

That is, only 10% of the space of the original data is needed!

2. Denoising

Real measured data $\mathbf{A}$ usually contain noise:

If the signal has low-rank structure, the signal-to-noise ratio is high enough, and the main signal modes are separable from the noise, then the large singular values may be dominated by the signal; small singular values may also contain important signal.

The truncated SVD $\mathbf{A}_k$ is singular-mode filtering, which in general is not the same as low-pass filtering. Under the low-rank signal assumption above:

• keep the large singular values (possibly dominated by the signal)

• discard the small singular values (possibly dominated by noise, but weak signals may also be lost)

3. Latent Semantic Analysis (LSA)

In text mining, the document-term matrix $\mathbf{A}$ is usually high-dimensional and sparse ($m$ documents, $n$ terms; $mn$ is large but most entries are 0).

The truncated SVD maps the high-dimensional term space to a low-dimensional “concept space”:

• Left singular vectors $\mathbf{u}_i$: the coordinates of documents in concept space

• Right singular vectors $\mathbf{v}_i$: the coordinates of terms in concept space

• Singular values $\sigma_i$: the “importance” of the $i$-th concept

Keeping the first $k$ concepts (typically $k \sim 100-300$) makes it possible to:

• find semantically similar documents (even if they share no terms)

• reduce dimensionality and improve search efficiency

• remove interference from synonyms and polysemous words

For a fixed $m\times n$ image, the Frobenius energy is $\|\mathbf{A}\|_F^2$ and the mean squared error is $\|\mathbf{A}-\mathbf{A}_k\|_F^2/(mn)$. Minimizing the error norm, its square, or the mean squared error gives the same optimal matrix, although the values differ.

11.3 Application 1: Principal Component Analysis (PCA)¶

The Eckart–Young theorem of §11.2 (Theorem 8) is a purely matrix-optimization result: among all matrices of rank at most $k$, the truncated SVD $\tilde{\mathbf{X}}_k$ is the one nearest to the original matrix. This conclusion by itself carries no statistical meaning. But as soon as the matrix holds data—each row a sample, each column a feature—the “best low-rank approximation” immediately takes on a statistical identity: it is precisely principal component analysis (PCA).

So what exactly does PCA do? The best-known answer is “find the directions of greatest variation in the data,” which is the viewpoint of Hotelling (1933). This section, however, adopts a more fundamental viewpoint that fits §11.2 more closely—that of Pearson (1901): “the low-dimensional plane that best fits the data.” In other words, the essence of PCA is not “finding the directions of greatest variation” but “finding the best low-dimensional projection”; the directions of greatest variation are merely a by-product of this best projection. This order matters, because it determines how we understand PCA; we return in §11.3.2 to discuss “the directions of greatest variation” as an observation.

The measure used throughout this section is variance. Experiment 5 (Blocks C and D) has prepared the tool: the total variance $V_{\text{total}} = \operatorname{tr}(\mathbf{S}) = \frac{1}{m-1}\|\tilde{\mathbf{X}}\|_F^2$. With it, the argument of §11.3.1 proceeds in two steps. First, the Frobenius Pythagorean decomposition guarantees that any orthogonal projection conserves the total variance—dimensionality reduction neither creates nor erases variance out of nothing; it only divides it between the approximation and the error. Second, the Eckart–Young theorem picks out the best one among all projections of the same dimension, namely the one given by the right singular vectors of the SVD. Conservation is the universal background; optimization is the real protagonist of PCA.

This section proceeds in two steps. §11.3.1 establishes the core logic (conservation + best projection); §11.3.2 first points out that the best directions of the SVD are exactly “the directions of greatest variation in the data” (an incidental observation, together with a variational characterization by variance maximization), and then, from this nature of “looking only at linear variance,” marks the limits of PCA—it fails on curved manifolds. The rank-by-rank transfer of variance, and the application of PCA to real data to see it reveal cluster structure, are left to the PCA experiment notebook of this chapter for a full demonstration.

11.3.1 From Low-Rank Approximation to Variance Decomposition¶

The Eckart–Young theorem of §11.2.2 (Theorem 8) is a matrix-optimization result. This section translates it into the language of statistics and clarifies the true logical order of PCA: first comes the invariant “every projection conserves the total variance,” and then Eckart–Young picks the best projection from among them.

We continue with the centered data matrix $\tilde{\mathbf{X}} \in \mathbb{R}^{m\times n}$ of Experiment 5 (samples as rows), the sample covariance matrix $\mathbf{S} = \frac{1}{m-1}\tilde{\mathbf{X}}^{\top}\tilde{\mathbf{X}}$, and the total variance

(For readers who have not done Experiment 5, one line suffices: $\operatorname{tr}(\mathbf{S}) = \sum_j \frac{1}{m-1}\sum_i \tilde{x}_{ij}^2 = \frac{1}{m-1}\|\tilde{\mathbf{X}}\|_F^2$.) Take the SVD of $\tilde{\mathbf{X}}$ and write it as an outer-product expansion (§11.1):

Step 1: Every Orthogonal Projection Conserves the Total Variance¶

Dimensionality reduction means projecting the data onto some low-dimensional subspace. Take any orthogonal projection matrix $\boldsymbol{\Pi}\in\mathbb{R}^{n\times n}$ (not yet specified); the projected data are $\tilde{\mathbf{X}}\boldsymbol{\Pi}$ and the residual is $\tilde{\mathbf{X}}(\mathbf{I}-\boldsymbol{\Pi})$. How do the variances of these two divide the total variance?

The whole point of this equality lies in the word “conservation”: no matter which subspace we project onto, the total variance splits exactly into two pieces, approximation and error, and the total does not change one bit. So dimensionality reduction is never a question of “how much variance can be preserved”—that holds automatically; the real question is “within the fixed total, how to keep as much variance as possible in the approximation and lose as little as possible to the error.” This leads to Step 2.

Step 2: Eckart–Young Selects the Best Projection¶

The choice made by PCA is thus justified: we choose the right singular vectors of the SVD not because these directions are special in themselves, but because projecting onto the subspace they span preserves the most variance among all projections of the same dimension. As for the geometric appearance of these best projection directions—they happen to be “the directions of greatest variation in the data”—we leave that to §11.3.2 as an observation.

11.3.2 Directions of Greatest Variation, and the Limitations of PCA¶

Look back at the best projection $\boldsymbol{\Pi}_k$ chosen in §11.3.1. Taking $k=1$, its direction is $\mathbf{v}_0$, which preserves the most variance among all unit directions; in turn, $\mathbf{v}_0, \mathbf{v}_1, \ldots$ are “the directions of variation in the data from largest to smallest”—the 0-th, the 1st, ... principal components. This is also the best-known description of PCA: finding the directions of greatest variation in the data. But we must stress its place in the logic of this chapter: it is an observation that falls out of the “best projection” incidentally, not the goal we are pursuing. What we want throughout is the best low-rank approximation; the direction of greatest variation is just how it looks when $k=1$.

The Variance-Maximization Viewpoint¶

Let us write the observation above as a precise step-by-step definition. The 0-th principal component $\mathbf{v}_0$ is the unit vector that maximizes the variance of the projected data:

(Since $\tilde{\mathbf{X}}$ is already centered, the projection scores $\tilde{\mathbf{X}}\mathbf{v}$ have mean zero, so their variance is $\frac{1}{m-1}\|\tilde{\mathbf{X}}\mathbf{v}\|^2$.) This is exactly the variational characterization of the SVD in §11.1.3 (Theorem 4)! The 1st principal component $\mathbf{v}_1$ is the vector that maximizes the variance among directions orthogonal to $\mathbf{v}_0$:

Continuing in this way, the $i$-th principal component $\mathbf{v}_i$ is the direction of maximal variance orthogonal to the first $i$ principal components $\mathbf{v}_0, \ldots, \mathbf{v}_{i-1}$. This characterization—“successively maximize the variance in the orthogonal complement”—arrives at the same place as the “best projection” of §11.3.1: the former picks greedily one direction at a time, the latter chooses the whole subspace at once, and both select the right singular vectors of the SVD.

The two are not equally “strong,” however. Eckart–Young is a global result: it weighs all $k$-dimensional subspaces simultaneously and determines the entire optimal subspace in one stroke. The greedy algorithm is local (greedy): it picks one direction after another, at each step performing only a one-dimensional maximization under the constraint “orthogonal to the directions already chosen,” and never goes back. In general, local greed is a weaker, operational result—stepwise optimality does not guarantee global optimality. In this problem, however, the two are exactly equivalent: each step maximizes the variance in the orthogonal complement of the first $i-1$ principal components and selects precisely the next singular direction, and the accumulated steps span exactly the globally optimal subspace of Eckart–Young. This coincidence of “greedy equals global” stems from the nested structure of the eigenspaces of a symmetric positive semidefinite matrix and is not a property enjoyed by optimization in general—and precisely because of it, PCA can reach the global optimum by the plainest step-by-step procedure.

The Limitations of PCA¶

It is precisely this nature of “looking only at linear variance” that marks the limits of PCA. Once the structure of the data is not linear, the conservation of §11.3.1 still holds (the total variance is always conserved), but “keeping the directions of greatest variance” need not mean “keeping the most useful structure.” Here are several limitations to remember:

1. Linearity assumption: PCA can capture only linear relationships. For nonlinear structures (spirals, rings), kernel PCA or manifold learning (t-SNE, UMAP) is needed.

2. Global method: PCA looks for the globally optimal linear subspace and is insensitive to local structure.

3. Sensitivity to outliers: variance maximization can be dominated by outliers. Robust alternatives include Robust PCA and L1-PCA.

4. Interpretability: principal components are linear combinations of the original features and often have a less clear meaning than the original features.

5. Unsupervised: PCA does not use labels, so in classification/regression tasks it may discard directions useful for prediction. In that case, supervised dimensionality reduction (such as Linear Discriminant Analysis, LDA) should be considered.

The first of these is the most fatal and deserves to be seen with your own eyes. Consider the classic Swiss roll: the data lie on a rolled-up two-dimensional sheet in three-dimensional space. Below, the height direction $y$ is deliberately stretched so that it becomes the direction of greatest variance (the colors mark four segments cut along $t$). The left panel shows the original 3D structure—the spiral is clearly visible, and the height axis is just a noise direction superimposed on it. The middle panel takes PCA’s default first two principal components PC0–PC1, but PC0 is exactly the height axis that has nothing to do with the structure, so the four segments are stirred into overlapping bands and no spiral can be seen. The right panel takes PC1–PC2 (the spiral plane), which PCA regards as secondary, and there the four segments separate clearly along the spiral. The lesson is sharp: PCA always picks the directions of greatest variance, but the greatest variance is not necessarily the most meaningful direction—when an irrelevant high-variance axis gets mixed in, PCA’s default projection actually hides the structure. Remedy: at the linear level, one can look at the secondary principal components instead or perform feature selection first; to truly “unroll” the spiral (the spiral in the right panel is still curled), one needs manifold learning (Isomap, LLE) or kernel PCA, which measure distance along the manifold rather than along Euclidean straight lines.

11.4 Application 2: The Moore–Penrose Pseudoinverse¶

In the discussion of systems of linear equations in Chapter 6, we saw that for a square matrix $\mathbf{A} \in \mathbb{R}^{n \times n}$, if $\mathbf{A}$ is invertible, then the equation $\mathbf{A}\mathbf{x} = \mathbf{b}$ has the unique solution $\mathbf{x} = \mathbf{A}^{-1}\mathbf{b}$. Real-world problems, however, often involve non-square or rank-deficient matrices:

• Overdetermined systems: $m > n$, more equations than unknowns, usually no exact solution

• Underdetermined systems: $m < n$, fewer equations than unknowns, usually infinitely many solutions

• Rank-deficient square matrices: $\text{rank}(\mathbf{A}) < n$, not invertible even though square

Question. Can we define a “generalized inverse” that exists for every matrix $\mathbf{A} \in \mathbb{R}^{m \times n}$ and reduces to the ordinary inverse when $\mathbf{A}$ is invertible?

This section takes a path from the concrete to the abstract: first we let the simplest case, diagonal matrices, reveal the rudiments of the answer (§11.4.1); then we use the SVD, a “universal rotation,” to carry the same trick over unchanged to an arbitrary matrix, obtaining the general least-squares solution (§11.4.2); only at the end do we go back and verify that the matrix so constructed satisfies the four Penrose conditions, and cite its uniqueness theorem—this is the formal theory of the Moore–Penrose pseudoinverse (§11.4.3).

11.4.1 Starting from Diagonal Matrices: The Rudiments of the Pseudoinverse¶

In a system of equations with an arbitrary matrix, the variables are entangled with one another, and it is not easy to see at a glance how “inversion” should be generalized. But in a system with a diagonal matrix, the components do not interfere with one another, and the answer can be read off term by term.

The trivial case. If $\boldsymbol{\Sigma} = \text{diag}(\sigma_0,\ldots,\sigma_{n-1})$ is an $n\times n$ square matrix and all $\sigma_i \neq 0$, then $\boldsymbol{\Sigma}\mathbf{x}=\mathbf{b}$ can be solved term by term: $x_i = b_i/\sigma_i$, that is, $\boldsymbol{\Sigma}^{-1}=\text{diag}(1/\sigma_0,\ldots,1/\sigma_{n-1})$—the familiar inverse matrix, with nothing new.

The really interesting case is the rank-deficient one:

This example states the rule clearly: the pseudoinverse $\boldsymbol{\Sigma}^{+}$ of a diagonal matrix (including a rectangular diagonal matrix) is the matrix obtained by taking reciprocals of the nonzero diagonal entries, transposing the shape, and filling all other positions with zeros. A column corresponding to a zero diagonal entry (a free variable) naturally becomes a zero row in $\boldsymbol{\Sigma}^{+}$—because setting that variable to zero is exactly the choice that minimizes the norm; a row corresponding to a zero diagonal entry (an equation that cannot be satisfied) is ignored, contributes nothing to the solution, but leaves an irreducible residual.

This intuition holds for every diagonal matrix, square or rectangular, full rank or rank-deficient. The problem is that a general matrix $\mathbf{A}$ is not diagonal, and its components are entangled with one another. The next subsection uses the SVD to “untangle” them, reducing the general case exactly to the diagonal case already solved here.

11.4.2 Generalizing to Least Squares with the SVD¶

The key idea is to use the SVD to “rotate” the problem for an arbitrary matrix into the diagonal problem already solved in §11.4.1. Let $\mathbf{A}=\mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{\top}\in\mathbb{R}^{m\times n}$ and consider

Change of coordinates. $\mathbf{U}$ and $\mathbf{V}$ are both orthogonal matrices, and multiplying by an orthogonal matrix does not change the length of a vector (for any $\mathbf{w}$, $\|\mathbf{U}\mathbf{w}\|=\|\mathbf{w}\|$). Therefore

Let $\mathbf{y}:=\mathbf{V}^{\top}\mathbf{x}$ (again, since $\mathbf{V}$ is orthogonal, $\|\mathbf{y}\|=\|\mathbf{x}\|$) and $\mathbf{c}:=\mathbf{U}^{\top}\mathbf{b}$. The problem becomes

—exactly the diagonal problem already solved in §11.4.1! Moreover, $\|\mathbf{x}\|=\|\mathbf{y}\|$ guarantees that the minimum-norm solution found in $\mathbf{y}$-coordinates is still the minimum-norm solution after changing back to $\mathbf{x}$-coordinates. By §11.4.1, this solution is $\hat{\mathbf{y}}=\boldsymbol{\Sigma}^{+}\mathbf{c}$ (reciprocals of the nonzero singular values, zeros elsewhere). Changing back to the original coordinates:

This is the construction formula for the pseudoinverse—it is not assumed out of thin air, but is the inevitable result of carrying the solution method for diagonal matrices over unchanged to general matrices along the “universal rotation” of the SVD. Let us state this argument as a formal theorem:

11.4.3 The Moore–Penrose Pseudoinverse in Theory¶

We have now constructed the same matrix $\mathbf{A}^{+}=\mathbf{V}\boldsymbol{\Sigma}^{+}\mathbf{U}^{\top}$ concretely twice—first seeing its logic in the rudimentary diagonal case of §11.4.1, then carrying it over to general matrices with the SVD in §11.4.2 and verifying that it does give the least-squares, minimum-norm solution. This subsection closes the theory: we verify that the constructed matrix satisfies the four Penrose conditions; for uniqueness we cite Penrose’s theorem without proving it here—this is precisely the characterization Penrose gave in 1955, and it is the origin of the name “Moore–Penrose pseudoinverse.”

Definition of the Moore–Penrose Pseudoinverse¶

Constructing the Pseudoinverse from the SVD¶

Through the SVD, the pseudoinverse has an extremely concise expression (this is precisely the construction formula we derived in §11.4.2; here we take a different angle and formally prove that it satisfies the Penrose axioms):

Outer-product expansion:

This formula shows clearly that the pseudoinverse takes the reciprocal of the singular value in each term of the outer-product expansion of $\mathbf{A}$ and swaps the left and right singular vectors.

The Geometric Action of the Pseudoinverse¶

The geometric meaning of $\mathbf{A}^{+}$ can be understood through the four fundamental subspaces (recall §11.1.2):

How $\mathbf{A}^{+}$ acts:

Given a vector $\mathbf{b} \in \mathbb{R}^m$, compute $\mathbf{x} = \mathbf{A}^{+}\mathbf{b}$:

1. Filter: project $\mathbf{b}$ onto the image $\text{Col}(\mathbf{A})$, ignoring its component in the left null space $\ker(\mathbf{A}^{\top})$

   $$\mathbf{b}_{\text{proj}} = \mathbf{A}\mathbf{A}^{+}\mathbf{b} = \sum_{i=0}^{r-1}⟨ \mathbf{u}_i|\mathbf{b}⟩\mathbf{u}_i$$

   (180)

2. Inverse mapping: map the components in the image back to the row space $\text{Row}(\mathbf{A})$ through the reciprocals of the singular values

   $$\mathbf{x} = \mathbf{A}^{+}\mathbf{b} = \sum_{i=0}^{r-1}\frac{⟨ \mathbf{u}_i|\mathbf{b}⟩}{\sigma_i}\mathbf{v}_i$$

   (181)

3. Output: we obtain the unique vector in the row space

Diagram:

In the diagram, the null space $\ker(\mathbf{A})$ is not connected to any arrow: the range of $\mathbf{A}^{+}$ is the row space $\text{Row}(\mathbf{A})$, and the solution $\mathbf{A}^{+}\mathbf{b}$ never has a null-space component (it is orthogonal to the null space).

Key properties:

1. $\mathbf{A}\mathbf{A}^{+}$ is the orthogonal projection onto $\text{Col}(\mathbf{A})$

2. $\mathbf{A}^{+}\mathbf{A}$ is the orthogonal projection onto $\text{Row}(\mathbf{A})$

3. $\mathbf{A}^{+}\mathbf{b}$ always lies in the row space (that is, $\mathbf{A}^{+}\mathbf{b} \perp \ker(\mathbf{A})$)

11.4.4 Numerical Instability and Ill-Conditioned Problems¶

Although the pseudoinverse solves the least-squares problem perfectly in theory, in numerical computation it can be extremely unstable!

The Problem of Small Singular Values in Ill-Conditioned Matrices¶

Recall the formula for the pseudoinverse:

The problem. When some singular values $\sigma_i$ are very small, $1/\sigma_i$ becomes very large, which leads to:

1. Sensitivity to noise: tiny perturbations in $\mathbf{b}$ are amplified by a factor of $1/\sigma_i$

2. Numerical overflow: an extremely small $\sigma_i$ may make $1/\sigma_i$ exceed the range of floating-point numbers

3. Unstable solutions: different numerical algorithms may give completely different solutions

Truncated SVD¶

The remedy. Ignore singular values that are too small and keep only the first $k$:

Given a tolerance threshold $\tau\geq0$, let $k=\#\{i:0\leq i<r,\ \sigma_i>\tau\sigma_0\}$. If $0<k<r$, this is equivalent to $\sigma_{k-1}>\tau\sigma_0\geq\sigma_k$; when $k=0$ we take the zero matrix, and when $k=r$ all positive singular values are kept. For the zero matrix we set $k=0$ separately.

Effect:

• Advantage: lowers the noise sensitivity of the truncated directions, but introduces truncation bias

• Disadvantage: it is not a strict least-squares solution, and it introduces bias

This leads to a more systematic regularization method—ridge regression.

11.4.5 Ridge Regression: Tikhonov Regularization¶

From Ill-Conditioned Problems to Regularization¶

The least-squares problem is

When $\mathbf{A}$ is ill-conditioned (has small singular values), the norm of the solution $\hat{\mathbf{x}} = \mathbf{A}^{+}\mathbf{b}$ may be very large and extremely sensitive to noise.

Tikhonov regularization (also called ridge regression) stabilizes the solution by adding a norm penalty term:

where $\lambda > 0$ is the regularization parameter.

Key Properties of Ridge Regression¶

1. Numerical stability

Compare the coefficients:

• Pseudoinverse: $\frac{1}{\sigma_i}$ ($\to \infty$ as $\sigma_i \to 0$)

• Ridge regression: $\frac{\sigma_i}{\sigma_i^2 + \lambda}$ ($\to 0$ as $\sigma_i \to 0$)

The ridge regression coefficients are bounded:

The amplification along small-singular-value directions is suppressed, trading bias for lower noise sensitivity.

2. The bias-variance trade-off

• $\lambda\to0^+$: the solution approaches the minimum-norm least-squares solution $\mathbf{A}^+\mathbf{b}$; in the rank-deficient case one cannot simply substitute $\lambda=0$ into the inverse-matrix formula

• $\lambda \to \infty$: $\hat{\mathbf{x}}_{\lambda} \to \mathbf{0}$, with small variance and large bias

• Moderate $\lambda$: strikes a balance between bias and variance

3. The filter viewpoint

Define the filter factor

Then

Properties of the filter factor:

• when $\sigma_i \gg \sqrt{\lambda}$, $f_i \approx 1$ (almost no attenuation)

• when $\sigma_i \ll \sqrt{\lambda}$, $f_i \approx \sigma_i^2/\lambda$ (strong attenuation)

Ridge regression is singular-mode filtering that suppresses the amplification along small-singular-value directions; the size of a singular value is in general not a temporal or spatial frequency, so it cannot be called low-pass filtering across the board.

Choosing the Regularization Parameter $\lambda$¶

In practice, the choice of $\lambda$ is crucial.

Common methods:

1. Cross-validation: choose the $\lambda$ that minimizes the prediction error on a validation set

2. The L-curve method: plot $\log\|\mathbf{A}\hat{\mathbf{x}}_{\lambda} - \mathbf{b}\|$ against $\log\|\hat{\mathbf{x}}_{\lambda}\|$ and choose the $\lambda$ at the “corner”

3. Generalized cross-validation (GCV): a statistical method that balances bias and variance automatically

4. Mode-transition scale: when $\lambda=\sigma_i^2$, $f_i=1/2$. To approximately keep the first $k$ modes ($0<k<r$), one may consider $\sigma_k^2\lesssim\lambda\lesssim\sigma_{k-1}^2$, but the parameter must be chosen according to the noise and the data; this is not an optimality theorem

The code below actually scans $\lambda$ from 10^-18 to 10^-2 for the ill-conditioned system of Example 14 and draws the L-curve: the horizontal axis is the norm of the solution $\|\hat{\mathbf{x}}_{\lambda}\|$ and the vertical axis is the residual $\|\mathbf{A}\hat{\mathbf{x}}_{\lambda}-\mathbf{b}\|$ (log-log axes). The corner of the L-curve is a heuristic candidate for the parameter; the marker in the figure only demonstrates a preset parameter and does not compute the corner or the optimum automatically.

11.5 ◆Application 3: The Schmidt Decomposition and Quantum Entanglement¶

Quantum entanglement is one of the most mysterious and most important phenomena of quantum mechanics—Einstein called it “spooky action at a distance”—and it is also the theoretical foundation of technologies such as quantum computing, quantum communication, and quantum sensing. The SVD provides a precise mathematical tool for quantifying entanglement—the Schmidt decomposition: given the state vector of a composite quantum system, the Schmidt decomposition reveals the structure and strength of the entanglement between the subsystems.

The role of mathematics. The SVD/Schmidt decomposition provides:

1. A quantitative tool: how can the “strength” of entanglement be measured?

2. Structural understanding: what is the internal structure of an entangled state?

3. A computational method: how can high-dimensional entangled states be handled efficiently?

11.5.1 Matricizing the State Vector and the Schmidt Decomposition¶

The Tensor Product Structure of Composite Quantum Systems¶

Two views of the tensor product, revisited In Chapter 4 we introduced the abstract definition of the tensor product space. Before entering this section, let us become familiar with it again through a concrete example and emphasize a key idea that recurs throughout quantum mechanics: the state vector of a composite system can be viewed either as a “long vector” or as a “matrix.” These two views do not contradict each other; rather, they are where the SVD comes into play.

Consider two quantum systems:

• System A: Hilbert space $\mathcal{H}_A$, dimension $d_A$, orthonormal basis $\{|i⟩_A\}_{i=0}^{d_A-1}$

• System B: Hilbert space $\mathcal{H}_B$, dimension $d_B$, orthonormal basis $\{|j⟩_B\}_{j=0}^{d_B-1}$

The state space of the composite system A⊗B is the tensor product space

of dimension $d_{AB} = d_A \times d_B$, with orthonormal basis

abbreviated as $\{|ij⟩\}$.

A general pure state can be expanded as

where the $c_{ij} \in \mathbb{C}$ are complex amplitudes satisfying the normalization condition

Key observation. The coefficients $c_{ij}$ can be arranged into a $d_A \times d_B$ matrix!

The Coefficient Matrix and the SVD¶

Define the coefficient matrix $\mathbf{C} \in \mathbb{C}^{d_A \times d_B}$ by

Then the pure state can be written as

and the normalization condition becomes

that is, $\mathbf{C}$ has Frobenius norm 1.

The core idea. Take the SVD of $\mathbf{C}$:

where:

• $\mathbf{U} \in \mathbb{C}^{d_A \times d_A}$ is a unitary matrix

• $\mathbf{V} \in \mathbb{C}^{d_B \times d_B}$ is a unitary matrix

• $\boldsymbol{\Sigma} \in \mathbb{R}^{d_A \times d_B}$ is a diagonal matrix whose diagonal entries are the singular values $\sigma_0 \geq \sigma_1 \geq \cdots \geq 0$

Let the rank be $r = \text{rank}(\mathbf{C})$, so that the first $r$ singular values are nonzero.

The meaning of the Schmidt decomposition:

1. Conciseness: whereas $d_A \times d_B$ coefficients were originally needed, the Schmidt decomposition needs only $r$ coefficients ($r \leq \min(d_A, d_B)$)

2. Symmetry: the descriptions of the two subsystems A and B are completely symmetric, each with its own orthonormal basis

3. Uniqueness: the Schmidt coefficients are unique and directly describe the structure of the entanglement

Proof. This is immediate. Schmidt rank 1 means there is only one nonzero Schmidt coefficient, in which case $|\psi⟩ = \lambda_0|\alpha_0⟩_A|\beta_0⟩_B$ is separable. Conversely, if the Schmidt rank is $\geq 2$, there are at least two terms, and the state cannot be written as a single tensor product. $\square$

Schmidt Rank and the von Neumann Entanglement Entropy¶

The Schmidt rank $r$ is the coarsest measure of entanglement:

But the Schmidt rank is not fine enough. For example:

Both have Schmidt rank 2, but their degrees of entanglement are clearly different!

A finer measure is the von Neumann entanglement entropy:

Properties:

1. Nonnegativity: $S_{\text{ent}} \geq 0$

2. Zero for separable states: $r = 1 \implies S_{\text{ent}} = 0$

3. Maximally entangled states have the largest entropy: for $d_A = d_B = d$, if $\lambda_k = 1/\sqrt{d}$ for all $k$, then

   $$S_{\text{ent}} = -d \cdot \frac{1}{d}\log_2(1/d) = \log_2 d$$

   (262)

4. Symmetry between subsystems: $S_{\text{ent}}$ describes the entanglement between A and B and does not depend on which one is the “main” system

11.5.2 Reduced Density Matrices and Purification¶

The Partial Trace and Quantum Measurement¶

The partial trace and the projective measurement of §10.4 are related but different operations, and they deserve to be distinguished carefully.

Measuring B and the conditional pure state of A

For a composite pure state $|\psi⟩_{AB}=\sum_{i,j}c_{ij}|i⟩_A|j⟩_B$, perform a rank-1 projective measurement of B in any orthonormal basis; if outcome $j_0$ has positive probability, the conditional state of A is the definite pure state obtained by normalizing $\sum_i c_{i,j_0}|i⟩_A$. The A states corresponding to different outcomes need not be orthogonal.

When measuring in the Schmidt basis $\{|\beta_k⟩_B\}$, $|\psi⟩_{AB}=\sum_k\lambda_k|\alpha_k⟩_A|\beta_k⟩_B$, and each outcome $k$ of positive probability corresponds to $|\alpha_k⟩_A$, giving orthogonal conditional states in one-to-one correspondence. A definite conditional pure state does not mean that measuring any observable on it shows no statistical fluctuation.

The partial trace corresponds to the average with “outcomes discarded”

Still measuring B in the basis $\{|j⟩_B\}$, the probability of outcome $j$ is

Only when $P(j)>0$ is the corresponding conditional pure state of A defined, namely $|\phi_j⟩_A = \dfrac{1}{\sqrt{P(j)}}\sum_i c_{ij}|i⟩_A$, with density matrix $\rho_A^{(j)} = |\phi_j⟩_A⟨ \phi_j|$.

If the outcome is not read after the measurement (or B simply cannot be accessed), the effective state of A is the probability-weighted average of all the conditional states:

This weighted average does not depend on the choice of measurement basis for B and equals the partial trace $\operatorname{Tr}_B(\rho_{AB})$. We verify this explicitly below in the language of matrix indices.