Sources

The concepts and notation used in this post are inspired from the following sources, which I encourage the reader to check out:

• Data-Driven Science and Engineering, Steven L. Brunton and J. Nathan Kutz, Cambridge University Press (2021)
• Measure, Integration & Real Analysis (MIRA), Sheldon Axler, open publication (2024)
• Regularization of Inverse Problems, Christian Clason, arXiv (2021)

Basic definitions

Linear maps

A linear map is a function that preserves the zero vector, vector addition and scalar multiplication between vector spaces. Geometrically, straight lines through the origin stay straight and the origin stays put.

Definition: Linear map

Let \(V\) and \(W\) be vector spaces over the scalar field \(\mathbf{F}\) and \(T:V\to W\) be a function. Then \(T\) is called a linear function or a linear map if \[T(ax+by)=aTx+bTy\] for all \(a,b\in\mathbf{F}\) and \(x,y\in V\).

Linear maps are themselves elements of a vector space. This space is often denoted by \(\mathcal{L}(V,W)\). To make this space a normed space, we need some way to measure the “size” of the linear maps. For this we can use the operator norm.

Definition: Operator norm, bounded linear map

Let \(V\) and \(W\) be Hilbert spaces and \(T:V\to W\) be a linear map. Then \(T\) is called bounded if the operator norm, defined by \[\|T\|:=\sup_{x\in V\setminus\lbrace0\rbrace}\frac{\|Tx\|_W}{\|x\|_V}\] is finite.

Note that the operator norm is equivalently expressed as \[\|T\|=\sup_{\|x\|_V=1}\|Tx\|_W.\] This follows by writing \(x=\|x\|_V\frac{x}{\|x\|_V}\) to get \[\sup_{x\in V\setminus\lbrace0\rbrace}\frac{\|Tx\|_W}{\|x\|_V}=\sup_{x\in V\setminus\lbrace0\rbrace}\frac{\|T(\|x\|_V\frac{x}{\|x\|_V})\|_W}{\|x\|_V}=\sup_{x\in V\setminus\lbrace0\rbrace}\|T\frac{x}{\|x\|_V}\|_W=\sup_{\|v\|_V=1}\|Tv\|_W,\] where the last part of transforming the supremum set is just another way to state that the vector \(x/\|x\|_V\) always has norm/length \(1\).

By definition, a linear map \(T\) is continuous if for every \(\epsilon>0\) there exists some \(\delta>0\) so that \[\|Tx-Ty\|_W\le\epsilon,\quad\forall x,y\in V\text{ with }\|x-y\|_V\le\delta.\] But since \(T\) is linear, we have by definition \[\|Tx-Ty\|_W=\|T(x-y)\|_W.\] If \(T\) is bounded, there exists some constant \(C>0\) such that \[\frac{\|Tx\|_W}{\|x\|_V}\le C\] for all \(x\in V\setminus\lbrace 0\rbrace\), which implies \[\|Tx\|_W\le C\|x\|_V.\] Clearly now \(\|Tx\|_W\le\epsilon\) by choosing any \(x\in V\) such that \(\|x\|_V\le\epsilon/C\), so that \(T\) is continuous. On the other hand, if \(T\) is continuous, it is clearly bounded. Therefore, a linear map is bounded if and only if it is a continuous function.

Another way of stating this is that a continuous linear function maps “small enough” closed balls to small closed balls: for every \(\epsilon>0\) there exists a \(\delta>0\) so that \[T\overline{B_V(0,\delta)}\subseteq\overline{B_W(0,\epsilon)}.\] This tells us that if a linear map is not bounded, it is also not continuous.

Adjoint operators

Definition: Adjoint operator

Let \(V\) and \(W\) be Hilbert spaces equipped with the inner products \(\langle\cdot,\cdot\rangle_V\) and \(\langle\cdot,\cdot\rangle_W\). Let \(T:V\to W\) be a linear map. The adjoint of \(T\) is the map \(T^*:W\to V\) such that \[\langle Tf,g\rangle_W=\langle f,T^*g\rangle_V\] for all \(f\in V\) and \(g\in W\).

Definition: Self-adjoint operator

If \(T=T^*\), we call \(T\) self-adjoint.

Intuitively, in certain cases the adjoint operator can be seen as a “reverse-direction” operator of the forward operator. For example, consider the transport PDE operator \[T=\partial_t+c\partial_x\] defined on a suitable domain. This PDE is used to describe e.g. fluid flow. Informally, the adjoint operator’s action (ignoring the details regarding the adjoint operator’s domain for the sake of brevity) is \[T^* = -\partial_t – c\partial_x.\] Physically, if \(T\) moves elements forward in time and space, \(T^*\) moves them backward, which justifies the intuition of the adjoint being a reverse operator.

Compact operators

Intuitively, a compact set is a set without much “room” so that any sequence in such a set must have values repeatedly getting arbitrarily close to each other. With this intuition in mind, it makes sense that in the finite-dimensional case, a set is compact if and only if it is closed and bounded. A compact operator maps bounded sets to relatively compact sets (i.e. the closure of the mapped set is compact).

Definition: Compact operator

A linear map \(T:X\to Y\) defined on Hilbert spaces \(X\) and \(Y\) is called compact if for every bounded sequence \(f_1,f_2,\dots\) in \(X\) the sequence \(Tf_1,Tf_2,\dots\) in \(Y\) has a convergent subsequence.

We state some useful results without proof:

• If \(T\) is compact, then \(T^*\) is compact.
• The composition of compact operators is compact, so that, in particular, \(T^*T\) is compact as well.

Spectral theorem

The spectral theorem is an important result when deriving the singular value decomposition:

Theorem: Spectral theorem

Let \(T\) be a self-adjoint compact operator on a Hilbert space \(V\). Then \(V\) has an orthonormal basis consisting of the eigenvectors of \(T\), and there is a countable set \(\Omega\), an orthonormal set \(\lbrace e_k\rbrace_{k\in\Omega}\) in \(V\), and a set \(\lbrace \alpha_k\rbrace_{k\in\Omega}\) in \(\mathbf{R}\setminus\lbrace0\rbrace\) such that \[Tf=\sum_{k\in\Omega}\alpha_k\langle f,e_k\rangle e_k\] for all \(f\in V\).

The set \(\lbrace e_k\rbrace_{k\in\Omega}\) forms an orthonormal basis for the closure of the range of \(T\), \(\overline{R(T)}\).

Singular value decomposition

The singular value decomposition (SVD) now follows directly from these preliminaries. While the spectral theorem requires an operator to be self-adjoint and map a space back to itself (\(V \to V\)), SVD generalizes this to any compact operator acting between two different spaces \(X \to Y\). It allows us to decompose the operator into a weighted sum of orthonormal modes. Furthermore, because every linear operator on a finite-dimensional space has finite rank and is automatically compact, every matrix possesses an SVD.

Theorem: Singular value decomposition

Let \(V\) be a Hilbert space and \(T:X\to Y\) be a compact linear operator. Then there exists

• 1. a sequence \((s_n)\) such that \(s_1\ge s_2\ge\dots> 0\),
  2. an orthonormal basis \((u_n)\subseteq Y\) for the space \(\overline{R(T)}\), and
  3. an orthonormal basis \((v_n)\subseteq X\) for the space \(\overline{R(T^*)}\)

such that \[Tv_n=s_nu_n,\quad T^*u_n=s_nv_n\] for all \(n\in\mathbb{N}\) and \[Tx=\sum_{n\in\mathbb{N}}s_n\langle x,v_n\rangle_Xu_n\] for all \(x\in X\).

First, since \(T\) is compact, so is \(T^*\). Thus \(T^*T\) is a self-adjoint compact operator. Therefore, by the spectral theorem it has an eigen-decomposition \[(T^*T)x=\sum_{k\in\Omega}\alpha_k\langle x,v_k\rangle v_k.\] The eigenvalues of \(T^*T\) are real-valued and non-negative, as is shown in the following claims.

Claim 1: The \(\alpha_k\) are real-valued.

Proof: By definition we have \[(T^*T)v_k=\alpha_kv_k.\] Thus \[\langle (T^*T)v_k,v_k\rangle_X=\langle\alpha_kv_k,v_k\rangle_X=\alpha_k\|v_k\|_X^2.\] Using the self-adjoint property and the definition of the inner product (with the convention of complex conjugate holding for the second slot), we get \[\langle (T^*T)v_k,v_k\rangle_X=\langle v_k,(T^*T)v_k\rangle_X=\overline{\alpha_k}\|v_k\|_X^2.\] Therefore \(\alpha_k=\overline{\alpha_k}\) and \(\alpha_k\) is a real number. \(\blacksquare\)

Claim 2: The \(\alpha_k\) are non-negative.

Proof: \[0\le \|Tv_k\|_Y^2=\langle Tv_k,Tv_k\rangle_Y=\langle v_k,(T^*T)v_k\rangle_X=\alpha_k\|v_k\|_X^2.\] Since eigenvectors can not equal the zero vector, \(\|v_k\|_X^2>0\) and thus \(\alpha_k\ge 0\). \(\blacksquare\)

Based on claim 1 and claim 2, the value \(s_k=\sqrt{\alpha_k}\) is a well-defined real number. The positive \(s_k\) are called the singular values of \(T\) and we can sort them in decreasing order.

Next we look at the relationship between the ranges and null spaces of the maps \(T^*T\) and \(T\). First we show that the orthogonal complement of the null space of a compact operator is equivalent to the closure of the adjoint operator’s range.

Claim 3: \(N(K^*)^\perp=\overline{R(K)}\) and \(N(K)^\perp=\overline{R(K^*)}\) for any compact operator \(K:X\to Y\).

Proof: \[\begin{align*} u\in N(K^*)&\iff K^*u=0\\ &\iff \langle v,K^*u\rangle_X=0\quad\forall v\in X\\ &\iff\langle Kv,u\rangle_Y=0\quad \forall v\in X\\ &\iff u\in R(K)^\perp \end{align*}\] The claim follows from the fact that \((A^\perp)^\perp=\overline{A}\) for any set \(A\). The argument for \(N(K)^\perp\) is similar. \(\blacksquare\)

The null space of \(T^*T\) is equivalent to the null space of \(T\):

Claim 4: \(N(T^*T)=N(T)\).

Proof: Let \(x\in N(T^*T)\). Then \[T^*Tx=0\implies \langle (T^*T)x,x\rangle_X=0\implies\|Tx\|_Y^2=0\implies Tx=0.\] Thus \(N(T^*T)\subseteq N(T)\). The other direction holds trivially. \(\blacksquare\)

From claims 3 and 4 it follows that \[\overline{R(T^*T)}=N(T^*T)^\perp=N(T)^\perp=\overline{R(T^*)}\] Since the set \(\lbrace v_k\rbrace\) forms an orthonormal basis for \(\overline{R(T^*T)}\), this shows that it also forms an orthonormal basis for \(\overline{R(T^*)}\).

Now, let \[u_n=\frac{Tv_n}{\sqrt{\alpha_n}},\quad \alpha_n>0\] and consider the set \(\lbrace u_n\rbrace\subseteq Y\). Since \(\lbrace v_k\rbrace\subseteq N(T)^\perp\), \(Tv_n\ne 0\) for all \(n\) and therefore \(u_n\ne 0\) for each \(n\). The set \(\lbrace u_n\rbrace\) forms an orthonormal set: \[\begin{align} \langle u_n,u_m\rangle_Y&=\langle \frac{1}{\sqrt{\alpha_n}}Tv_n,\frac{1}{\sqrt{\alpha_m}}Tv_m\rangle_Y\\ &=\langle \frac{1}{\sqrt{\alpha_n}}v_n,\frac{1}{\sqrt{\alpha_m}}T^*Tv_m\rangle_X\\ &=\frac{1}{\sqrt{\alpha_n}\sqrt{\alpha_m}}\langle v_n,\alpha_mv_m\rangle_X\\ &=\delta_{nm}, \end{align}\] where \(\delta_{nm}\) is the Kronecker delta operator. Thus, \(\lbrace u_n\rbrace\subseteq R(T)\) is an orthonormal set. Next we want to show that it forms a basis for \(\overline{R(T)}\).

Since \(\lbrace v_n\rbrace\) is a basis of \(\overline{R(T^*)}=N(T)^\perp\), for any \(x\in X\) we can write (Hilbert space projection theorem) \[x=P_{N(T)}+\sum_{n\in\mathbb{N}}\langle x,v_n\rangle_Xv_n.\] Applying \(T\) to \(x\), we get \[Tx=\sum_{n\in\mathbb{N}}\langle x,v_n\rangle_XTv_n=\sum_{n\in\mathbb{N}}\langle x,v_n\rangle_X\sqrt{\alpha_n}u_n=\sum_{n\in\mathbb{N}}\langle x,v_n\rangle_Xs_nu_n.\] This shows that any \(Tx\) can be expressed as a linear combination of the set \(\lbrace u_n\rbrace\) and thus \(\lbrace u_n\rbrace\) forms an orthonormal basis for \(\overline{R(T)}\).

The finite-dimensional case

In practical applications, one deals with finite-dimensional spaces. It is not too difficult to show that any bounded linear operator defined on a Hilbert space with a finite-dimensional range is compact. This covers many “real-world” cases.

For example, let us consider any real-valued continuous-time process, defined on a reproducing kernel Hilbert space (RKHS) \(H\), observed at some discrete set of times \(\lbrace t_k\rbrace_{k=1}^N\). We define the observation operator \[O:H\to\mathbb{R}^N,\quad Of:=(f(t_1),\dots,f(t_N))^{\mathrm{T}}.\] Now \(O\) is linear and bounded. Therefore \(O\) is a compact linear operator and it has the SVD \[Ox=\sum_{n=1}^Ns_n\langle x,v_n\rangle_Hu_n,\] where the sum is only from \(1\) to \(N\), since \(u_n\) spans \(\mathbb{R}^N\). We can write this in the familiar matrix form: \[Ox=USV^{\mathrm{T}}x,\] where \(V^{\mathrm{T}}:H\to\mathbb{R}^N\) is the operator \[V^{\mathrm{T}}x=(\langle x,v_1\rangle_H,\dots,\langle x,v_N\rangle_H)^{\mathrm{T}},\] which intuitively quantifies how much of each basis vector \(v_1,\dots,v_N\) is present in the function \(x\in H\).

If \(H\) is a finite-dimensional space, say \(\mathbb{R}^M\), the representation becomes even simpler. In that case, using the typical choice \(\langle a,b\rangle_H=a^{\mathrm{T}}b\), we can define \(V\) as the matrix \(N\times M\) with the vectors \(v_1,\dots,v_N\) in its columns.

Principal component analysis (PCA) via SVD

Principal component analysis (PCA) identifies the orthogonal directions in feature space that capture the maximum variance of the data. Using SVD, we can compute PCA without explicitly constructing the covariance matrix.

Let \(X \in \mathbb{R}^{n \times m}\) be a data matrix containing \(n\) observations across \(m\) features stored in its rows. First, we mean-center \(X\) by subtracting the column-wise feature means from each row to obtain the centered matrix \(Y \in \mathbb{R}^{n \times m}\).

The \(m \times m\) sample covariance matrix \(C\) of the features is then given by: \[C = \frac{1}{n-1} Y^{\mathrm{T}} Y.\] Now, consider the SVD of \(Y\): \[Y = U S V^{\mathrm{T}},\] where \(U \in \mathbb{R}^{n \times n}\) and \(V \in \mathbb{R}^{m \times m}\) are orthogonal matrices, and \(S \in \mathbb{R}^{n \times m}\) contains the singular values. The covariance matrix is then given by \[C = \frac{1}{n-1} (V S^{\mathrm{T}} U^{\mathrm{T}})(U S V^{\mathrm{T}}) = \frac{1}{n-1} V (S^{\mathrm{T}}S) V^{\mathrm{T}} = V \left( \frac{S^2}{n-1} \right) V^{\mathrm{T}}\]

Since \(V\) is orthogonal and \(S^2/(n-1)\) is diagonal, this is just the eigen-decomposition of \(C\):

• The columns of \(V\) are the eigenvectors of \(C\), representing the principal directions in feature space \(\mathbb{R}^m\).
• The diagonal entries \(\lambda_i = \frac{s_i^2}{n-1}\) are the variance values along each principal component.
• The product \(YV = US\) gives the projected coordinates of the data points along the principal axes.