Random Matrix Theory

Motivation

Random matrices with independent entries appear whenever a system is measured through many noisy, roughly independent channels: shot noise in a sensor array, edges in an Erdős–Rényi-type random graph, or the design matrix of a linear model with independent covariates. A basic question about any such matrix $A$ is how far it can stretch a vector — its operator norm $\|A\|$ — since this single number controls the stability of every linear statistic computed from $A$: least-squares estimates, spectral clustering, covariance estimation, and random projections all reduce, at some point, to bounding $\|A\|$.

The theory traces to Marchenko and Pastur's 1967 asymptotic law for the spectrum of large random matrices, and to Bai and Yin's 1988 almost-sure limit $\|A\|/\sqrt n \to 2$ for $n \times n$ matrices with i.i.d. mean-zero, unit-variance entries. Those results are asymptotic: they say what happens as the dimension $n \to \infty$, for a fixed matrix shape. The result formalized here, Theorem 4.4.5 of Vershynin's High-Dimensional Probability (2018) (DOI 10.1017/9781108231596), belongs to the more recent non-asymptotic strand of the theory: it gives an explicit, dimension-free bound that holds at every fixed $m, n$, with an explicit failure probability — the form of statement needed for finite-sample guarantees in statistics and data science, rather than limiting behavior.

Setting

Let $A$ be an $m \times n$ real matrix. Equip $\mathbb R^n$ and $\mathbb R^m$ with the Euclidean norm $\|\cdot\|_2$; $A$ acts as a linear map $\ell_2^n \to \ell_2^m$. Its operator norm (§4.1.2) is

$$\|A\| \;:=\; \max_{x \in S^{n-1}} \|Ax\|_2 ,$$

the largest factor by which $A$ can stretch a unit vector; equivalently, the largest singular value of $A$.

A real random variable $X$ is sub-gaussian with the mission's own convention (matching the Orlicz $\psi_2$ norm this series already carries as a published definition, HighDimProb.Concentration.subgaussianNorm) if

$$\|X\|_{\psi_2} \;:=\; \inf\{t > 0 : \mathbb E \exp(X^2/t^2) \le 2\} \;<\; \infty .$$

Bounded random variables and Gaussians are sub-gaussian; a Bernoulli($p$) variable and a $\pm 1$-valued coin flip both qualify, which is why the theorem below directly covers random matrices with i.i.d. Rademacher or Gaussian entries as special cases.

The mission's proof technique is the ε-net argument, developed in §4.2 and used nowhere before this chapter of the book: a metric space $(T,d)$, a subset $K \subseteq T$, and $\varepsilon > 0$ give rise to an ε-net $N \subseteq K$ — a finite set such that every point of $K$ is within $\varepsilon$ of some point of $N$ — and a covering number $N(K, d, \varepsilon)$, the smallest cardinality of such a net. The technique reduces a statement that must hold uniformly over an infinite (compact) set to a statement about finitely many points, paid for by a union bound whose cost is controlled by the covering number.

Formalization targets

Goal (Theorem 4.4.5)

$$\exists\, C > 0 \;\text{such that}\; \forall\, t > 0,\quad \mathrm{Prob}\bigl\{\, \|A\| \le CK(\sqrt m + \sqrt n + t) \,\bigr\} \;\ge\; 1 - 2\exp(-t^2),$$

for any $m \times n$ random matrix $A$ with independent, mean-zero, sub-gaussian entries $A_{ij}$ and $K = \max_{i,j} \|A_{ij}\|_{\psi_2}$. This is the weakest stable form of the bound — it fixes no numerical value for $C$, only its existence and absoluteness (independence from $m$, $n$, $A$, $t$), so later refinements of the constant do not invalidate it.

Significance

The result itself. The bound $\|A\| \lesssim \sqrt m + \sqrt n$ is sharp up to the constant: for entries of unit variance, $\mathbb E\|A\| \ge \tfrac14(\sqrt m + \sqrt n)$ for large $m, n$ (the book's Exercise 4.4.7), so no non-asymptotic bound of this shape can be improved beyond constants. It is the entry point to the rest of the book's random matrix theory: Corollary 4.4.8 specializes it to symmetric matrices, and it underlies the community-detection (§4.5) and covariance-estimation (§4.7) applications later in the same chapter, neither of which is part of this mission.

Formalizing it. The theorem is a classical, fully proved result; nothing about its truth is open. What this mission contributes is a machine-checked formal statement — together with the two pieces of chapter infrastructure its own textbook proof names by number (Corollary 4.2.13, Exercise 4.4.3(a)) — and a third, self-contained application of the same covering-number machinery (Theorem 4.3.5) that exercises the shared IsEpsNet/coveringNumber definitions on a different metric space (the Hamming cube), independently of the Euclidean case. Mathlib and the Prove2Me platform currently have no ε-net, covering-number, or packing-number infrastructure (checked by q=random matrix, q=operator norm, q=covering number, q=net on the platform, and by filename search in Mathlib): this mission is the first to introduce it, restated inside its own namespace since it is not otherwise available to build on.

Difficulty

The obvious first approach is to bound $\|A\| = \max_{x \in S^{n-1}} \|Ax\|_2$ directly by union-bounding a concentration inequality over the sphere $S^{n-1}$. This fails outright: $S^{n-1}$ is infinite (indeed uncountable) for $n \ge 2$, so no union bound over its points can converge — the naive approach gives $\infty \cdot (\text{tail probability})$. The ε-net argument is the fix, but it is not just "discretize and hope": the reduction from the sphere to a finite net (quadratic_form_on_net, Exercise 4.4.3(a)) loses a multiplicative factor $1/(1-2\varepsilon)$ that must be tracked, and the net's cardinality (covering_numbers_of_euclidean_ball_and_sphere, Corollary 4.2.13) is exponential in the dimension ($9^n$ at $\varepsilon = 1/4$) — so the per-point tail probability from Hoeffding-type concentration must itself decay fast enough (quadratically in the exponent) to survive multiplying by $9^{m+n}$ many points. Getting the union bound to close requires choosing the threshold $u$ in the tail bound proportionally to $\sqrt m + \sqrt n + t$, not to $t$ alone — the $\sqrt m + \sqrt n$ term is exactly what pays for the net's exponential size.

Formalization scope

$A$ is represented as Ω → Matrix (Fin m) (Fin n) ℝ; its entries A ω i j are the individual real random variables. Independence of the $mn$ entries is iIndepFun over the index type Fin m × Fin n; mean-zero is the vanishing of each entry's Bochner integral. The operator norm is the norm of the associated continuous linear map between EuclideanSpace ℝ (Fin n) and EuclideanSpace ℝ (Fin m) (matrixOpNorm, every linear map between finite-dimensional normed spaces being automatically continuous), matching the book's maxₓ∈Sⁿ⁻¹ ‖Ax‖₂ exactly. The sub-gaussian norm $K$ reuses this series' own published definition, HighDimProb.Concentration.subgaussianNorm, rather than a re-derivation. Covering numbers (coveringNumber) are restricted to finite (Finset) ε-nets, the only kind this chapter uses; this is a deliberate restriction, not a general-purpose covering-number formalization, and is disclosed as such. A hard-coded numeral for $C$, or an unquantified "with high probability" in place of the explicit failure probability $2\exp(-t^2)$, would each trivialize the statement and is ruled out: $C$ is existentially bound ahead of every other quantifier, and $t > 0$ is a free parameter with its own explicit bound, exactly as the book states it.

Definitions reusable beyond this mission: IsEpsNet and coveringNumber are stated for a general PseudoMetricSpace and apply unchanged to any later chapter's covering-number needs (e.g. Chapter 8's VC-dimension covering numbers), though per this series' rule that drafts cannot import drafts, a later chunk would restate rather than import them until this mission is published. matrixOpNorm is likewise chapter-agnostic. Contributions completing the sorry proofs of any of the four theorem items are welcome and independent of one another; the covering number and net-reduction items (covering_numbers_of_euclidean_ball_and_sphere, quadratic_form_on_net) are the standard prerequisites for the goal's own volumetric/ε-net proof.

Selected references

• Vershynin, R. High-Dimensional Probability: An Introduction with Applications in Data Science. Cambridge University Press, 2018. DOI 10.1017/9781108231596
• Bai, Z. D., Yin, Y. Q. "Necessary and sufficient conditions for almost sure convergence of the largest eigenvalue of a Wigner matrix." Annals of Probability 16 (1988), 1729–1741.
• Marchenko, V. A., Pastur, L. A. "Distribution of eigenvalues for some sets of random matrices." Mathematics of the USSR-Sbornik 1 (1967), 457–483.

Motivation

Estimating a $d\times d$ covariance matrix $\Sigma$ from $n$ samples is a canonical high-dimensional problem: the natural estimator, the sample covariance $\hat\Sigma$, is consistent in operator norm only when $n\gtrsim d$, which fails outright in the regime $d\gg n$ common to modern applications. When $\Sigma$ is additionally known to be sparse — few nonzero entries per row — a simple fix restores consistency even when $d\gg n$: threshold every entry of $\hat\Sigma$ below a data-dependent level to zero. This mission formalizes the matrix concentration machinery behind that fix and the resulting guarantee, following Wainwright, High-Dimensional Statistics: A Non-Asymptotic Viewpoint (Cambridge University Press, 2019), Chapter 6.

Setting

For a random matrix $Q$, the matrix variance is $\mathrm{var}(Q):=\mathbb E[Q^2]-(\mathbb E[Q])^2$, and the Loewner order $A\preceq B$ on symmetric matrices means $B-A$ is positive semidefinite. A zero-mean symmetric random matrix $Q$ satisfies Bernstein's condition (Definition 6.10) with parameter $b>0$ if $\mathbb E[Q^j]\preceq\tfrac12 j!\,b^{j-2}\mathrm{var}(Q)$ for $j=3,4,\dots$ — the matrix analogue of the scalar Bernstein condition of Chapter 2. The operator (spectral) norm $|\!|\!|M|\!|\!|_2$ is $M$'s largest singular value.

Given a threshold $\lambda>0$, the hard-thresholding operator is $T_\lambda(u):=u\cdot \mathbb 1[|u|>\lambda]$, extended entrywise to matrices (Eq. (6.52)). A covariance matrix $\Sigma$'s sparsity pattern is captured by its adjacency matrix $A_{j\ell}:=\mathbb 1[\Sigma_{j\ell}\neq0]$ (p. 181); $|\!|\!|A|\!|\!|_2\le d$ always, with equality only when $\Sigma$ has no zero entries, and $|\!|\!|A|\!|\!|_2\le s$ whenever $\Sigma$ has at most $s$ nonzero entries per row.

Formalization targets

Goal — Theorem 6.23 (thresholding-based covariance estimation)

Let $\{x_i\}_{i=1}^n$ be i.i.d. zero-mean random vectors with covariance $\Sigma$, each coordinate sub-Gaussian with parameter at most $\sigma$. If $n>\log d$, then for any $\delta>0$, the thresholded sample covariance $T_{\lambda_n}(\hat\Sigma)$ with $\lambda_n/\sigma^2 = 8\sqrt{\log d/n}+\delta$ satisfies

$$\mathbb P\big[\,|\!|\!|T_{\lambda_n}(\hat\Sigma)-\Sigma|\!|\!|_2 \ge 2|\!|\!|A|\!|\!|_2\lambda_n\,\big] \;\le\; 8e^{-\frac{n}{16}\min\{\delta,\delta^2\}}.$$

The error scales with the graph sparsity $|\!|\!|A|\!|\!|_2$, not with $d$ directly — the whole point of thresholding when the ambient dimension is much larger than the sample size.

Milestone — Theorem 6.17 (the matrix Bernstein bound)

For independent, zero-mean, symmetric random matrices $\{Q_i\}$ satisfying Bernstein's condition with parameter $b$,

$$\mathbb P\Big[\frac1n\Big|\!\Big|\!\Big|\sum_{i=1}^n Q_i\Big|\!\Big|\!\Big|_2 \ge \delta\Big] \;\le\; 2\,\mathrm{rank}\Big(\sum_{i=1}^n\mathrm{var}(Q_i)\Big)\exp\Big(-\frac{n\delta^2}{2(\sigma^2+b\delta)}\Big).$$

This is the general matrix concentration tool the whole chapter builds toward; Theorem 6.23 is one of its corollaries.

Milestone — Eq. (6.54) (the deterministic thresholding bound)

For any $\lambda_n$ with $\|\hat\Sigma-\Sigma\|_{\max}\le\lambda_n$, $|\!|\!|T_{\lambda_n}(\hat\Sigma)-\Sigma|\!|\!|_2\le2|\!|\!|A|\!|\!|_2\lambda_n$ — a purely deterministic fact, with no probability involved, that reduces Theorem 6.23's proof to a single probabilistic input: controlling $\|\hat\Sigma-\Sigma\|_{\max}$.

Significance

Theorem 6.17 is the workhorse of the whole chapter: besides Theorem 6.23, it also underlies Corollary 6.20's operator-norm bound for the unstructured sample covariance matrix (used in turn for the two example ensembles in Section 6.4.5), by taking $Q_i:=x_ix_i^T-\Sigma$. Theorem 6.23 itself is the standard justification for thresholding-based covariance estimators used throughout high-dimensional statistics whenever the sparsity pattern of $\Sigma$, though unknown, is believed to be structured.

Formalizing it. No faithful prior art exists on the platform for the matrix Bernstein bound or covariance thresholding (a fresh search for "matrix Bernstein," "Bernstein condition matrix," "thresholding covariance," and "sparse covariance" returned no hits). The existing RademacherWigner.* items formalize spectral-edge and empirical-spectral-distribution bounds for the Wigner ensemble specifically (a symmetric matrix with i.i.d. entries) — a narrower, different random-matrix model from the general Bernstein-condition matrices this chapter treats, and not reused here. All three theorems are drafted as open goals (:= by sorry).

Difficulty

The naive approach to a matrix tail bound — apply the scalar Chernoff/Bernstein technique directly to the operator norm — fails because the operator norm is not a linear functional of $Q$, so the scalar moment generating function bound does not translate directly. The resolution (Lemma 6.13, not itself part of this mission) instead bounds the trace of the matrix exponential of the sum, which is linear-algebraically tractable via the Golden–Thompson-type inequality and the union bound over the (at most $\mathrm{rank}(\bar V)$) nonzero eigenvalue directions of the aggregate variance $\bar V=\sum_i\mathrm{var}(Q_i)$ — this is exactly why the bound's prefactor is $\mathrm{rank}(\bar V)$, not the ambient dimension $d$: a low-rank aggregate variance (e.g. from a highly structured collection of matrices) yields a much tighter bound than a naive union bound over all $d$ eigen-directions would. Theorem 6.23's own difficulty is entirely in the reduction to Theorem 6.17 plus Eq. (6.54): checking that $Q_i:=x_ix_i^T-\Sigma$ satisfies the Bernstein condition with the stated parameters, and that $\|\hat\Sigma-\Sigma\|_{\max}$ concentrates at the stated rate via a union bound over the $O(d^2)$ matrix entries.

Formalization scope

"Symmetric" is realized via Mathlib's Matrix.IsSymm; the Loewner order via (B-A).PosSemidef. The matrix moment generating function itself (ΨQ(λ):=E[e^{λQ}]) is not formalized, since none of this mission's three theorem statements need it directly — Bernstein's condition (Definition 6.10) is stated purely in terms of polynomial moments $\mathbb E[Q^j]$, avoiding the substantial extra machinery (Mathlib's NormedSpace.exp for matrices) that a faithful matrix-exponential definition would require, at no cost to faithfulness for the three statements actually drafted.

Independence of $\{Q_i\}$/$\{x_i\}$ is realized via Mathlib's iIndepFun; this required adding a local MeasurableSpace (Matrix (Fin d) (Fin d) ℝ) instance in the workspace, since Matrix is a def, not an abbrev, over the underlying Pi type, so Lean's instance search does not automatically find the Pi-type MeasurableSpace instance for it. "Identically distributed" is realized via Mathlib's IdentDistrib against a fixed reference index. "Each component sub-Gaussian" is restated locally as an explicit MGF bound at the reference index, per this book's cross-chapter rule against importing another chapter's draft definitions.

If a statement admits a trivializing formalization: the rank(...) prefactor of Theorem 6.17 is kept exactly, not replaced by the ambient dimension $d$ (which would be a strictly weaker, non-trivializing but unfaithful simplification, since the bound's whole point is that rank(...) ≤ d can be much smaller). ‖·‖_max and opNorm at d=0 return Mathlib's junk value 0 (harmless: no matrix entries or singular values exist at that degenerate size either).

Out of scope for this mission: Theorem 6.15 (the sub-Gaussian-tail matrix concentration result that precedes Theorem 6.17 in the same section) and Corollary 6.20 (the unstructured covariance-estimation corollary) — both natural follow-on work, cut for this chunk's time budget; see STATUS.md.

Selected references

• M. J. Wainwright, High-Dimensional Statistics: A Non-Asymptotic Viewpoint, Cambridge University Press, 2019. DOI: 10.1017/9781108627771. Chapter 6.
• J. A. Tropp, "User-friendly tail bounds for sums of random matrices," Foundations of Computational Mathematics, 12(4):389–434, 2012.

Motivation

Many questions in high-dimensional probability reduce to bounding the expected supremum of a random process $(X_t)_{t\in T}$ — the maximum, over an entire indexed family of random variables, of how large any one of them can get. When $T$ is finite this is routine (a union bound over $|T|$ terms suffices), but the interesting cases have $T$ infinite, even uncountable: a supremum over a continuum of test functions, a norm expressed as a supremum over a sphere, or an empirical process indexed by a whole class of functions. A naive union bound is unusable here, since $|T|$ is infinite.

R. M. Dudley's 1967 entropy bound (R. M. Dudley, The sizes of compact subsets of Hilbert space and continuity of Gaussian processes, Journal of Functional Analysis 1 (1967), 290–330) resolved this for Gaussian processes, controlling the expected supremum purely in terms of the metric entropy of $T$ — how many balls of radius $\varepsilon$ are needed to cover $T$, at every scale $\varepsilon$. The technique behind the proof, chaining, builds a sequence of increasingly fine finite approximations to $T$ and telescopes the resulting bounds; it is one of the central tools of the field, reused throughout empirical process theory, statistical learning theory (via Vapnik-Chervonenkis theory), and non-asymptotic random matrix theory. This mission formalizes the chapter's generalization of Dudley's bound beyond Gaussian processes, to any process with sub-gaussian increments, together with the purely combinatorial Sauer-Shelah lemma that the chapter's applications to statistical learning theory build on.

Setting

Fix a probability space $(\Omega,\mathcal F,P)$. A random process is a family $(X_t)_{t\in T}$ of real random variables on $(\Omega,\mathcal F,P)$ indexed by an arbitrary set $T$, with no independence or measurability-of-the-supremum assumed between different $t$'s. Since $\sup_{t\in T}X_t(\omega)$ need not be measurable in $\omega$ for a general index set $T$, its expectation is understood — following the book's own convention, set once in Chapter 7 and reused throughout — through the process's finite-dimensional marginals:

$$\mathbb E\sup_{t\in T}X_t \;:=\; \sup_{T_0\subseteq T\text{ finite, nonempty}}\ \mathbb E\max_{t\in T_0}X_t.$$

Now fix a metric $d$ on $T$, making $(T,d)$ a metric space. The covering number $N(T,d,\varepsilon)$, for $\varepsilon>0$, is the smallest cardinality of a finite $\varepsilon$-net of $T$: a finite set $N\subseteq T$ such that every point of $T$ lies within distance $\varepsilon$ of some point of $N$ (or $N(T,d,\varepsilon):=\infty$ if no finite $\varepsilon$-net exists). The quantity $\log N(T,d,\varepsilon)$ is the metric entropy of $T$ at scale $\varepsilon$: it measures how large $T$ looks when resolved only down to scale $\varepsilon$.

A process $(X_t)_{t\in T}$ has sub-gaussian increments with parameter $K\ge0$ if

$$\|X_t-X_s\|_{\psi_2}\;\le\;K\,d(t,s)\qquad\text{for all }t,s\in T,$$

where $\|\cdot\|_{\psi_2}$ is the sub-gaussian (Orlicz) norm of Chapter 2: the smallest $u>0$ with $\mathbb E\exp((X_t-X_s)^2/u^2)\le2$. This says the increments of the process are controlled by the metric $d$ the way a Gaussian process's increments are controlled by its own canonical metric $d(t,s):=\|X_t-X_s\|_{L^2}$ — but without assuming $(X_t)_{t\in T}$ is Gaussian.

A class of Boolean functions $F$ on a set $\Omega$ shatters a subset $\Lambda\subseteq\Omega$ if every function $g:\Lambda\to\{0,1\}$ arises as the restriction to $\Lambda$ of some $f\in F$. The VC (Vapnik-Chervonenkis) dimension $\mathrm{vc}(F)$ is the largest cardinality of a subset of $\Omega$ shattered by $F$ (or $\infty$ if arbitrarily large finite subsets, or an infinite one, are shattered) — a purely combinatorial measure of how rich the class $F$ is.

Formalization targets

Goal (Theorem 8.1.3, Dudley's integral inequality)

$$\exists\,C>0:\quad\mathbb E\sup_{t\in T}X_t\;\le\;CK\int_0^\infty\sqrt{\log N(T,d,\varepsilon)}\;d\varepsilon$$

for every mean-zero random process $(X_t)_{t\in T}$ on a metric space $(T,d)$ with sub-gaussian increments parameter $K\ge0$, whenever the integral is finite. $C$ is the book's own unnamed absolute constant, never depending on $T$, $K$, or the process. This is the weakest stable form of the claim: no numeral is hard-coded for $C$, and the statement asks only for the shape of the bound, matching what the book actually proves.

Milestone (Theorem 8.3.16, Sauer-Shelah lemma)

$$|F|\;\le\;\sum_{k=0}^{d}\binom nk\;\le\;\left(\frac{en}{d}\right)^{d},\qquad d:=\mathrm{vc}(F),$$

for every class $F$ of Boolean functions on a finite $n$-point set $\Omega$. This is a purely combinatorial fact, with no probability involved, but it is the bridge (via the covering-number bound Theorem 8.3.18, outside this mission's scope) between the chapter's Dudley-inequality engine and its statistical-learning applications — a bound on how large a finite class of Boolean functions can be, in terms of a single combinatorial complexity parameter.

Significance

Dudley's inequality is, in the book's own words, "the main result" of the chaining chapter: it converts a purely geometric quantity — the metric entropy of an index set, computable in many cases from covering-number estimates already available for balls, ellipsoids, and other convex bodies — into a probabilistic control on the size of a random process indexed by that set. This is what lets later chapters (uniform laws of large numbers over function classes, the matrix deviation inequality, the Dvoretzky-Milman theorem on almost-spherical sections of convex bodies) bound suprema over infinite, even uncountable, index sets without ever performing a union bound. The bound is also known to be tight only up to a logarithmic factor in general — Sudakov's minoration inequality (Chapter 7) gives a matching lower bound for Gaussian processes, and the book's own Exercise 8.1.12 exhibits a set where the two bounds genuinely diverge — so the constant $C$ here cannot in general be sharpened away.

The Sauer-Shelah lemma is one of the two founding results of VC theory (together with the Glivenko-Cantelli-type uniform convergence it feeds into), independently discovered by Vapnik and Chervonenkis, Sauer, and Shelah in the early 1970s; it underlies the sample-complexity bounds of statistical learning theory (a hypothesis class with finite VC dimension is PAC-learnable) and, through Theorem 8.3.18, gives one of the two standard routes (the other being direct combinatorial counting) to bounding covering numbers of infinite function classes.

Both results are decades old and have long-established, standard proofs; no open mathematical question is being formalized. What this mission contributes is the machine-checked statement infrastructure — the goal and the Sauer-Shelah milestone, together with the definitions (CoveringNumber, ProcessESup, Shatters, VcDim) a faithful Lean rendering of either result needs — for a solver to close with a proof. No formalization of Dudley's inequality or the Sauer-Shelah lemma is known to exist on the platform prior to this mission.

Difficulty

The natural first idea for bounding $\mathbb E\sup_{t\in T}X_t$ is a single-scale $\varepsilon$-net argument: replace $T$ by a finite $\varepsilon$-net, bound the maximum over the (finite) net by a union bound using the sub-gaussian tail, and separately bound the error of replacing $T$ by the net using the Lipschitz-in-probability control the sub-gaussian-increments hypothesis gives. This works, but it only ever sees $T$ at one fixed resolution $\varepsilon$, and optimizing over $\varepsilon$ afterward gives a bound with an extra $\sqrt{\log(1/\varepsilon)}$-type loss that does not match Dudley's inequality. The actual difficulty is genuinely multi-scale: chaining builds a whole sequence of nets at dyadic scales $\varepsilon=2^{-k}$ simultaneously, connects each point of $T$ to its nearest net point at every scale to form a "chain" of successive approximations back to a single fixed basepoint, and telescopes the resulting sum of increments — turning $X_t$ itself into a sum of differences between successive links of the chain, each individually well controlled by the sub-gaussian hypothesis at its own scale. Passing from the resulting discrete sum over dyadic scales (Theorem 8.1.4) to the continuous integral of the goal is a further, separate technical step.

For the Sauer-Shelah lemma, the natural first idea — bound $|F|$ directly by counting — has no obvious purchase on an arbitrary class of Boolean functions. The actual argument goes through Pajor's lemma, which reduces bounding $|F|$ to counting the shattered subsets of $\Omega$ instead of the functions in $F$ themselves; only then does the cardinality bound $d=\mathrm{vc}(F)$ on shattered sets become directly usable, via a binomial-sum estimate.

Formalization scope

CoveringNumber T ε is ℕ∞-valued (ℕ∞ = WithTop ℕ), defined as the infimum, over the subtype of finite $\varepsilon$-nets of the whole type T (an instance of MetricSpace T), of their cardinality; the infimum of the empty family in this complete lattice is ⊤, reproducing "$N:= \infty$ if no finite net exists" with no case split. ProcessESup is EReal-valued, defined as the supremum over finite nonempty $T_0\subseteq T$ of the Bochner integral of the finite max — EReal, not ℝ, because a real-valued supremum would silently return the junk value $0$ if the family of marginal expectations were unbounded above. Shatters and VcDim are direct transcriptions of Definition 8.3.1, with VcDim valued in ℕ∞ via a supremum of Set.encard over the (always-nonempty, since $\varnothing$ is trivially shattered) subtype of shattered subsets. The goal's mean-zero hypothesis is stated as Integrable (X t) P ∧ ∫ X t = 0 rather than the bare equation, since a non-integrable variable's Bochner integral is 0 in Mathlib by convention regardless of its true mean — a bare-equation hypothesis would let a non-mean-zero, non-integrable process satisfy the theorem vacuously. Two further hypotheses make explicit what the book's own displayed statement treats as understood without spelling out: that $N(T,d,\varepsilon)$ is finite for every $\varepsilon>0$ (total boundedness of $T$), and that the resulting integrand is integrable on $(0,\infty)$ — both hold whenever $T$ is totally bounded, since the integrand vanishes once $\varepsilon\ge\mathrm{diam}(T)$, so neither hypothesis excludes any case the book's own proof does not also need. [Nonempty T] excludes the degenerate empty index set. The absolute constant $C$ is existentially quantified ahead of every type, instance, and hypothesis it is uniform over, and pinned to no numeral, matching "$C$ is an absolute constant" — a formalization hard-coding a specific numeral for $C$ would be invalidated by the next sharper constant in the literature and would not match what the book proves.

A trivializing formalization of the Sauer-Shelah lemma would fix $\mathrm{vc}(F)$ at a hard-coded small value, or drop the second (exponential) inequality in favor of the weaker first one; this mission's statement keeps both inequalities, with $d$ genuinely computed from VcDim, and handles the $d=0$ boundary (where the exponential bound's base involves a division by zero under Lean's x/0=0 convention) explicitly rather than excluding it, since x^0=1 still recovers the book's correct bound $|F|\le1$ there.

This mission covers Theorem 8.1.3 and Theorem 8.3.16 only; Theorem 8.3.18 (covering numbers via VC dimension) and Theorem 8.2.3 (the uniform law of large numbers, the chapter's direct application of Dudley's inequality) are left out, not approximated, for lack of the additional empirical- process measurability machinery — the class of Lipschitz functions of Eq. (8.22), measurability of the resulting empirical process — that a faithful statement of either would need beyond what this mission's items already provide. CoveringNumber and ProcessESup are reusable by any later chapter needing a metric space's covering numbers or a general random process's expected supremum (this book's own Chapters 7, 9, and 11 all use one or both); Shatters and VcDim are reusable by any later development of VC theory or statistical learning theory. Solvers' contributions are welcome on: the chaining argument itself (the mission's hardest open leaf, via the discrete dyadic form of Theorem 8.1.4), Pajor's lemma underlying Sauer-Shelah, and the binomial-sum estimate closing its second inequality.

Selected references

• R. M. Dudley, The sizes of compact subsets of Hilbert space and continuity of Gaussian processes, Journal of Functional Analysis 1 (1967), 290–330. https://doi.org/10.1016/0022-1236(67)90017-1
• N. Sauer, On the density of families of sets, Journal of Combinatorial Theory, Series A 13 (1972), 145–147. https://doi.org/10.1016/0097-3165(72)90019-2
• V. N. Vapnik, A. Ya. Chervonenkis, On the uniform convergence of relative frequencies of events to their probabilities, Theory of Probability & Its Applications 16 (1971), 264–280. https://doi.org/10.1137/1116025
• R. Vershynin, High-Dimensional Probability: An Introduction with Applications in Data Science, Cambridge University Press, 2018, Chapter 8. https://doi.org/10.1017/9781108231596

Motivation

Sums of independent random variables are well understood: Bernstein's inequality and its relatives give sharp, non-asymptotic tail bounds for $\sum_i a_i X_i$ whenever the $X_i$ are independent and light-tailed. Many quantities that arise in high-dimensional statistics and random matrix theory, however, are not linear but quadratic in an independent sample — the squared norm of a random vector after a linear transformation, a quadratic-form test statistic, the diagonal of a sample covariance matrix, or the number of edges cut by a random partition in a random graph. A quadratic form $X^\top A X = \sum_{i,j} A_{ij} X_i X_j$ is a sum with dependent terms: $X_iX_j$ and $X_iX_k$ share the factor $X_i$, so classical sum-of-independent-variables tools do not apply directly.

The Hanson-Wright inequality, first obtained by Hanson and Wright (1971) for sub-gaussian variables and later sharpened and popularized in this form by Rudelson and Vershynin (2013, "Hanson-Wright inequality and sub-gaussian concentration," Electronic Communications in Probability), closes this gap: it gives a concentration inequality for $X^\top A X$ around its mean with the same two-regime (sub-gaussian near the center, sub-exponential in the tail) shape as Bernstein's inequality for linear sums. It is now a standard tool wherever quadratic statistics of independent data are analyzed: covariance estimation, compressed sensing, randomized numerical linear algebra, and the analysis of random matrices more broadly draw on it routinely.

Setting

Fix a probability space and let $X = (X_1, \dots, X_n)$ be a random vector whose coordinates $X_1, \dots, X_n$ are independent, mean zero, and sub-gaussian: each $X_i$ has a finite sub-gaussian (Orlicz $\psi_2$) norm $\|X_i\|_{\psi_2}$, the smallest $t > 0$ with $\mathbb E \exp(X_i^2/t^2) \le 2$. Write $K = \max_i \|X_i\|_{\psi_2}$.

Let $A = (A_{ij})_{i,j=1}^n$ be an $n \times n$ real matrix, with no constraint on its diagonal, and form the quadratic form

$$X^\top A X = \sum_{i,j=1}^n A_{ij} X_i X_j.$$

Two matrix norms measure the size of $A$: the Frobenius norm $\|A\|_F = \bigl(\sum_{i,j} A_{ij}^2\bigr)^{1/2}$ (the Euclidean norm of $A$'s entries) and the operator (spectral) norm $\|A\| = \sup_{\|x\|_2=1} \|Ax\|_2$ (the largest singular value of $A$). Always $\|A\| \le \|A\|_F \le \sqrt{n}\,\|A\|$, so the two norms can differ by a factor as large as $\sqrt n$ — the gap between them is exactly what produces the inequality's two regimes below.

Formalization targets

Goal — Theorem 6.2.1 (Hanson-Wright inequality)

$$P\bigl\{\, |X^\top A X - \mathbb E\, X^\top A X| \ge t \,\bigr\} \;\le\; 2 \exp\!\left[-c \min\!\left(\frac{t^2}{K^4 \|A\|_F^2},\ \frac{t}{K^2 \|A\|}\right)\right] \qquad \text{for every } t \ge 0,$$

where $c > 0$ is an absolute constant, not depending on $n$, $X$, $A$, or $t$. Stating the constant only as "some absolute $c$" (rather than pinning it to a numeral) is deliberate: the book's own proof does not track a sharp value, and a goal that only asserts the shape of the bound survives any later improvement to $c$.

Significance

The result itself. Hanson-Wright turns a two-dimensional (in $i,j$) dependency structure into a one-dimensional concentration statement controlled by two scalar quantities, $\|A\|_F$ and $\|A\|$. This is what makes it usable: a practitioner bounding a quadratic statistic need only compute these two norms, not analyze the joint dependency structure of $\{X_iX_j\}$ directly. It specializes to Bernstein's inequality (Chapter 2 of this book) when $A$ is diagonal, and it underlies non-asymptotic guarantees for covariance estimation, the Johnson-Lindenstrauss lemma via a different route, and the concentration of Lipschitz functions of sub-gaussian vectors.

Formalizing it. The published proof of Hanson-Wright is not a single argument but a chain of four steps: a decoupling reduction (Section 6.1), a direct computation for Gaussian chaos (Lemma 6.2.2), a comparison lemma extending the Gaussian bound to general sub-gaussian vectors via a replacement trick (Lemma 6.2.3), and a final assembly that separates the diagonal part (handled by Bernstein's inequality) from the off-diagonal part (handled by decoupling and comparison). This mission formalizes the goal theorem's statement and the first, most reusable link in that chain — the decoupling machinery of Section 6.1, which reduces the analysis of the dependent chaos $X^\top A X$ to the independent-once-conditioned bilinear form $X^\top A X'$ — together with the chapter's separate contraction principle (Section 6.7), a general comparison tool for Rademacher-weighted sums used repeatedly in the book's later chaining chapters. The Gaussian MGF computation and the replacement-trick comparison lemma (Lemmas 6.2.2–6.2.3) are left as future milestones on top of this mission: they require Gaussian rotation invariance and the singular value decomposition of $A$, substantially more machinery than the milestones included here.

Difficulty

The obvious first idea — treat $X^\top A X = \sum_{i,j} A_{ij}X_iX_j$ as if it were a sum of independent terms and apply Bernstein's inequality termwise — fails immediately: the terms $A_{ij}X_iX_j$ for fixed $i$ are not independent across $j$, since they all share the factor $X_i$. Decoupling (Theorem 6.1.1) is the non-obvious fix: it replaces the off-diagonal chaos by a bilinear form $X^\top A X'$ in an independent copy $X'$, which genuinely does become a sum of independent terms once one of the two vectors is conditioned on. The price is a universal constant factor of $4$ and the restriction to diagonal-free matrices, which is exactly why the full Hanson-Wright proof must separate the diagonal contribution to $\mathbb E\,X^\top A X$ (handled directly by Bernstein's inequality, Chapter 2) before decoupling can be applied to what remains.

Formalization scope

Random variables and vectors are real-valued on an explicit probability space $(\Omega, \mathcal F, P)$. The sub-gaussian norm is HighDimProb.Concentration.subgaussianNorm, the Orlicz-$\psi_2$-norm definition already published for this series (01-concentration), reused here as a reference item rather than redefined. $K = \max_i \|X_i\|_{\psi_2}$ is written as a finite supremum over the coordinate index, ⨆ i, subgaussianNorm P (X i); because the index type is always a Fintype (Fin n), this supremum is well-defined and, at the degenerate index $n=0$, reduces to a true (if content-free) instance of the inequality rather than a vacuous or false one. The Frobenius and operator norms of $A$ are this mission's own frobeniusNorm and opNorm, stated directly from their defining formulas rather than through Mathlib's scoped matrix-norm typeclass instances, which are deliberately not global defaults (to avoid a diamond between the two norms) and so are unsuitable for a statement that needs both simultaneously. Every place the goal or a milestone integrates a quantity, that quantity is required Integrable, guarding against Mathlib's convention of returning 0 for the Bochner integral of a non-integrable function — without these hypotheses, a mean-zero or expectation hypothesis could hold vacuously, or a conclusion could hold trivially, for reasons having nothing to do with the book's mathematics.

The formalization deliberately does not restrict $A$'s diagonal in the goal theorem: doing so would collapse Hanson-Wright to a restatement of Bernstein's inequality for the special case of a diagonal matrix, discarding the chapter's actual content, which is handling the off-diagonal, genuinely quadratic dependence between coordinates. The diagonal-free restriction does appear, correctly, in the Decoupling theorem (6.1.1), whose proof needs it.

Reusable beyond this mission: frobeniusNorm and opNorm are needed by any future chapter using matrix norms (Chapter 4's random matrix norms, Chapter 9's matrix deviation inequality); the decoupling theorem and convex decoupling lemma are the standard entry point for any later formalization of chaos concentration; the contraction principle is reused throughout the book's chaining chapters (7 and 8). Welcome contributions include the Gaussian MGF and comparison lemmas (6.2.2–6.2.3) needed to complete a full proof of the goal theorem, and the two-sided version of Bernstein's inequality needed for the diagonal part of that proof.

Selected references

• R. Vershynin, High-Dimensional Probability: An Introduction with Applications in Data Science, Cambridge University Press, 2018. DOI: 10.1017/9781108231596.
• D. L. Hanson, F. T. Wright, "A bound on tail probabilities for quadratic forms in independent random variables," Annals of Mathematical Statistics 42 (1971), 1079–1083.
• M. Rudelson, R. Vershynin, "Hanson-Wright inequality and sub-gaussian concentration," Electronic Communications in Probability 18 (2013), no. 82, 1–9. https://arxiv.org/abs/1306.2872

Motivation and context

Large matrices are often reduced by retaining only some rows or columns. The reduced matrix is cheaper to store and apply, but its spectral norm can be distorted by a single unusually large row, column, or entry. A useful estimate must show both the average effect of sampling and the cost of these local concentrations. In §5.2 of Tropp's monograph, matrix Chernoff inequalities give numerical estimates for independently sampled columns and for independently sampled rows and columns. The latter is the capstone of this mission.

The general matrix Chernoff inequality, Theorem 5.1.1, is already an Open goal in tc's Chapter 5 mission on Prove2Me. The intrinsic-dimension Chernoff inequality, Theorem 7.2.1, is likewise Open in tc's Chapter 7 mission. This mission concerns the application statements that those chapter missions have not posed: the numerical bound (5.2.2), its three displayed intermediate estimates, and companion consequences in §§5.1.1, 5.2.1 and 7.2.2. The motivation for choosing this target is that it connects a general concentration result to a concrete random matrix model; its terms can be read directly from a deterministic input matrix.

Setting

Fix a complex matrix $B\in\mathbb C^{d\times n}$ with $d,n\ge1$. For each row $j$, let $\delta_j$ equal one when that row is retained and zero otherwise. For each column $k$, let $\xi_k$ have the analogous meaning. The whole family is independent, with $\delta_j\sim\operatorname{Bernoulli}(p/d)$ and $\xi_k\sim\operatorname{Bernoulli}(r/n)$, where $0\le p\le d$ and $0\le r\le n$. Thus $p$ and $r$ are the expected numbers of retained rows and columns. Put

$$P=\operatorname{diag}(\delta_1,\ldots,\delta_d),\qquad R=\operatorname{diag}(\xi_1,\ldots,\xi_n),\qquad Z=PBR.$$

The zeros remain in the matrix: $Z$ has the same dimensions as $B$, but entries in discarded rows or columns vanish. Write $b_{j:}$ for row $j$, $b_{:k}$ for column $k$, and $b_{jk}$ for an entry. Their vector norms are Euclidean, while $\|B\|$ and $\|Z\|$ denote the spectral operator norm. The source uses $\mathbb E$ for expectation over the selections. Its maxima range over all rows, columns or entries of $B$, so a highly uneven matrix is distinguished from one whose mass is spread evenly.

The first milestone is the streamlined upper matrix Chernoff expectation bound. For a finite independent family of positive semidefinite Hermitian matrices $X_k$ with eigenvalues at most $L$, let $Y=\sum_kX_k$ and $\mu_{\max}=\lambda_{\max}(\mathbb EY)$. The bound uses the book's numerical constant $1.72$. The subsequent milestones specialize this setting to the row and column selectors above. A separate companion states the lower bound with $0.63$ and the two weakened tail inequalities of §5.1.1.

Formalization targets

The goal is the exact displayed estimate (5.2.2), with the book's rounded constants:

$$\begin{aligned} \mathbb E\|PBR\|^2\le{}&3\frac pd\frac rn\|B\|^2 +2\frac{p\log d}{d}\max_k\|b_{:k}\|^2\\ &+2\frac{r\log n}{n}\max_j\|b_{j:}\|^2 +(\log d)(\log n)\max_{j,k}|b_{jk}|^2. \end{aligned}$$

Its milestone chain is (5.1.8), followed by (5.2.3), (5.2.4), and (5.2.5). These are the claims the text places immediately before the final estimate. They bound, respectively, the first column-selection stage, the norm of the row-selected matrix, and the largest column norm after row selection. The companion (5.2.1) gives both the top and least squared singular-value bounds for column selection alone. The companion from §7.2.2 replaces $\log d$ by $\log(2\operatorname{srank}(B))$, where $\operatorname{srank}(B)=\|B\|_F^2/\|B\|^2$ is stable rank.

Significance

The goal separates four sources of size in a sampled matrix. The first term scales the squared spectral norm by the expected row and column retention fractions. The remaining terms record concentrations in a column, a row, and an individual entry. Consequently, the estimate can be evaluated from $B$ and the two sampling rates without specifying a distribution for $B$ itself. It is an expectation bound, not a claim that every realized submatrix is close to $B$.

The formal development would make the application layer of Tropp's Chapter 5 available alongside tc's general, currently Open matrix Chernoff theorem. A proof of that general theorem does not itself identify the random projectors, their joint independence, or the deterministic quantities on the right side of (5.2.2). The stable-rank companion also links the Chapter 7 notion of intrinsic dimension to a concrete rectangular matrix statistic. Both the chapter's general theorems and this mission's statements remain Open on the platform; none is being presented as already machine proved. Tropp, §§5.2 and 7.2.2 gives the mathematical results.

Difficulty

Sampling both rows and columns introduces two interacting random projectors. Merely applying a one-stage column estimate to $B$ leaves the row selection uncontrolled; merely applying a one-stage row estimate leaves the sampled column norms random. The two families must be independent of each other as well as independent internally. That distinction is substantial: separate independence claims for the row family and the column family do not imply independence between $P$ and $R$. The final estimate must also retain the book's exact $3$ and $2$ coefficients rather than replacing the right side by an unspecified constant times a simpler expression.

Formalization scope

Matrices are Matrix (Fin d) (Fin n) ℂ, matching the book's general complex-matrix convention. The spectral norm, extreme eigenvalues, Borel measurable matrix space, and intrinsic dimension come from tc's published definitions, which this mission references. Random selectors are measurable real-valued functions on a probability space; their Bernoulli law specifies almost-sure values zero or one and the probability of one. The row and column selectors form one jointly independent family. Products with diagonal selector matrices define $PB$, $BR$ and $PBR$; a deterministic diagonal or an unconstrained zero-one function would not be the sampled model.

Every sampling rate is within its Bernoulli domain. Nonzero dimensions keep the maxima over finite, inhabited sets and give the logarithms their intended arguments. Lean's Bochner integral otherwise returns zero on a nonintegrable function; each expectation inequality records integrability in its conclusion, which follows from measurable, bounded Bernoulli selectors and Tropp's standing regularity convention in §2.2.1. The maxima of column, row, and entry magnitudes are real finite suprema, equal to maxima for these finite nonempty index sets. The intrinsic-dimension identity assumes $p>0$ and $B\ne0$ for scale invariance; the final stable-rank bound itself includes zero sampling and the zero matrix.

The definitions of Bernoulli selection, random diagonal projections, sampled matrices, finite norm maxima and stable rank are reusable outside this mission. Contributions can address their structural identities, integrability, the four milestones, the goal, or the independent column-selection companions. The book's constants $1.72$ and $0.63$ are retained as decimals; $3$ and $2$ are the rounded coefficients printed in (5.2.2), as the text explains on p. 66.

Selected references

• Joel A. Tropp, An Introduction to Matrix Concentration Inequalities, Foundations and Trends in Machine Learning 8(1–2), 2015, arXiv:1501.01571v1, §§2.1, 2.2, 5.1, 5.2, 7.2.
• Joel A. Tropp, An Introduction to Matrix Concentration Inequalities, Foundations and Trends in Machine Learning 8(1–2), 2015, DOI:10.1561/2200000048. The formalization cites the pinned arXiv version above for page and equation indices.

Motivation

Random matrices arise when observations, signs, or model coefficients fluctuate. Their spectral norm measures the largest amplification they can cause to a vector. A bound on that norm can turn entrywise randomness into a statement about the behavior of an entire linear operator. Tropp's chapter develops a common framework for Gaussian and Rademacher sums of fixed matrices, then applies it to Wigner matrices, rectangular Gaussian matrices, random signs, Toeplitz matrices, and a rounding scheme for a quadratic optimization relaxation (Tropp, Chapter 4). The applications show what the general estimate says when the matrices have concrete structure.

The chapter's general matrix-series theorem, Theorem 4.1.1, is already an open statement in tc's Prove2Me Chapter 4 mission. This mission targets the squared-norm comparison in §4.1.2 and the worked consequences that are not part of that existing goal. They are useful benchmarks because their dimensions, coefficients, and variance statistics are explicit. The existing platform status is a formalization status: Tropp proves these mathematical bounds in the cited book; the Chapter 4 theorem and the results proposed here do not yet have machine-checked proofs.

Setting

Fix positive integers $d_1,d_2$ and a finite family of complex $d_1\times d_2$ matrices $B_k$. Let $\gamma_k$ be independent standard Gaussian real variables, each with mean zero and variance one. Their matrix Gaussian series is $Z=\sum_k\gamma_k B_k$. For a complex matrix $M$, write $M^*$ for its conjugate transpose and $\|M\|$ for its Euclidean operator norm. The matrix variance statistic is

$$v(Z)=\max\{\|\mathbb E(ZZ^*)\|,\ \|\mathbb E(Z^*Z)\|\}.$$

For this Gaussian series, the two expected products equal the corresponding sums of coefficient products, so $v(Z)$ can be computed from the fixed $B_k$. This is the same statistic and norm as Theorem 4.1.1, equations (4.1.2)–(4.1.4) (Tropp, p. 42). Expectations are entrywise Bochner integrals for matrix-valued functions and ordinary integrals for real quantities.

A Rademacher variable takes values $-1$ and $1$ with equal probability. Some companion examples replace the Gaussian coefficients by independent Rademacher variables. For a Hermitian square matrix $Y$, $\lambda_{\min}(Y)$ denotes its smallest real eigenvalue. In the signed-entry example, the relevant matrix statistic is the larger of the greatest squared row norm and the greatest squared column norm of a fixed matrix.

Formalization targets

Squared norm of a Gaussian series

The goal is the pair of inequalities in (4.1.7), together with explicit integrability of the random squared norm:

$$v(Z)\le\mathbb E\|Z\|^2\le2v(Z)\bigl(1+\log(d_1+d_2)\bigr).$$

The constants and the logarithm are exactly those printed on pp. 43–44. The left side identifies a necessary scale; the right side gives a dimension-dependent upper estimate. The two milestones are the Jensen lower-bound display and the tail-integration display in §4.1.2. The latter is stated for any measurable matrix with the chapter's tail bound and positive variance, because the displayed integration step itself uses only those properties. It retains the exact split-point expression and the factor $2$.

Worked series

Companion statements record the minimum-eigenvalue consequences (4.6.4)–(4.6.5), the Wigner estimate (4.2.3), the rectangular Gaussian estimate (4.2.6), the randomly signed-entry estimate (4.3.3), the Toeplitz estimate (4.4.3), and the MaxQP rounding scale in §4.5. Their deterministic variance identities are part of the statements, rather than assumptions. In particular the Wigner bound keeps $\log(2d)$ and the Toeplitz bound keeps the exact coefficient sum $dI_d$ (Tropp, pp. 46–52).

Significance

The squared-norm comparison says that $v(Z)$ is the correct baseline for the size of a Gaussian series in mean square, while allowing a logarithmic penalty for matrix dimension. It also shows why an expectation estimate and a tail estimate answer different questions: the chapter discusses examples where a tail bound can be loose even when the expected norm is controlled (Tropp, pp. 44–45). The example statements make the abstract variance statistic usable for standard matrix ensembles and for signed matrices with prescribed entries.

A machine-checked development would connect the existing spectral and probability definitions to explicit examples, finite-sum variance identities, and integrability of Gaussian matrix norms. Those components can be reused in later sampling and numerical linear algebra missions. The general Chapter 4 matrix-series theorem remains an open platform target; these application statements are not claimed to be already formalized or proved there.

Difficulty

The central obstacle is moving from a tail estimate for $\|Z\|$ to a finite, sharp expectation for $\|Z\|^2$. The tail estimate has a dimension factor, while the desired bound has its logarithm; an unspecialized bound at a single threshold does not provide the stated constant. At the modeling level, matrix variance uses expectations of matrix products, so it must coincide with the finite sums of coefficient squares for each example. A definition that silently returns zero for a nonintegrable Bochner integral would make the upper inequality misleading. The Gaussian law, measurability, and moment behavior therefore matter to the statement as well as to eventual proofs.

Formalization scope

Matrices use Matrix (Fin d₁) (Fin d₂) ℂ with $d_1,d_2\ge1$ and tc's spectralNorm and rectSecondMoment definitions. The probability space is arbitrary; coefficient functions are measurable and independent, and their pushforward laws are exactly standard Gaussian or Rademacher as the page requires. Families may be empty. Real random coefficients scale complex matrices. Integrability of $\|Z\|^2$ is a conclusion of the goal and tail-integration milestone, while the minimum-eigenvalue expectation companion likewise asserts integrability. The tail-integration milestone assumes $v(Z)>0$ because its displayed denominator divides by $v(Z)$; the Gaussian goal itself includes $v(Z)=0$.

The Wigner coefficients are indexed by pairs $i<j$; rectangular and signed-entry coefficients by row-column pairs; Toeplitz coefficients by diagonals $-(d-1),\ldots,d-1$. The Toeplitz index in Lean is shifted to Fin (2*d-1). Row and column maxima range over the nonempty finite row and column sets. The signed-entry theorem permits complex fixed entries and uses $|b_{ij}|^2$; the printed example starts with real entries. MaxQP uses the positive scale $\alpha=1/\sqrt{2\log(d_1+d_2)}$ and tc's positive-semidefinite order for the two relaxation constraints. None of the coefficient-square identities is supplied as a hypothesis; each is a conclusion to establish.

Context. The existing TroppMatrixConcentration.ch4_matrix_series formalizes Theorem 4.1.1 and is Open. Its Gaussian/Rademacher variance identity and norm/tail bounds provide the chapter context. The present goal specializes the same Gaussian model and asks for the separate second-moment result (4.1.7). Reusable contributions include finite-matrix identities, Gaussian moment integrability, and tail-integration lemmas. A valid proof must preserve the actual standard laws and independence, and cannot redefine $v(Z)$ as the expected squared norm.

Selected references

• Joel A. Tropp, An Introduction to Matrix Concentration Inequalities, Foundations and Trends in Machine Learning 8(1–2), 2015, arXiv:1501.01571v1, Chapter 4, pp. 42–52.

Motivation

Many computations need the trace of a matrix $A$ that is never formed explicitly and can only be applied to vectors: the trace of a matrix function such as $\mathrm{trace}(A^{-1})$ or $\log\det A = \mathrm{trace}(\log A)$ in statistics and lattice QCD, the Frobenius norm $\|B\|_F^2 = \mathrm{trace}(B^TB)$ of an operator, or the number of triangles of a graph. The standard tool is Monte-Carlo estimation, introduced by M. F. Hutchinson (Hutchinson 1989): average $M$ quadratic forms $z_i^TAz_i$ over random sign vectors $z_i$. Each sample costs one matrix–vector product, uses one random bit per entry, and needs only additions and subtractions.

Before Avron and Toledo (2011), only the variance of such estimators had been analysed. A small variance does not say how many samples guarantee a given relative error with a given probability. Avron and Toledo gave the first bounds of this kind for several estimators. This mission covers the bound for Hutchinson's estimator, their Theorem 7.1.

Setting

A Rademacher random variable takes the values $+1$ and $-1$, each with probability $1/2$. Let $A \in \mathbb{R}^{n\times n}$ be symmetric positive semi-definite. Draw $M \ge 1$ independent random vectors $z_1, \ldots, z_M \in \mathbb{R}^n$ whose $Mn$ entries are independent Rademacher variables. Hutchinson's trace estimator is

$$H_M = \frac{1}{M}\sum_{i=1}^{M} z_i^TAz_i .$$

A single sample $z^TAz$ is an unbiased estimator of $\mathrm{trace}(A)$ (Lemma 2.1 of the paper, due to Hutchinson). For symmetric $A$ its variance is $2\bigl(\|A\|_F^2 - \sum_i A_{ii}^2\bigr)$, twice the squared Frobenius mass of $A$ off the diagonal.

Given $\epsilon > 0$ and $\delta \in (0,1)$, a random estimator $T$ is an $(\epsilon,\delta)$-approximator of $\mathrm{trace}(A)$ if

$$\Pr\bigl(|T - \mathrm{trace}(A)| \le \epsilon\,\mathrm{trace}(A)\bigr) \ge 1 - \delta .$$

The rank $\mathrm{rank}(A)$ is the number of nonzero eigenvalues $\lambda_j$ of $A$, counted with multiplicity.

Formalization targets

Goal: Theorem 7.1, sample bound for $H_M$

For every symmetric positive semi-definite $A$, every $0 < \epsilon \le 1/2$ and every $0 < \delta < 1$,

$$M \ge 6\epsilon^{-2}\ln\!\left(\frac{2\,\mathrm{rank}(A)}{\delta}\right) \quad\Longrightarrow\quad H_M \text{ is an } (\epsilon,\delta)\text{-approximator of } \mathrm{trace}(A).$$

The bound depends on $A$ only through its rank. It does not depend on the dimension $n$, on the condition number, or on how the trace is spread over the diagonal.

Milestones

1. Lemma 7.2 (Achlioptas 2001, Lemma 5). For a unit vector $\alpha \in \mathbb{R}^n$ and $S = \frac1M\sum_{i=1}^M(\alpha^Tz_i)^2$, for every $\epsilon > 0$,

$$\Pr(|S - 1| \ge \epsilon) \le 2\exp\!\left(-\frac{M}{2}\left(\frac{\epsilon^2}{2} - \frac{\epsilon^3}{3}\right)\right).$$

2. Per-direction bound (proof of Theorem 7.1, p. 8:11). Let $r \ge 1$ and $0 < \epsilon \le 1/2$. If $M \ge 6\epsilon^{-2}\ln(2r/\delta)$, then $\Pr(|S - 1| \ge \epsilon) \le \delta/r$.
3. From directions to the trace (proof of Theorem 7.1, p. 8:11). Write $A = U\Lambda U^T$ and $y_i = U^Tz_i$. If $|\frac1M\sum_i y_{ij}^2 - 1| \le \epsilon$ for every $j$ with $\lambda_j \ne 0$, then $|H_M - \mathrm{trace}(A)| \le \epsilon\,\mathrm{trace}(A)$. This step is deterministic.
4. Lemma 2.1 (Hutchinson). $\mathrm{E}(z^TAz) = \mathrm{trace}(A)$, and for symmetric $A$, $\mathrm{Var}(z^TAz) = 2(\|A\|_F^2 - \sum_i A_{ii}^2)$.

Significance

Theorem 7.1 gives a practitioner an explicit number of matrix–vector products after which Hutchinson's method is guaranteed to reach relative accuracy $\epsilon$ with confidence $1-\delta$. No bound was available before, although the method had been in wide use for two decades. The bound exceeds the one the same paper proves for Gaussian test vectors ($20\epsilon^{-2}\ln(2/\delta)$) by a $\ln(\mathrm{rank}(A))$ factor. The authors conjecture that this factor is not needed. Later work removed it: Roosta-Khorasani and Ascher (2015) proved a rank-free bound for the Rademacher case, and Cortinovis and Kressner (2022) extended sample bounds to indefinite matrices. Sample bounds of this type underlie variance-reduced estimators such as Hutch++ (Meyer et al. 2021).

The theorem is proved in the literature. To the platform's knowledge it has no machine-checked proof. Formalizing it produces a checked Rademacher concentration inequality for averages of squared linear forms (Achlioptas' lemma, which is also a core lemma of database-friendly Johnson–Lindenstrauss projections), and a checked spectral reduction from a quadratic-form estimator to its eigen-directions. Neither is currently in Mathlib.

Difficulty

The estimator is not an average of independent copies of a bounded variable with a small range: a single sample $z^TAz$ can deviate from $\mathrm{trace}(A)$ by an amount comparable to $\|A\|_F\sqrt n$. Chebyshev's inequality with the variance of Lemma 2.1 only gives a polynomial dependence on $1/\delta$. Hoeffding's inequality applied to $z^TAz$ directly gives a dependence on $n$ and on the size of the entries of $A$. Unlike the Gaussian case, the estimator cannot be written as a weighted sum of independent chi-squared variables, because rotating a Rademacher vector does not give another Rademacher vector. The central difficulty is Lemma 7.2: a tail bound for $(\alpha^Tz)^2$ that holds uniformly over every unit direction $\alpha$, including directions in which $\alpha^Tz$ is far from Gaussian (for $\alpha = e_1$ it is a constant).

Formalization scope

Matrices are Matrix (Fin n) (Fin n) ℝ, and positive semi-definiteness is Matrix.PosSemidef. The Rademacher law is $\tfrac12(\delta_1 + \delta_{-1})$ on ℝ. The sample space is Fin M → Fin n → ℝ with the product of $Mn$ copies of this law (hutchinsonSampleMeasure), so the law of the estimator is constructed, not assumed. $H_M(\omega) = (M:\mathbb{R})^{-1}\sum_i \omega_i\cdot(A\omega_i)$. Probabilities are Measure.real, and the $(\epsilon,\delta)$-approximator is Definition 4.1 verbatim. $\mathrm{rank}(A)$ is Matrix.rank. In Lemma 7.2 the i.i.d. copies $Q_i$ are the functions $(\alpha\cdot z_i)^2$ on this product space.

Corrections of printed statements, each recorded in the item's Formalization Note:

• Theorem 7.1 is printed without a range on $\epsilon$. Its proof needs $\frac{M}{2}(\frac{\epsilon^2}{2} - \frac{\epsilon^3}{3}) \ge \frac{M\epsilon^2}{6}$, which holds exactly when $\epsilon \le 1/2$. Without a range the printed statement is false: take $A = \frac1n\mathbf 1\mathbf 1^T$, $n = 10^4$, $\epsilon = 10$ and $\delta = 10^{-4}$. The condition admits $M = 1$, while $\Pr(H_1 > 11) \approx 9\cdot10^{-4} > \delta$. The goal and milestone 2 are therefore stated for $0 < \epsilon \le 1/2$.
• Lemma 2.1 is printed for an arbitrary $n\times n$ matrix. The variance formula fails for non-symmetric $A$: for $A = \begin{pmatrix}0&1\\0&0\end{pmatrix}$ the variance is $1$, not $2$. Symmetry is assumed for the variance part only.
• Proof of Theorem 7.1. The proof writes $\Lambda = UAU^T$ together with $y_i = U^Tz_i$. These fit together only for $A = U\Lambda U^T$, which is the convention of milestone 3.

For $A = 0$ the threshold involves $\ln 0$. Lean's Real.log 0 = 0 turns the condition into $M \ge 0$, which agrees with the paper's reading $\ln 0 = -\infty$. The conclusion then holds because $H_M = \mathrm{trace}(A) = 0$, so no extra hypothesis is added. A formalization that took the law of the samples as a hypothesis could make that hypothesis unsatisfiable and the theorem vacuous; the constructed product space rules this out.

A complete development needs:

• the product Rademacher measure and moment generating functions of Rademacher sums;
• a Chernoff bound for averages of i.i.d. bounded variables on a product space;
• the spectral theorem for real symmetric matrices, with rank equal to the number of nonzero eigenvalues;
• a finite union bound.

The Rademacher concentration results are reusable beyond this mission. Proofs of any milestone are welcome, as are alternative proofs of Lemma 7.2.

Selected references

• H. Avron and S. Toledo, Randomized algorithms for estimating the trace of an implicit symmetric positive semi-definite matrix, Journal of the ACM 58(2), Article 8, 2011. https://doi.org/10.1145/1944345.1944349
• M. F. Hutchinson, A stochastic estimator of the trace of the influence matrix for Laplacian smoothing splines, Communications in Statistics – Simulation and Computation 18(3), 1059–1076, 1989. https://doi.org/10.1080/03610918908812806
• D. Achlioptas, Database-friendly random projections, Proceedings of PODS 2001, 274–281. https://doi.org/10.1145/375551.375608
• F. Roosta-Khorasani and U. Ascher, Improved bounds on sample size for implicit matrix trace estimators, Foundations of Computational Mathematics 15, 1187–1212, 2015. https://doi.org/10.1007/s10208-014-9220-1
• A. Cortinovis and D. Kressner, On randomized trace estimates for indefinite matrices with an application to determinants, Foundations of Computational Mathematics 22, 875–903, 2022. https://doi.org/10.1007/s10208-021-09525-9
• R. A. Meyer, C. Musco, C. Musco and D. P. Woodruff, Hutch++: Optimal stochastic trace estimation, SOSA 2021, 142–155. https://doi.org/10.1137/1.9781611976496.16

Motivation

Many matrices in scientific computing are available only implicitly: one can multiply a vector by $A$, for instance by solving a linear system or applying a matrix function, but the entries of $A$ are never formed. Quantities such as the trace must then be estimated from a small number of matrix–vector products. Randomized trace estimators do exactly this: draw random vectors $z$, compute the quadratic forms $z^TAz$, and average.

Avron and Toledo (J. ACM 58(2), 2011) gave the first sample bounds of the form "$M$ samples suffice for relative error $\epsilon$ with probability $1-\delta$" for the standard estimators. Among their results is a case in which the estimate is not merely approximate but exact: when $A$ is a projection matrix, its trace equals its rank, an integer, and rounding a Gaussian trace estimate recovers that integer with high probability. Computing the rank of a projection arises, for example, when charge densities are computed in electronic structure calculations without diagonalization (Bekas, Kokiopoulou and Saad 2007).

Setting

Let $n \ge 1$ and $A \in \mathbb{R}^{n\times n}$. A projection matrix here is an orthogonal projection: $A$ is symmetric and $A^2 = A$. Equivalently, $A$ is symmetric and every eigenvalue of $A$ is $0$ or $1$; its rank $\mathrm{rank}(A)$ is the number of eigenvalues equal to $1$.

Fix a number of samples $M \ge 1$. Let $z_1,\ldots,z_M \in \mathbb{R}^n$ be random vectors whose $Mn$ entries are independent standard normal random variables. The Gaussian trace estimator (Definition 3.1 of the paper) is

$$G_M = \frac{1}{M}\sum_{i=1}^{M} z_i^T A z_i .$$

Its expectation is $\mathrm{trace}(A)$. For $x \in \mathbb{R}$, $\mathrm{round}(x)$ denotes the nearest integer to $x$. For $k \ge 1$, a random variable $X$ has the $\chi^2$ distribution with $k$ degrees of freedom, $X \sim \chi^2(k)$, if it has the law of $g_1^2 + \cdots + g_k^2$ for independent standard normal $g_1, \ldots, g_k$.

Formalization targets

Goal: Lemma 5.3 (p. 8:9)

For a projection matrix $A$, a failure probability $\delta > 0$, and every integer $M \ge 1$ with $M \ge 24\,\mathrm{rank}(A)\ln(2/\delta)$,

$$\Pr\bigl(\mathrm{round}(G_M) \ne \mathrm{rank}(A)\bigr) \le \delta .$$

The number of samples depends on the rank and on $\delta$ only; there is no accuracy parameter, and none on $n$.

Milestones, in the order the proof uses them

1. Law of $MG_M$. $MG_M \sim \chi^2(M\,\mathrm{rank}(A))$.
2. $\chi^2$ tail bound (cited by the paper from Li, Hastie and Church 2007). For $X \sim \chi^2(k)$, $k\ge1$ and $0 < \epsilon \le \tfrac12$,

$$\Pr(|X - k| \ge \epsilon k) \le 2\exp(-k\epsilon^2/6).$$

3. Tail of $G_M$. For $\mathrm{rank}(A) \ge 1$ and $0 < \epsilon \le \tfrac12$,

$$\Pr(|G_M - \mathrm{rank}(A)| \ge \mathrm{rank}(A)\epsilon) \le 2\exp(-M\,\mathrm{rank}(A)\epsilon^2/6).$$

4. Eq. (2). If moreover $M \ge 6\,\mathrm{rank}(A)^{-1}\epsilon^{-2}\ln(2/\delta)$, then $\Pr(|G_M - \mathrm{rank}(A)| \ge \mathrm{rank}(A)\epsilon) \le \delta$.
5. Trace of a projection. $\mathrm{trace}(A) = \mathrm{rank}(A)$.

Significance

The lemma turns a randomized estimator into an exact algorithm with a controlled failure probability: the rank of an implicitly given projection is obtained from $O(\mathrm{rank}(A)\log(1/\delta))$ matrix–vector products, independently of the dimension $n$. It is also an instance where the paper's general relative-error bound for the Gaussian estimator (Theorem 5.2, whose sample count grows like $\epsilon^{-2}$) is improved by exploiting the spectrum of $A$: an absolute error below $\tfrac12$ is a relative error $\epsilon = 1/(2\,\mathrm{rank}(A))$, for which the general bound would require a number of samples quadratic in the rank, whereas Lemma 5.3 needs only a linear number.

The result is proved in the paper. The mission produces a machine-checked version of it, together with two pieces of reusable substrate: the exact $\chi^2$ law of a Gaussian quadratic form in an orthogonal projection, and a two-sided $\chi^2$ tail bound with explicit constant $1/6$ on the range $0<\epsilon\le 1/2$. To the best of available knowledge, none of these statements is formalized in Mathlib; a platform mission on the Johnson–Lindenstrauss lemma states a $\chi^2$ concentration bound with a different exponent, $(\epsilon^2-\epsilon^3)/4$, on the open range $0<\epsilon<1/2$, which does not cover the value $\epsilon = 1/2$ needed here when $\mathrm{rank}(A)=1$.

Difficulty

The deterministic part is short; the probabilistic part is not. The step "$y = Uz$ has independent standard normal entries because $U$ is orthogonal" is the rotation invariance of the standard Gaussian measure on $\mathbb{R}^n$, and it must be combined with the independence of the $M$ samples to identify the law of a sum of $M\,\mathrm{rank}(A)$ squares; this is a statement about product measures and pushforwards, not about moments. The $\chi^2$ tail bound is quoted by the paper without proof. Its constant $1/6$ is not the constant of the usual textbook $\chi^2$ estimates, and the bound is false outside a restricted range of $\epsilon$ (see below), so a generic sub-exponential concentration inequality with unspecified constants does not deliver it. Finally, the rounding step requires matching Mathlib's round with the event $|G_M - \mathrm{rank}(A)| < \tfrac12$.

Formalization scope

All declarations live in the namespace TraceEstimation.ProjectionRank. Matrices are Matrix (Fin n) (Fin n) ℝ. The sample space of $G_M$ is Fin M → Fin n → ℝ with the product measure Measure.pi (fun _ => Measure.pi (fun _ => gaussianReal 0 1)), so the law of the estimator is constructed, not assumed. Probabilities are Measure.real of events. "Projection matrix" is A.IsHermitian ∧ A * A = A (orthogonal projection); a non-symmetric idempotent also has eigenvalues $0$ and $1$, but the paper's proof diagonalizes $A$ by a unitary matrix, which requires symmetry. round is Mathlib's round : ℝ → ℤ; the rank is Matrix.rank. The $\chi^2$ distribution is not defined: its role is played by the pushforward of a product of standard normals under $g \mapsto \sum_l g_l^2$.

Corrections of the printed text, both in the $\chi^2$ tail bound (milestone 2):

• The paper prints $\Pr(|X - k| \le \epsilon k) \le 2\exp(-k\epsilon^2/6)$. The inner $\le$ is a misprint for $\ge$; the next display applies the bound with $\ge$.
• The paper gives no range for $\epsilon$. The bound fails for $\epsilon = 1$ and large $k$, since $\Pr(X \ge 2k)$ decays like $e^{-k(1-\ln 2)/2}$, slower than $e^{-k/6}$. The mission states it for $0<\epsilon\le 1/2$, and milestones 3 and 4 inherit that range; the proof uses only $\epsilon = 1/(2\,\mathrm{rank}(A)) \le 1/2$.

The goal keeps the paper's hypotheses: $\delta > 0$ with no upper bound, and $\mathrm{rank}(A) = 0$ allowed (both are true cases). A statement about the event $|G_M - \mathrm{rank}(A)| < \tfrac12$, or Eq. (2) alone, is not the goal; the goal is about $\mathrm{round}(G_M)$. The sample measure is a probability measure, so the bound $\le \delta$ is not vacuous.

Contributions welcome beyond the milestones: rotation invariance of the standard Gaussian vector under orthogonal matrices, the moment generating function of $\chi^2(k)$, and the spectral fact $\mathrm{trace}(A) = \mathrm{rank}(A)$ for idempotents, each reusable well outside this mission.

Selected references

• H. Avron and S. Toledo, Randomized algorithms for estimating the trace of an implicit symmetric positive semi-definite matrix, Journal of the ACM 58(2), Article 8, 2011. https://doi.org/10.1145/1944345.1944349
• P. Li, T. Hastie and K. Church, Nonlinear estimators and tail bounds for dimension reduction in $l_1$ using Cauchy random projections, in Learning Theory (COLT 2007), Lecture Notes in Computer Science 4539, Springer, 514–529, 2007 (the version the paper cites); journal version in Journal of Machine Learning Research 8, 2497–2532, 2007, https://jmlr.org/papers/v8/li07b.html
• C. Bekas, E. Kokiopoulou and Y. Saad, An estimator for the diagonal of a matrix, Applied Numerical Mathematics 57(11–12), 1214–1229, 2007. https://doi.org/10.1016/j.apnum.2007.01.003
• M. F. Hutchinson, A stochastic estimator of the trace of the influence matrix for Laplacian smoothing splines, Communications in Statistics – Simulation and Computation 19(2), 433–450, 1990. https://doi.org/10.1080/03610919008812866