A. Linear Algebra II: Norms, Positive Definiteness & the SVD

From a first linear algebra course to the matrix tools every later module uses

Book contents · Apply this chapter to a real decision · Glossary

Almost every certificate in the modules is a statement about a matrix: this matrix is positive semidefinite (an LMI), this matrix stretches no vector by more than $L$ (a Lipschitz bound), this Gram matrix has these eigenvalues (a GP confidence bound). The modules need to compare sizes: of vectors in several norms, of matrices as maps, of quadratic forms in different directions. This primer builds those tools from $2\times2$ examples and says where each is used.

1. Vectors, Inner Products & Norms

A first course measures length by $\sqrt{x_1^2+\dots+x_n^2}$. The modules need several notions of size, inequalities that convert between them, and distances between points and sets.

First-course refresher: what a matrix calculation means

A vector is a list of coordinates. A matrix combines them by linear rules. Recover the following operations before learning the new norms; each calculation below has a direct substitution check.

Check before continuing: multiply A by $(2,-1)^\top$ to obtain $(0,-1)^\top$, solve $Ax=(4,1)^\top$ to obtain $(2,1)^\top$, and verify the eigenvector equations above. If the steps are clear, continue to the inner product; fluency returns with the local Easy exercises.

1.1 Inner products and the Cauchy–Schwarz inequality

Intuition and example. In the plane $x^\top y=\|x\|\,\|y\|\cos\theta$: the inner product measures alignment. For $x=(1,2)$, $y=(3,-1)$: $x^\top y=1$, $\|x\|=\sqrt5$, $\|y\|=\sqrt{10}$, $\theta\approx81.9^\circ$. With $M=\mathrm{diag}(1,4)$ (coordinate 2 counts double): $\langle x,y\rangle_M=3-8=-5$, $\|x\|_M=\sqrt{17}$, $\|y\|_M=\sqrt{13}$, an obtuse angle. Which vectors are long or perpendicular depends on the inner product.

Proof. If $y=0$ both sides are $0$. Otherwise, for every real $t$, $0\le\|x-ty\|^2=\|x\|^2-2t\langle x,y\rangle+t^2\|y\|^2$. This parabola in $t$ is smallest at $t^\star=\langle x,y\rangle/\|y\|^2$, where it equals $\|x\|^2-\langle x,y\rangle^2/\|y\|^2\ge0$. Equality forces $x=t^\star y$. $\square$ Only (i)–(iii) were used, so it holds for weighted inner products and for functions (Section 8).

Uses. $|a^\top d|\le\|a\|\bar d$ for a disturbance with $\|d\|\le\bar d$, and $\max_{\|u\|\le u_{\max}}b^\top u=u_{\max}\|b\|$. A logit difference is $(e_y-e_j)^\top f$ with $\|e_y-e_j\|=\sqrt2$, so an $L$-Lipschitz $f$ moves the margin by at most $\sqrt2L\|\delta\|$ (the $\sqrt2$ of Module 12's certified radius). For $H\succ0$, $u^\top H^{-1}v$ is an inner product, and for CPO's $q=g^\top H^{-1}g$, $r=g^\top H^{-1}b$, $s=b^\top H^{-1}b$ Cauchy–Schwarz reads $qs\ge r^2$.

1.2 Norms and their unit balls

Example. $x=(3,-4)$: $\|x\|_1=7$, $\|x\|_2=5$, $\|x\|_\infty=4$. The unit balls are a diamond, a disc and a square, nested as $\{\|x\|_1\le1\}\subseteq\{\|x\|_2\le1\}\subseteq\{\|x\|_\infty\le1\}$; all are convex and symmetric.

Triangle inequality for $\ell_2$ from Cauchy–Schwarz: $\|x+y\|^2=\|x\|^2+2x^\top y+\|y\|^2\le(\|x\|+\|y\|)^2$. The reverse triangle inequality $\big|\,\|x\|-\|y\|\,\big|\le\|x-y\|$ follows from $\|x\|\le\|x-y\|+\|y\|$ and symmetry: every norm is 1-Lipschitz. Stacked vectors $z=(x,u)$: $\|z\|_2^2=\|x\|_2^2+\|u\|_2^2$, $\|z\|_1=\|x\|_1+\|u\|_1$, $\|z\|_\infty=\max\{\|x\|_\infty,\|u\|_\infty\}$. A punctured ball $\{u:0\lt\|u\|_2\le1\}$ excludes the centre, where $u/\|u\|_2$ is undefined.

1.3 Hölder's inequality and dual norms

Proof of the $(1,\infty)$ case (used in Module 9): $|\sum_iu_iv_i|\le\sum_i|u_i||v_i|\le\max_j|v_j|\sum_i|u_i|$. The general case follows from Young's inequality $ab\le a^p/p+b^q/q$ (see Further Reading). Example: $w=(2,-1,3)$, $\varepsilon=0.1$: $w^\top\delta$ ranges over $\pm0.6$ on the $\ell_\infty$ ball, $\pm0.374$ on the $\ell_2$ ball, $\pm0.3$ on the $\ell_1$ ball. Since the box $[-1,1]^n$ is the $\ell_\infty$ unit ball, the hyperplane $w^\top a=b$ meets it iff $|b|\le\|w\|_1$.

1.4 Norm equivalence: the dimension factors

Proof. $\max_ix_i^2\le\sum_ix_i^2\le(\sum_i|x_i|)^2$ (cross terms are $\ge0$); $\sum_i|x_i|\cdot1\le\|x\|_2\|\mathbf 1\|_2$ by Cauchy–Schwarz; $\sum_ix_i^2\le n\max_ix_i^2$. $\square$ E.g. $x=(3,-4,12)$: $12\le13\le19\le13\sqrt3\approx22.5$.

Balls and certificates. $\{\|\delta\|_\infty\le\varepsilon/\sqrt n\}\subseteq\{\|\delta\|_2\le\varepsilon\}\subseteq\{\|\delta\|_\infty\le\varepsilon\}$, so an $\ell_2$ certificate of radius $\varepsilon$ certifies only the $\ell_\infty$ radius $\varepsilon/\sqrt n$: on MNIST ($\sqrt{784}=28$) an $\ell_2$ radius $1.58$ gives $0.056$. Lipschitz constants convert likewise, $\|f(x)-f(y)\|_\infty\le\|f(x)-f(y)\|_2\le L\|x-y\|_2\le L\sqrt n\|x-y\|_\infty$: changing the input norm costs a factor in the input dimension $n$, changing the output norm one in the output dimension $m$.

Covering numbers. A cube of side $s$ lies in the $\ell_2$ ball of radius $s\sqrt d/2$ about its centre. Cutting $[0,1]^d$ into $m^d$ cubes with $m=\lceil\sqrt d/(2\tau)\rceil$ covers it by $m^d$ balls of radius $\tau$; their centres form a $\tau$-net, and the fewest balls that suffice is the covering number $M(\tau,[0,1]^d)$. For $\tau=0.1$: $64$ balls in $d=2$, about $1.1\times10^{12}$ in $d=10$. If an $L$-Lipschitz function has margin $m_0$ at every net point, its margin is at least $m_0-L\tau$ everywhere (Primer B).

1.5 Metrics, balls and distances to sets

Examples. $\mathrm{dist}(x,\{x_1\le0\})=\max\{x_1,0\}$; for the unit disc, $\max\{\|x\|_2-1,0\}$. The $\ell_\infty$ ball $\{\|x-x_0\|_\infty\le\varepsilon\}$ is the box $|x_i-x_{0,i}|\le\varepsilon$, handled coordinate by coordinate by interval methods. Shortest-path length on a grid is a metric too. Certified balls: if $g$ is $L$-Lipschitz with $L\gt0$ and $g(z)=c\gt0$, then $g\ge0$ on $B_{c/L}(z)$, since $g(y)\ge c-L\,d(y,z)$; a union of such balls is certified. If $L=0$, the Lipschitz inequality makes $g$ constant on its domain, so the positive value at $z$ certifies the whole domain; no division by zero is needed.

Signed distance and Hausdorff distance. A signed distance to a closed set equals $\mathrm{dist}(x,C)$ on one side of the boundary and $-\mathrm{dist}(x,\text{complement})$ on the other; which side is positive is a declared convention (reachability modules put the safe side positive). It vanishes on the boundary, so "$\ell\gt0$" and "$\ell\ge0$" differ exactly there. The Hausdorff distance $d_H(A,B)=\max\{\sup_{a\in A}\mathrm{dist}(a,B),\sup_{b\in B}\mathrm{dist}(b,A)\}$ compares nonempty compact sets, e.g. $d_H([0,1],[0,2])=1$.

Isometric embeddings. $E:D\to\mathbb R^m$ is isometric if $\|E(x)-E(x')\|=d(x,x')$; e.g. a line $E(t)=x_0+tu$ with $\|u\|=1$. A kernel that depends only on distance, $k=\kappa(d(x,x'))$, is unchanged under such a map; a learned nonlinear embedding is generally not isometric and changes which points count as close.

Practice: Vectors, norms and worst-case directions

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

2. Matrix Norms & Operator Norms

A matrix is a map $x\mapsto Ax$; its most useful size is the largest factor by which it stretches a vector, i.e. the Lipschitz constant of the map.

Proofs. $|(Ax)_i|\le\sum_j|a_{ij}|\,\|x\|_\infty$, attained by $x_j=\mathrm{sign}(a_{i^\star j})$ for the worst row $i^\star$. $\|Ax\|_1\le\sum_j|x_j|\sum_i|a_{ij}|$, attained at $x=e_{j^\star}$. Each output coordinate $a_{i,:}x$ is bounded by $\|a_{i,:}\|_q\|x\|_p$ (Hölder, attainable). $\square$ Running example $A=\begin{bmatrix}2 & 2\\ -1 & 1\end{bmatrix}$: $\|A\|_1=3$, $\|A\|_\infty=4$ (at $x=(1,1)$), $\|A\|_{1\to\infty}=2$, $\|A\|_{2\to\infty}=\sqrt8$, $\|A\|_2=2\sqrt2\approx2.83$: the longest semi-axis of the ellipse $A$ makes of the unit circle (Section 5, explorer).

Examples. A column-stochastic $P$ (entries $\ge0$, columns summing to $1$) has $\|P\|_1=1$; for $0\le\gamma\lt1$, $G=(I-\gamma P)^{-1}=\sum_k\gamma^kP^k\ge0$ and $\mathbf 1^\top G=\mathbf 1^\top/(1-\gamma)$, so $\|G\|_1=1/(1-\gamma)$ exactly; for a substochastic (transient) chain $(I-P)^{-1}=\sum_tP^t$ counts expected visits. For $K\succeq0$ with $\|K\|\lt\lambda$: $(K+\lambda I)^{-1}=\lambda^{-1}I-\lambda^{-2}K+\lambda^{-3}K^2-\dots$, with remainder after the $\lambda^{-2}$ term of norm $\le\lambda^{-1}(\|K\|/\lambda)^2/(1-\|K\|/\lambda)=O(\lambda^{-3})$.

The image of a box is a zonotope. $W$ maps $\{\|\delta\|_\infty\le\varepsilon\}$ to a sum of segments $\varepsilon[-w_j,w_j]$ over its columns (a parallelogram for invertible $2\times2$ $W$). The tightest enclosing box has half-widths $\varepsilon\|w_{i,:}\|_1$, which is what interval bound propagation keeps: for $W=\begin{bmatrix}1 & 1\\ 1 & -1\end{bmatrix}$, $\varepsilon=1$, the image is the square with corners $(\pm2,0),(0,\pm2)$, area $8$, inside the box $[-2,2]^2$ of area $16$.

Practice: Matrices as maps and norm bounds

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

3. Symmetric Matrices & the Spectral Theorem

In a first course eigenvalues can be complex (a rotation has no real ones) and eigenvectors can be missing (the shear $\begin{bmatrix}1 & 1\\ 0 & 1\end{bmatrix}$ has one direction). Symmetric matrices, $A=A^\top$, have neither problem, and every certificate matrix in the modules is symmetric.

Two short proofs. Real eigenvalues: if $Av=\lambda v$, $v\ne0$ complex, then $v^*Av=\lambda v^*v$; the number $v^*Av$ equals its own conjugate $v^*A^\top v=v^*Av$, so it is real, hence so is $\lambda$ (and the real singular matrix $A-\lambda I$ has a real null vector). Orthogonality: if $Av=\lambda v$, $Aw=\mu w$, $\lambda\ne\mu$, then $\lambda v^\top w=(Av)^\top w=v^\top Aw=\mu v^\top w$, so $v^\top w=0$.

Example. $A=\begin{bmatrix}2 & 1\\ 1 & 2\end{bmatrix}$: $(2-\lambda)^2-1=0$ gives $\lambda_1=3$, $\lambda_2=1$, $q_1=(1,1)/\sqrt2$, $q_2=(1,-1)/\sqrt2$, and $3q_1q_1^\top+q_2q_2^\top=\tfrac32\begin{bmatrix}1 & 1\\ 1 & 1\end{bmatrix}+\tfrac12\begin{bmatrix}1 & -1\\ -1 & 1\end{bmatrix}=A$. In words: a symmetric matrix is a pure stretch along perpendicular axes (here by $3$ along the diagonal, by $1$ across it), with no rotation or shear.

Functions of a symmetric matrix. $f(A)=Qf(\Lambda)Q^\top=\sum_if(\lambda_i)q_iq_i^\top$ agrees with powers, polynomials and $A^{-1}$ (all $\lambda_i\ne0$), gives $(A+\eta I)^{-1}=\sum_i(\lambda_i+\eta)^{-1}q_iq_i^\top$ (all $\lambda_i+\eta\ne0$) and square roots (Section 4). Functions of one $A$ share its eigenvectors, so comparing them reduces to comparing scalars eigenvalue by eigenvalue: for $K\succeq0$ and $\eta\gt0$, $K(K+\eta I)^{-1}$ has eigenvalues $\lambda_i/(\lambda_i+\eta)\in[0,1)$.

Proof. With $y=Q^\top x$ ($\|y\|=\|x\|$): $x^\top Ax=\sum_i\lambda_iy_i^2\le\lambda_{\max}\|x\|^2$, and similarly from below. $\square$ For the example, $x=(1,0)$ gives $2$ and $x=(1,1)/\sqrt2$ gives $3=\lambda_{\max}$. A Lyapunov function $x^\top Px$ is thus squeezed between $\lambda_{\min}(P)\|x\|^2$ and $\lambda_{\max}(P)\|x\|^2$.

Complex matrices. $M$ is Hermitian if $M^*=M$. Everything above holds with $^\top$ replaced by $^*$: real eigenvalues, $M=U\Lambda U^*$ with $U$ unitary ($U^*U=I$), and $v^*Mv$ real for all complex $v$. Example: $\begin{bmatrix}2 & i\\ -i & 2\end{bmatrix}$ has eigenvalues $3,1$, and $v=(1,i)$ gives $v^*Mv=2$. Hermitian semidefiniteness ($v^*Mv\ge0$ for all $v\in\mathbb C^n$) is the complex version of Section 4.

Practice: Eigenvectors, quadratic forms and spectral limits

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

4. Positive (Semi)Definite Matrices & Quadratic Forms

A quadratic form is $q(x)=x^\top Px=\sum_{ij}p_{ij}x_ix_j$; only the symmetric part of $P$ matters, so take $P\in\mathbb S^n$. In the plane, if $p_{11}\gt0$, completing the square gives $$q(x)=p_{11}\Big(x_1+\frac{p_{12}}{p_{11}}x_2\Big)^2+\Big(p_{22}-\frac{p_{12}^2}{p_{11}}\Big)x_2^2,$$ so $q\gt0$ for all $x\ne0$ iff $p_{11}\gt0$ and $p_{11}p_{22}-p_{12}^2\gt0$. The second coefficient is our first Schur complement (Section 6).

Examples. $\begin{bmatrix}1 & 2\\ 2 & 1\end{bmatrix}$ has positive entries, yet $x=(1,-1)$ gives $-2$: indefinite. $\begin{bmatrix}1 & -1\\ -1 & 1\end{bmatrix}$ gives $(x_1-x_2)^2\ge0$: PSD and singular. For diagonal matrices the notions agree: $T=\mathrm{diag}(\lambda_i)\succeq0$ iff every $\lambda_i\ge0$. The level sets $q=1$ are ellipses (PD), parallel lines (singular PSD) or hyperbolas (indefinite): see the explorer.

Minors. A principal minor is $\det P[I,I]$, the determinant of the submatrix that keeps the same index set $I$ of rows and columns; the leading ones use $I=\{1,\dots,k\}$. A $3\times3$ matrix has $3$ leading principal minors ($p_{11}$, the top-left $2\times2$ determinant, $\det P$) but $7$ nonempty principal minors (also $p_{22}$, $p_{33}$ and the determinants for $I=\{1,3\}$ and $\{2,3\}$).

Proof of (a), (d). With $y=Q^\top x$, $x^\top Px=\sum_i\lambda_iy_i^2$: positive for all $y\ne0$ iff all $\lambda_i\gt0$ (take $x=q_i$). If $P=R^\top R$, $x^\top Px=\|Rx\|^2$, positive for $x\ne0$ iff $R$ has independent columns; conversely take $R=\Lambda^{1/2}Q^\top$. $\square$ Module 5's Gram matrix $K=\begin{bmatrix}1 & 0 & 0\\ 0 & 1 & 2\\ 0 & 2 & 5\end{bmatrix}$ has leading minors $1,1,1$, so $K\succ0$ (eigenvalues $3\mp2\sqrt2$, $1$).

4.1 Gram matrices, outer products and regularisation

For any $R$, $R^\top R\succeq0$ and $RR^\top\succeq0$ ($x^\top R^\top Rx=\|Rx\|^2$); $R^\top R\succ0$ iff $R$ has independent columns. For $\epsilon\gt0$, $x^\top(\epsilon I+R^\top R)x=\epsilon\|x\|^2+\|Rx\|^2\gt0$ ($x\ne0$), so $\epsilon I+R^\top R\succ0$ is invertible: a neural certificate uses this to get $\|(\epsilon I+R^\top R)(\xi-\xi^\ast)\|_1\gt0$ for $\xi\ne\xi^\ast$. Likewise $K+\lambda I\succ0$ for a kernel matrix $K\succeq0$ and $\lambda\gt0$, so a regularised GP never fails Cholesky in exact arithmetic. The outer product $aa^\top$ is PSD, of rank one when $a\ne0$ and $0$ when $a=0$ ($v^\top aa^\top v=(a^\top v)^2$; for $a\ne0$, eigenvalue $\|a\|^2$ on $a$, $0$ on $a^\perp$), and adding $\beta aa^\top$, $\beta\ge0$, is a PSD increase: the quadratic form grows by $\beta(a^\top v)^2$ and decreases nowhere (whether a damped loop then converges faster is a separate question).

4.2 The Loewner order, congruence and the PSD cone

Meaning and facts. $\mathrm{Var}(u^\top X)=u^\top\Sigma u$, so $\Sigma_1\succeq\Sigma_2$ means at least as much variance in every direction (a "consistent" uncertainty estimate never under-reports). $\lambda_{\min}(A)I\preceq A\preceq\lambda_{\max}(A)I$; $A\succeq B$ implies $\lambda_i(A)\ge\lambda_i(B)$ for each ordered eigenvalue (Courant–Fischer) and $\mathrm{tr}A\ge\mathrm{tr}B$; and, for $\gamma\ge0$, $\|W\|_2\le\gamma\iff W^\top W\preceq\gamma^2I$ (both say $\|Wx\|^2\le\gamma^2\|x\|^2$; squaring loses the sign of $\gamma$).

Example. $T=\begin{bmatrix}1 & -2\\ 0 & 1\end{bmatrix}$ turns $M=\begin{bmatrix}1 & 2\\ 2 & 1\end{bmatrix}$ into $T^\top MT=\mathrm{diag}(1,-3)$: it completed the square, changed the eigenvalues ($3,-1$) but not their signs (Sylvester's law of inertia). The PSD cone: $A,B\succeq0$, $\alpha,\beta\ge0\Rightarrow\alpha A+\beta B\succeq0$, so PSD matrices form a convex cone (Primer B) with the PD matrices as interior; for $2\times2$, $\begin{bmatrix}a & b\\ b & c\end{bmatrix}\succeq0\iff a,c\ge0,\ ac\ge b^2$. A linear matrix inequality $F_0+\sum_iz_iF_i\succeq0$ is a convex constraint on $z$ (Module 2).

4.3 Square roots, weighted norms and ellipsoids

Example. $\begin{bmatrix}2 & 1\\ 1 & 2\end{bmatrix}^{1/2}=\sqrt3\,q_1q_1^\top+q_2q_2^\top=\tfrac12\begin{bmatrix}\sqrt3+1 & \sqrt3-1\\ \sqrt3-1 & \sqrt3+1\end{bmatrix}\approx\begin{bmatrix}1.366 & 0.366\\ 0.366 & 1.366\end{bmatrix}$. Weighted norms. For $H\succ0$, $\|x\|_H=\sqrt{x^\top Hx}=\|H^{1/2}x\|$, and $|g^\top x|=|(H^{-1/2}g)^\top(H^{1/2}x)|\le\|g\|_{H^{-1}}\|x\|_H$, with equality for $x\propto H^{-1}g$. Steps live in the $H$-norm, gradients in the $H^{-1}$-norm: $\max\{g^\top x:\tfrac12x^\top Hx\le\delta\}=\sqrt{2\delta\,g^\top H^{-1}g}$ for $\delta\ge0$, at $x^\star\propto H^{-1}g$ if $g\ne0$ (the natural-gradient step), and $\varphi^\top V^{-1}\varphi$ is a "squared $V^{-1}$-norm". For $a\ne0$ the $H$-norm projection of $u_0$ onto $\{a^\top u\ge b\}$ moves along $H^{-1}a$: $u^\star=u_0+[b-a^\top u_0]_+H^{-1}a/(a^\top H^{-1}a)$, with the positive part $[t]_+=\max\{t,0\}$ (no move if $u_0$ is already feasible). For $a=0$ the set is all of $\mathbb R^n$ if $b\le0$ and empty if $b\gt0$. A merely PSD $H$ may be singular; the damped $H+\epsilon I\succ0$ is not.

Inversion reverses the order. If $A\succeq B\succ0$, congruence gives $B^{-1/2}AB^{-1/2}\succeq I$; its eigenvalues are $\ge1$, those of its inverse $B^{1/2}A^{-1}B^{1/2}$ are $\le1$, and congruence again gives $A^{-1}\preceq B^{-1}$.

Generalised eigenvalues: for $A\succeq0$, $M\succ0$, the nonsymmetric $AM^{-1}=M^{1/2}(M^{-1/2}AM^{-1/2})M^{-1/2}$ is similar to a symmetric PSD matrix, so its eigenvalues are real and $\ge0$, $\lambda_{\max}(AM^{-1})=\max_{x\ne0}x^\top Ax/x^\top Mx$, and for $F\ge0$: $M-FA\succeq0\iff F\lambda_{\max}(AM^{-1})\le1$.

Reading off sizes ($c=0$). Cauchy–Schwarz with $y=P^{1/2}x$ gives $\max_{\mathcal E}a^\top x=\sqrt{a^\top Sa}$, so the enclosing box has half-widths $\sqrt{S_{ii}}$ and the farthest point lies at $1/\sqrt{\lambda_{\min}(P)}=\sqrt{\lambda_{\max}(S)}$. For $P=\begin{bmatrix}2 & 1\\ 1 & 2\end{bmatrix}$: semi-axes $1/\sqrt3$ along $(1,1)$ and $1$ along $(1,-1)$, area $\pi/\sqrt3\approx1.81$. The slice $x_2=0$ has half-width $1/\sqrt{p_{11}}\approx0.707$, the projection (shadow) on $x_1$ has $\sqrt{S_{11}}=\sqrt{2/3}\approx0.816$: minimising over $x_2$ leaves the Schur complement $p_{11}-p_{12}^2/p_{22}=1.5$. Slices underestimate projections unless $p_{12}=0$.

4.4 Cholesky factorisation and linear solves

Worked example. $P=\begin{bmatrix}4 & 2 & -2\\ 2 & 5 & 1\\ -2 & 1 & 3\end{bmatrix}$: $L_{11}=2$, $L_{21}=1$, $L_{31}=-1$, $L_{22}=\sqrt{5-1}=2$, $L_{32}=(1+1)/2=1$, $L_{33}=\sqrt{3-1-1}=1$. To solve $Px=b$, $b=(2,7,5)$, never form $P^{-1}$: solve $Ly=b$ forwards, $y=(1,3,3)$, then $L^\top x=y$ backwards, $x=(2,0,3)$ (each $\approx n^2$ operations). $\log\det P=2\sum_j\log L_{jj}=\log16$.

Setting for Algorithm 1. Given observations $y_i=f(x_i)+\varepsilon_i$ at inputs $x_1,\dots,x_t$ ($y=(y_1,\dots,y_t)$), take a zero-mean GP prior on $f$ with kernel $k$ and noise $\varepsilon_i\sim\mathcal N(0,\sigma_n^2)$, $\sigma_n^2\gt0$, independent of each other and of $f$. Write $K_{ij}=k(x_i,x_j)$ for the prior covariances among the data and $k(x)=(k(x_i,x))_{i\le t}$ for those between the data and a new input $x$. The algorithm returns the posterior mean and variance of the latent function value at $x$; a new noisy observation there has variance $\sigma_t^2(x)+\sigma_n^2$. The probabilistic reading is in Primer C; here it is linear algebra, two triangular solves instead of an inverse.

Orientation. NumPy's cholesky returns lower $L$ ($P=LL^\top$), MATLAB's chol upper $R=L^\top$ ($P=R^\top R$). Module 13 writes $F=L^\top L$ with an upper-type $L$; a library's lower factor silently transposes its formulas. Conjugate gradients solve $Hx=g$ for $H\succ0$ using only products $v\mapsto Hv$ (Fisher-vector products), by minimising $\tfrac12x^\top Hx-g^\top x$ over growing subspaces until the residual $r_k=g-Hx_k$ has $\|r_k\|\le\mathrm{tol}\|g\|$; the output exactly solves $Hx=g-r_k$, and the $H$-norm error shrinks at least like $2\big((\sqrt\kappa-1)/(\sqrt\kappa+1)\big)^k$ (Section 5).

Sum-of-squares certificates. If $p(x)=z(x)^\top Qz(x)$ with a vector of monomials $z(x)$ and $Q=\sum_il_il_i^\top\succeq0$, then $p=\sum_i(l_i^\top z)^2\ge0$; e.g. $z=(1,x,x^2)$, $Q=\mathbf 1\mathbf 1^\top$ gives $x^4+2x^3+3x^2+2x+1=(1+x+x^2)^2$. Finding such a $Q$ is an SDP; the test is only sufficient (some nonnegative polynomials are not SOS).

Practice: Definiteness, factors and ellipsoid support

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

5. The Singular Value Decomposition

Eigenvalues can be complex and say little about stretching. The SVD exists for every matrix, rectangular or not, and describes exactly how it stretches.

Geometry. Read $U\Sigma V^\top x$ right to left: rotate or reflect by $V^\top$, stretch coordinate $i$ by $\sigma_i$, rotate or reflect by $U$. The unit ball becomes a filled, possibly flat ellipsoid with semi-axes $\sigma_iu_i$ ($i\le r$), and $v_1$ is the most-stretched input direction. For an injective $A$ the unit sphere maps onto the boundary of this ellipsoid within $\mathrm{range}(A)$ (an invertible $2\times2$ matrix turns the circle into an ellipse); if $A$ has a nontrivial null space, the sphere covers the whole filled ellipsoid ($A=\begin{bmatrix}1 & 0\end{bmatrix}$ maps the unit circle onto $[-1,1]$).

Consequences. (1) $\|A\|_2=\sigma_1$, because $\max_{\|x\|=1}\|Ax\|^2=\max_{\|x\|=1}x^\top A^\top Ax=\lambda_{\max}(A^\top A)$ (Rayleigh). Hence $\|A^\top\|_2=\|A\|_2$ and $\|A\|_F^2=\sum_i\sigma_i^2$; for $n\le m$, $\sigma_n=\min_{\|x\|=1}\|Ax\|$. (2) $\mathrm{range}(A)=\mathrm{span}(u_1,\dots,u_r)$ and $\mathrm{null}(A)=\mathrm{span}(v_{r+1},\dots,v_n)$: orthonormal bases of both, e.g. the nullspace basis of elimination arguments. (3) Keeping $k$ terms gives the best rank-$k$ approximation, with spectral error $\sigma_{k+1}$ (Eckart–Young).

5.1 Conditioning and the pseudo-inverse

Why Gram matrices become ill-conditioned. $K^{-1}y=\sum_i(q_i^\top y/\lambda_i)q_i$ amplifies components along small-eigenvalue directions by $1/\lambda_i$. Two inputs $0.01$ apart under a squared-exponential kernel with length scale $1$ give $K=\begin{bmatrix}1 & 0.99995\\ 0.99995 & 1\end{bmatrix}$, eigenvalues $1.99995$ and $5\times10^{-5}$, $\kappa\approx4\times10^4$; adding $\lambda=0.01$ lifts the small one to $0.01005$ and $\kappa\approx200$. Sixty close grid points are far worse.

Range conditions. $b\in\mathrm{range}(A)\iff(I-AA^+)b=0$. Module 2 writes $C^\dagger$ and $(I-CC^\dagger)B^\top=0$ (the columns of $B^\top$ lie in $\mathrm{range}(C)$), Module 9 $G^+$, Module 13 $(A^\top)^+$: the same object. Example: $\begin{bmatrix}1 & 2\\ 2 & 4\end{bmatrix}=5uu^\top$, $u=(1,2)/\sqrt5$, has rank 1 and pseudo-inverse $\tfrac15uu^\top$.

5.2 Power iteration and matrix factors

Since $\|Ax\|\le\sigma_1$ for every unit $x$, every estimate is a lower bound. If $\sigma_1\gt\sigma_2$ and $v_1^\top x_0\ne0$, the angle to $v_1$ shrinks like $(\sigma_2/\sigma_1)^{2k}$, slowly when $\sigma_2\approx\sigma_1$; if $\sigma_1=\sigma_2\gt0$, there is no single dominant $v_1$; provided $x_0$ has a nonzero projection onto the dominant right-singular subspace, $x_k$ approaches that subspace while $\hat\sigma_k\to\sigma_1$. The step breaks down only if $A^\top Ax_k=0$, i.e. $Ax_k=0$; all iterates after the first lie in $\mathrm{range}(A^\top)=\mathrm{null}(A)^\perp$, so this can only happen at $x_0$, which a random start avoids unless $A=0$. For the walkthrough matrix from $(1,0)$: $2.236,\ 2.765,\ 2.824,\ 2.8282\to2\sqrt2$. Dividing $W$ by such an estimate can leave $\|W/\hat\sigma\|_2\gt1$; guaranteed upper bounds are $\|A\|_F$ and $\sqrt{\|A\|_1\|A\|_\infty}$ ($\sqrt{12}\approx3.46$ here).

Polar factors and factor recovery. $A=(UV^\top)(V\Sigma V^\top)$ splits a square $A$ into an orthogonal factor and a PSD stretch; $UV^\top$ is the orthogonal matrix nearest to $A$ in Frobenius norm, the usual projection back onto orthogonal matrices after a gradient step. To write an invertible nonsymmetric $N$ as $VU^\top$, take $N=\tilde U\tilde\Sigma\tilde V^\top$ and $V=\tilde U\tilde\Sigma^{1/2}$, $U=\tilde V\tilde\Sigma^{1/2}$: both factors are invertible with condition number $\sqrt{\kappa(N)}$ (any $VT$, $UT^{-\top}$ also works; the SVD choice is balanced). Complex matrices have $A=U\Sigma V^*$ with unitary $U,V$, and $\sigma_{\max}(G(e^{j\omega}))$ is the gain of a transfer matrix at frequency $\omega$ (Primer D).

Practice: Singular values, least squares and power iteration

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

6. Block Matrices, Determinants, Trace & log det

6.1 Blocks, stacked vectors and sparsity

Partition into blocks of compatible sizes and multiply blocks like numbers, keeping their order: $\begin{bmatrix}A & B\\ C & D\end{bmatrix}\begin{bmatrix}x\\ u\end{bmatrix}=\begin{bmatrix}Ax+Bu\\ Cx+Du\end{bmatrix}$ (a state-space model in one matrix), and $\begin{bmatrix}x\\ y\end{bmatrix}^\top\begin{bmatrix}A & B\\ B^\top & C\end{bmatrix}\begin{bmatrix}x\\ y\end{bmatrix}=x^\top Ax+2x^\top By+y^\top Cy$. $\mathrm{blkdiag}(M_1,\dots,M_k)\succeq0$ iff every $M_i\succeq0$, which is how several LMIs become one. A selection matrix $E_i=[0\ \cdots\ I\ \cdots\ 0]$ extracts block $i$ of a stacked vector, $E_iz=z_i$, and $E_i^\top z_i$ zero-pads it back. Blocks are often named by the variables they couple ($P_{x,u}$: rows of $x$, columns of $u$); zero blocks are left blank.

Sparsity. If $L$ is block-bidiagonal (blocks only on the diagonal and first subdiagonal), $(LL^\top)_{ij}=\sum_kL_{ik}L_{jk}^\top$ vanishes unless $|i-j|\le1$: $LL^\top$ is block-tridiagonal. Conversely the block Cholesky factor of a PD block-tridiagonal matrix is block-bidiagonal (no fill-in), so layer-wise LMIs are cheap to factor. More generally, a symmetric sparsity pattern is a graph (edge $\{i,j\}$ if entry $ij$ may be nonzero); if it is chordal (every cycle of length $\ge4$ has a chord), a PSD matrix with that pattern is a sum of PSD matrices each supported on one clique (Vandenberghe & Andersen 2014, Thm. 9.2), so a large LMI splits into small ones. A tridiagonal pattern is a path; its cliques are consecutive pairs.

6.2 Block elimination and the Schur complement

Example (multiply out the factors to verify the theorem). For $M=\begin{bmatrix}2 & 1 & 0\\ 1 & 2 & 1\\ 0 & 1 & 2\end{bmatrix}$ with $A=2$: $S=\begin{bmatrix}2 & 1\\ 1 & 2\end{bmatrix}-\tfrac12\begin{bmatrix}1 & 0\\ 0 & 0\end{bmatrix}=\begin{bmatrix}1.5 & 1\\ 1 & 2\end{bmatrix}$, $\det S=2$, $\det M=2\cdot2=4$.

Proof. Complete the square: $x^\top Ax+2x^\top By+y^\top Cy=x^\top(A-BC^{-1}B^\top)x+(y+C^{-1}B^\top x)^\top C(y+C^{-1}B^\top x)$. The second term is $\ge0$ and vanishes at $y=-C^{-1}B^\top x$. $\square$ (A congruence with a block-triangular $T$.)

Why GP regression needs it. With features $\varphi(x_i)^\top$ as rows of $\Phi$, $\Phi\Phi^\top=K_t$ is the $t\times t$ kernel matrix while $\Phi^\top\Phi=\sum_i\varphi(x_i)\varphi(x_i)^\top$ lives in feature space (possibly $p=\infty$). Push-through moves between them, and Woodbury gives, with $k_t(x)=\Phi\varphi(x)$, $$\varphi(x)^\top(\lambda I+\Phi^\top\Phi)^{-1}\varphi(x)=\lambda^{-1}\big(k(x,x)-k_t(x)^\top(K_t+\lambda I)^{-1}k_t(x)\big)=\lambda^{-1}\sigma_t^2(x).$$ Gaussian conditioning uses the same Schur complement (Primer C).

6.3 Determinants and log det

Meaning. $\{x:x^\top P^{-1}x\le1\}$ has volume $\propto\det(P)^{1/2}$, so $\tfrac12\log\det P$ is a log-volume, as in the Gaussian entropy $\tfrac12\log\det(2\pi e\Sigma)$ (Primer C); products of factors become sums of logs. Pitfall: $\det(2I_3)=8$, so $\det\big(\lambda^{-1}(\lambda I_t+K_t)\big)=\lambda^{-t}\det(\lambda I_t+K_t)$.

Proof. Eliminate either corner of $\begin{bmatrix}I_m & -A\\ B & I_n\end{bmatrix}$: the Schur complement of $I_m$ gives $\det(I_n+BA)$, that of $I_n$ gives $\det(I_m+AB)$. Then $\det(V+uv^\top)=\det V\det(I+V^{-1}uv^\top)$. $\square$ Uses. (1) $\det(I_t+\lambda^{-1}K_t)=\det(I_p+\lambda^{-1}\Phi^\top\Phi)$; in infinite-dimensional feature space the right side is a Fredholm determinant, which for the rank-$\le t$ operator $\Phi^\top\Phi$ is just $\prod_i(1+\mu_i/\lambda)$ over its nonzero eigenvalues, those of $K_t$. (2) Adding a point, $V_t=V_{t-1}+\varphi_t\varphi_t^\top$, raises $\log\det$ by $\log(1+\|\varphi_t\|_{V_{t-1}^{-1}}^2)$; the sum over $t$ is the "elliptical potential" (Exercise A.6).

6.4 Trace

Lifting. $y^\top My=\langle M,yy^\top\rangle$, and $Y=yy^\top$ is PSD, rank one, $\mathrm{diag}(Y)=(y_i^2)$: a nonconvex quadratic problem in $y$ becomes linear in $Y$, and dropping "rank one" but keeping $Y\succeq0$ gives an SDP whose value upper-bounds the original maximum (the Goemans–Williamson idea).

Total variance: $\mathrm{tr}\,\Sigma=\sum_i\mathrm{Var}(X_i)=\mathbb E\|X-\mathbb EX\|^2$, while $\mathrm{tr}(\Sigma^{1/2})=\sum_i\sqrt{\lambda_i}$ is the sum of marginal standard deviations only for diagonal $\Sigma$; for $\begin{bmatrix}1 & 0.8\\ 0.8 & 1\end{bmatrix}$ it is $\sqrt{1.8}+\sqrt{0.2}\approx1.79$, not $2$.

Tight frames: unit $x_1,\dots,x_T\in\mathbb R^d$ with $\sum_tx_tx_t^\top=cI$ force $c=T/d$ (take traces: $T=cd$); three unit vectors at $120^\circ$ give $\tfrac32I$.

Practice: Block algebra, Schur complements and updates

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

7. Orthogonal, Skew-Symmetric & Structured Matrices

7.1 Orthogonal matrices, reflections and projections

Householder reflections. For a unit vector $v$, $H=I-2vv^\top$ reflects across $v^\perp$ ($Hv=-v$, $Hw=w$ for $w\perp v$); it is symmetric, orthogonal and its own inverse, while $I-vv^\top$ is the (singular) orthogonal projector onto $v^\perp$. For $v=(1,1)/\sqrt2$, $H=\begin{bmatrix}0 & -1\\ -1 & 0\end{bmatrix}$ maps $(1,0)$ to $(0,-1)$. Since $\det$ is continuous with values $\pm1$ on orthogonal matrices, no continuous path joins a reflection to a rotation: orthogonal matrices form two separate pieces, and a continuous parameterisation from $\mathbb R^p$ reaches only one (Primer 0).

7.2 Skew-symmetric matrices and the Cayley transform

Proof of orthogonality. $I\pm A$ are polynomials in $A$, so they commute; and $XY=YX$ with $Y$ invertible gives $Y^{-1}X=XY^{-1}$ (multiply by $Y^{-1}$ on both sides), so $I-A$ also commutes with $(I+A)^{-1}$. With $(I+A)^\top=I-A$: $Q^\top Q=(I-A)^{-1}(I+A)(I-A)(I+A)^{-1}=(I-A)^{-1}(I-A)(I+A)(I+A)^{-1}=I$. $\square$ Example: $A=\begin{bmatrix}0 & -\tfrac12\\ \tfrac12 & 0\end{bmatrix}$ gives $Q=\begin{bmatrix}0.6 & 0.8\\ -0.8 & 0.6\end{bmatrix}$, a rotation by $-53.1^\circ$. Optimising over orthogonal matrices: differentiating $Q(t)^\top Q(t)=I$ gives $\dot Q^\top Q+Q^\top\dot Q=0$, so the tangent directions at $Q$ are $QA$ with $A$ skew; methods move along them ($Q\,\mathrm{Cayley}(tA)$, $Qe^{tA}$) or take a plain step and project back with the polar factor $UV^\top$ (Section 5).

7.3 Diagonal, permutation, Toeplitz and circulant matrices

$D=\mathrm{diag}(d)$ scales coordinates, $\|D\|_p=\max_i|d_i|$, and $D\succeq0$ iff $d\ge0$. A permutation matrix $\Pi$ is orthogonal, and $\Pi^\top M\Pi$ keeps eigenvalues and definiteness. $\ell_p$ norms are permutation invariant; GroupSort and MaxMin permute entries within groups, so where differentiable their Jacobian is a permutation matrix: orthogonal, norm-preserving.

Example. $h=(1,2)$: $T_4=\begin{bmatrix}1 & 0 & 0 & 0\\ 2 & 1 & 0 & 0\\ 0 & 2 & 1 & 0\\ 0 & 0 & 2 & 1\end{bmatrix}$. The matrix grows with the signal, the filter does not: $\|T_N\|_2=2.41,\,2.83,\,2.95,\,3.00$ for $N=2,4,8,32$, approaching the frequency-response gain $\max_\omega|1+2e^{-i\omega}|=3$. The $4\times4$ circulant has eigenvalues $3,\,1-2i,\,-1,\,1+2i$ and norm $3$; for real $h$, $\hat h_{N-k}=\overline{\hat h_k}$ (conjugate symmetry), so per-frequency operations can return real outputs. Vector-valued signals give block-Toeplitz matrices; finite memory makes them banded. Time invariance is the Toeplitz structure, finite memory the band.

7.4 Kronecker products and flattening

Examples. $I_2\otimes\begin{bmatrix}1 & 2\\ 3 & 4\end{bmatrix}=\begin{bmatrix}1 & 2 & 0 & 0\\ 3 & 4 & 0 & 0\\ 0 & 0 & 1 & 2\\ 0 & 0 & 3 & 4\end{bmatrix}$ is block diagonal; $\begin{bmatrix}1 & 2\\ 3 & 4\end{bmatrix}\otimes I_2=\begin{bmatrix}1 & 0 & 2 & 0\\ 0 & 1 & 0 & 2\\ 3 & 0 & 4 & 0\\ 0 & 3 & 0 & 4\end{bmatrix}$ interleaves. So $I_N\otimes Q$ applies $Q$ to each of $N$ pixels, $\mathrm{diag}(\lambda)\otimes I_g$ gives each group of $g$ neurons one multiplier, and $z^\top(\mathrm{diag}(\gamma)\otimes\mathbf 1\mathbf 1^\top)z=\sum_i\gamma_i(\mathbf 1^\top z_i)^2$ collects group sums. Flattening order matters: $c$ channels at $N$ pixels flattened pixel by pixel turn $\sum_px_p^\top Qx_p$ into $\mathrm{vec}^\top(I_N\otimes Q)\mathrm{vec}$; channel by channel it is $Q\otimes I_N$. Complex counterparts: unitary ($U^*U=I$) and skew-Hermitian ($A^*=-A$) matrices; the Cayley transform of the latter is unitary.

Practice: Orthogonal maps, projections and Cayley transforms

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

8. Functions as Vectors: Inner Product & Hilbert Spaces

Functions $f:D\to\mathbb R$ add and scale pointwise, so they are vectors. On a finite $D=\{x_1,\dots,x_n\}$, $f$ is the vector $(f(x_1),\dots,f(x_n))$; on an interval there are infinitely many independent functions ($1,x,x^2,\dots$). Kernel methods (Modules 3–6) do linear algebra in such infinite-dimensional spaces.

Examples. On $[0,1]$: $\langle1,x\rangle=\tfrac12$, $\|x\|=1/\sqrt3$, angle $30^\circ$. In the Fourier basis $1,\sqrt2\cos2\pi kx,\sqrt2\sin2\pi kx$, $3+4\sqrt2\cos2\pi x$ has coefficients $(3,4,0,\dots)$ and norm $5$: Pythagoras for functions.

What infinite dimensions add. (1) A Hilbert space is an inner-product space that is complete: every Cauchy sequence converges inside it. Polynomials are not complete in $L^2$ on $[0,1]$ (the Taylor polynomials of $e^x$ converge to $e^x$, which is no polynomial); the completion adds such limits, as the reals complete the rationals. (2) Norms are no longer equivalent: $x^n$ on $[0,1]$ has sup norm $1$ but $L^2$ norm $1/\sqrt{2n+1}\to0$. (3) Projection still works: the best approximation of $f$ in $\mathrm{span}(g_1,\dots,g_t)$ solves $G\alpha=b$, $G_{ij}=\langle g_i,g_j\rangle$, $b_i=\langle g_i,f\rangle$, and $\|f\|^2=\|f_\parallel\|^2+\|f_\perp\|^2$.

Functionals and adjoints. On $\mathbb R^n$ every linear functional is $x\mapsto a^\top x$; on a Hilbert space every bounded one ($|L(f)|\le C\|f\|$) is $\langle g_L,f\rangle$ for a unique $g_L$ (Riesz). Point evaluation $f\mapsto f(x)$ is not even defined on $L^2(dx)$ (changing $f$ at one point changes nothing), and on continuous functions with the $L^2$ norm it is unbounded: the spikes $f_n(y)=\max\{1-n|y-x|,0\}$ have $f_n(x)=1$ but $\|f_n\|_{L^2}^2\le\tfrac2{3n}\to0$. (For a discrete $\mu$ with a point mass $p\gt0$ at $x$ it is bounded, $|f(x)|\le\|f\|_{L^2}/\sqrt p$.) The spaces where it is bounded are the RKHSs below. The adjoint is defined by $\langle Tv,w\rangle=\langle v,T^*w\rangle$ ($T^*=T^\top$ on $\mathbb R^n$). For the design operator $\Phi g=(\langle\varphi(x_i),g\rangle)_{i\le t}$, $v^\top\Phi g=\langle\sum_iv_i\varphi(x_i),g\rangle$ gives $\Phi^*v=\sum_iv_i\varphi(x_i)$; $\Phi\Phi^*$ is the $t\times t$ matrix $[\langle\varphi(x_i),\varphi(x_j)\rangle]$, and $\Phi^*\Phi=\sum_i\varphi(x_i)\langle\varphi(x_i),\cdot\,\rangle$ is a positive operator of rank $\le t$, the operator version of $\sum_i\varphi_i\varphi_i^\top$.

8.1 Feature maps and kernels

Examples. $\varphi(x)=(1,\sqrt2x,x^2)$ gives $k(x,x')=(1+xx')^2$ (at $x=1$, $x'=2$ both sides are $9$): a kernel computes feature inner products without forming features. The squared-exponential kernel $\sigma_f^2\exp(-\|x-x'\|^2/(2\ell^2))$ has infinitely many features; $\sigma_f^2=k(x,x)$ is the prior variance and $\ell$ the length scale: $k/\sigma_f^2=0.61,\,0.14,\,0.011$ at distances $\ell,2\ell,3\ell$. Building kernels: sums, positive multiples and products of kernels are kernels. A product kernel on pairs, $k_1(\theta,\theta')k_2(t,t')$, has features $\varphi_1(\theta)\otimes\varphi_2(t)$, and its Gram matrix $K_1\circ K_2$ (entrywise product) is a principal submatrix of the PSD $K_1\otimes K_2$. The Wiener kernel $\min(t,t')=\int_0^\infty\mathbf 1[s\le t]\,\mathbf 1[s\le t']\,ds$ ($t,t'\ge0$) has the feature functions $\mathbf 1_{[0,t]}$.

8.2 The reproducing kernel Hilbert space (the bridge to Module 3)

The RKHS norm measures wiggliness. If $K\succ0$ (e.g. the SE kernel at distinct inputs), the smallest-norm $f$ with $f(x_i)=y_i$ has $\|f\|_k^2=y^\top K^{-1}y$. (For singular $K$ an interpolant exists iff $y\in\mathrm{range}(K)$, and the smallest squared norm is then $y^\top K^+y$: for the constant kernel $k\equiv1$ at two inputs, $K=\mathbf 1\mathbf 1^\top$, the data $(1,1)$ cost $1$ and $(1,-1)$ cannot be interpolated.) SE kernel, $\ell=1$, two points $0.1$ apart ($k=0.995$): $y=(1,1)$ needs $\|f\|_k\approx1.00$, $y=(1,-1)$ needs $\|f\|_k\approx20.0$. For a smooth kernel such as this one, a bound $\|f\|_k\le B$ rules out fast oscillation (the assumption behind GP confidence bounds); the alternating data sit on the small-eigenvalue direction of $K$ (Section 5).

Intuition. On a finite $D$ with $K=Q\Lambda Q^\top\succ0$, $\|f\|_k^2=f^\top K^{-1}f=\sum_j(q_j^\top f)^2/\lambda_j$, so small-eigenvalue components are expensive: the RKHS unit ball is an ellipsoid with semi-axes $\sqrt{\lambda_j}$ along the eigenvectors, and fast eigenvalue decay makes it thin in all but a few directions. For a smooth kernel such as the SE kernel the expensive components are the rough, oscillating ones; in general smoothness depends on the eigenfunctions as well as the eigenvalues. Stationary kernels are a separate whole-space case: if $k(x,x')=\kappa(x-x')$ on $\mathbb R^d$ with $\kappa$ continuous and integrable, and $\hat g(\omega)=\int g(x)e^{-i\omega^\top x}dx$, Fourier modes replace eigenfunctions: if $\hat\kappa\gt0$, $\mathcal H_k$ consists of the continuous $f\in L^2$ with $\|f\|_k^2=(2\pi)^{-d}\int|\hat f(\omega)|^2/\hat\kappa(\omega)\,d\omega\lt\infty$ (Kanagawa et al., arXiv 2018, Thm. 2.4; smooth kernels have $\hat\kappa$ tiny at high frequencies, which makes rough $f$ expensive). On a compact $D$ this does not apply directly: the RKHS norm of $f$ is the smallest whole-space norm among its extensions to $\mathbb R^d$. Operator eigenvalues are not those of a sampled $K_t$, though $\lambda_j(K_t)/t\approx\lambda_j$ for large $t$ if $\mu$ is a probability measure and the inputs are i.i.d. from $\mu$ (for a general finite $\mu$, sampling from $\mu/\mu(D)$ gives the limit $\lambda_j/\mu(D)$).

Distances between distributions: the mean embedding $m_P=\mathbb E_{X\sim P}\varphi(X)$ gives the maximum mean discrepancy $\mathrm{MMD}(P,Q)=\|m_P-m_Q\|_k$, with $\mathrm{MMD}^2=\mathbb Ek(X,X')-2\mathbb Ek(X,Y)+\mathbb Ek(Y,Y')$ ($X,X'\sim P$, $Y,Y'\sim Q$ independent); for the SE kernel it vanishes only if $P=Q$, the basis of kernel two-sample tests.

Practice: Functions as vectors, kernels and RKHS norms

Work from Easy to Medium before trying Hard. Use the review link for a missing definition, then try again with the explanation closed.

Interactive: Quadratic Forms, Eigenvectors & the SVD

Set the entries of $A=\begin{bmatrix}a_{11} & a_{12}\\ a_{21} & a_{22}\end{bmatrix}$. Left: $A$ maps the unit circle (grey) to an ellipse (blue); dashed teal: right singular vectors $v_1,v_2$; blue: their images $\sigma_iu_i$; purple: real eigen-directions; orange: the power iterate $x_k$ and $Ax_k$. Right: level sets $x^\top Sx=c$, $S=\tfrac12(A+A^\top)$, for $c=\pm\tfrac12,\pm1,\pm2$ (blue positive, orange negative, $c=1$ thick). The readouts use closed-form formulas evaluated in floating-point arithmetic (the determinant and eigenvalue discriminant come from the integer slider ticks, so singular matrices are recognised exactly); the level-set curves are approximated numerically on a grid.

What to look for. Walkthrough: no real eigenvectors ($\rho=2$), yet two singular directions with $\sigma_1=2\sqrt2$; the power estimate climbs to $2.83$ from below. Symmetric PD: eigen- and singular vectors coincide, the level sets are ellipses with semi-axes $1/\sqrt{\lambda_i}$, Cholesky succeeds. Indefinite: hyperbolas; Cholesky fails at pivot 2 although all entries are positive. Singular PSD: the ellipse collapses to a segment (rank 1), level sets become lines. Nilpotent: $\rho(A)=0$ but $\|A\|_2=2$, and power iteration from $(1,0)\perp v_1$ breaks down at once: $(1,0)\in\mathrm{null}(A)$, so $A^\top A(1,0)=0$ and there is no next iterate (the readout says so). Rotation: $\sigma_1=\sigma_2=1$ (every direction is a dominant one, so $v_1$ is not unique), and $S=0.6I$ is PD although $A$ is not symmetric.

Exercises

Ready to use these matrix tools?

Before Primer B, redo the Easy and Medium exercises on vector norms, eigenvectors, definiteness and solves, and singular values. You should be able to explain why a symmetric PSD matrix has a nonnegative quadratic form, and why an eigenvalue is different from a singular value.

Before Module 2, check the Schur complement practice. Before Module 3, check the feature and RKHS practice. A missed calculation points to that section's local refresher; a missed assumption points to its definition or pitfall box. Next: Primer B: derivatives and optimisation.

Further Reading

Flashcards