Motivation

When a system is controlled indefinitely and decisions are frequent — a router admitting packets, a queue accepting jobs, a machine being maintained — discounting future rewards is often unjustified, and what matters is the long-run average reward per period. Puterman's Chapter 8 (doi:10.1002/9780470316887) develops the theory of this criterion, and its central object is a single equation, the average reward optimality equation $0=\max_{a\in A_s}\{r(s,a)-g+\sum_jp(j\mid s,a)h(j)-h(s)\}$, whose unknowns are a scalar gain $g$ and a bias function $h$. For unichain models, in which every stationary policy generates a Markov chain with one recurrent class, this equation determines the optimal gain and an optimal stationary policy. The results go back to Howard (Dynamic Programming and Markov Processes, MIT Press, 1960) for the recurrent case and to Blackwell (Discrete dynamic programming, Annals of Mathematical Statistics 33, 1962, doi:10.1214/aoms/1177704593) and Derman for the general finite case; Puterman's Section 8.4 proves them through the discounted theory of mission II, by letting the discount factor tend to one.

Setting

The model is stationary (Assumption 8.0.1): a finite set $S$ of states, for each $s$ a finite nonempty set $A_s$ of actions, a reward $r(s,a)$ and transition probabilities $p(j\mid s,a)$, none depending on the decision epoch. A policy $\pi\in\Pi^{HR}$ may randomize and may depend on the whole history; the deterministic stationary policy $d^\infty$ applies the decision rule $d:S\to A$ at every epoch. Its transition matrix is $P_d(i,j)=p(j\mid i,d(i))$.

For a policy $\pi$, $v^\pi_{N+1}(s)=\mathbb E^\pi_s[\sum_{t=1}^Nr(X_t,Y_t)]$ is the expected reward over $N$ epochs. Since the limit of $N^{-1}v^\pi_{N+1}(s)$ need not exist (Example 8.1.1), the chapter works with the lim sup and lim inf average rewards $g^\pi_+(s)$ and $g^\pi_-(s)$, and with $g^*_\pm(s)=\sup_{\pi}g^\pi_\pm(s)$. A policy $\pi^*$ is average optimal when $g^{\pi^*}_-(s)\ge g^\pi_+(s)$ for all $s$ and $\pi$, the strongest of the three criteria of Section 8.1.2.

The optimality residual is $B(g,h)(s)=\max_{a\in A_s}\{r(s,a)-g+\sum_jp(j\mid s,a)h(j)-h(s)\}$, and the optimality equation is $B(g,h)=0$. A decision rule is $h$-improving when it attains $\max_{a\in A_s}\{r(s,a)+\sum_jp(j\mid s,a)h(j)\}$ at every state. A transition matrix is unichain when it consists of a single recurrent class plus a possibly empty set of transient states, and the MDP is unichain when $P_d$ is unichain for every deterministic decision rule.

Formalization targets

Goal — Theorem 8.4.5 (printed p. 361)

For a finite unichain model: (a) some deterministic stationary policy is average optimal; (b) the optimality equation $B(g^*,h^*)=0$ has a solution, and (d) its scalar satisfies $g^*_+(s)=g^*_-(s)=g^*$ for every $s$; (c) for every solution, every $h^*$-improving decision rule gives an average optimal stationary policy.

Theorem 8.4.1 (printed p. 356)

If $B(g,h)\le 0$ then $g\ge g^*_+$; if $B(g,h)\ge 0$ then $g\le\sup_{d}g^{d^\infty}_-\le g^*_-$; if $B(g,h)=0$ then $g^*_+=g^*_-=g$.

Theorem 8.4.3 (printed p. 358)

In a finite unichain model $B(g,h)=0$ has a solution, and every solution has the same $g$.

Theorem 8.4.4 (printed p. 361)

If $B(g^*,h^*)=0$ and $d^*$ is $h^*$-improving, then $(d^*)^\infty$ is average optimal.

Significance

Theorem 8.4.1(c) is what the source calls "one of the most important results for average reward models": a solution of the optimality equation with constant $g$ pins down the optimal gain under every criterion at once, so that in finite unichain models the three optimality criteria of Section 8.1.2 coincide. Theorem 8.4.3 guarantees such a solution exists, and Theorem 8.4.4 reads an optimal policy off it. Together, Theorem 8.4.5 reduces the infinite-horizon average reward problem over all history-dependent randomized policies to a finite system of equations in $(g,h)$, which is what policy iteration, value iteration and linear programming solve in Sections 8.5 to 8.8.

The results are classical and proved. Formalizing them fixes the chain-structure hypothesis in a checkable form and pins down which criterion "average optimal" means, two places where the literature is loose. The platform's MarkovDecisionProcesses series has the finite-horizon (mission I) and discounted (mission II) models; this mission adds the undiscounted stationary model, the gains, and the unichain classification, on which Chapter 9's multichain optimality equations and Chapter 10's sensitive discount optimality can be built.

Difficulty

The obvious argument for Theorem 8.4.3 is to take the discounted optimal value $v^*_\lambda$ of mission II and let $\lambda\uparrow 1$. It fails as stated because $v^*_\lambda$ blows up like $(1-\lambda)^{-1}$; what converges is the Laurent expansion $v^{d^\infty}_\lambda=(1-\lambda)^{-1}ge+h+o(1)$ of the value of a fixed stationary policy, Corollary 8.2.4, and that expansion needs the limiting matrix $P_d^*$ and the deviation matrix $H_{P_d}$ of a unichain chain. So the proof must first develop the Markov chain theory of Section 8.2 and Appendix A, choose a subsequence of discount factors along which one policy is discount optimal (possible because $D^{MD}$ is finite), and only then pass to the limit in the discounted optimality equation.

Theorem 8.4.1 looks elementary and hides the analytic step: iterating $ge\ge r_d+(P_d-I)h$ along an arbitrary history-dependent policy and dividing by $N$ requires the telescoping term $N^{-1}(P^\pi_N-I)h$ to vanish, which uses boundedness of $h$, and requires the reduction from history-dependent randomized to Markov randomized policies (Theorem 8.1.2). For Theorem 8.4.4 the step is Corollary 8.2.7, that $r_d-ge+(P_d-I)h=0$ forces the gain of $d^\infty$ to be $g$, which is the multiplication by $P_d^*$ that annihilates $(P_d-I)$.

The traps are in the definitions. Recurrence and the unichain property must be stated so that the source's Example 8.4.3 comes out as the book says — the policy using $a_{1,1}$ has the absorbing state $s_2$ as its single recurrent class — and the optimality residual must use the lim sup / lim inf gains, since a definition through a limit that need not exist would be a junk value on the policies of Example 8.1.1.

Formalization scope

State and action spaces are Fintypes and admissible actions are nonempty Finsets, as in missions I and II; the stationary model is a new structure because mission II's DiscountedMDP bundles a discount factor, and carries the same data otherwise. Policies are history-dependent and randomized, so "average optimal" has its full strength; a stationary policy is the deterministic one built from a decision rule. Expected total reward is defined by the policy evaluation recursion, as in the earlier missions, rather than through a measure on trajectories.

Gains are Filter.limsup and Filter.liminf of $N^{-1}v^\pi_{N+1}(s)$ on $\mathbb R$; these are the source's because the sequence is bounded by $\max|r|$, and the suprema $g^*_\pm$ over the nonempty family of policies are genuine real suprema for the same reason. The residual $B(g,h)$ is a Finset.sup' over the admissible actions. Recurrence is "every state reachable from $i$ reaches $i$" and unichain is "any two recurrent states communicate", the definitions of Appendix A for finite chains, applied to $P_d$ for every admissible deterministic decision rule.

Restrictions relative to the printed text, all noted in the items: Theorem 8.4.1 is stated for finite $S$ where the source says countable, since the chapter's standing assumption and the model are finite; the gain $g^{d^\infty}$ of a stationary policy in (8.4.5) is written as its lim inf gain, which equals it; and the chain $ge=g^*=g^*_+=g^*_-$ of (8.4.6) is stated through $g^*_+$ and $g^*_-$, since $g^*$ presupposes existing limits. Nothing is trivialized: the existential in Theorem 8.4.5(a) has to produce a decision rule, and $B(g,h)=0$ with a junk maximum is impossible since every $A_s$ is nonempty. Welcome contributions beyond the milestones: Theorem 8.1.2 (reduction to Markov policies), Corollary 8.2.7 (the gain of a stationary policy from the evaluation equations), and the equivalence of the three optimality criteria in finite models.

Selected references

• Martin L. Puterman, Markov Decision Processes: Discrete Stochastic Dynamic Programming, Wiley, 1994, Chapter 8. doi:10.1002/9780470316887
• Ronald A. Howard, Dynamic Programming and Markov Processes, MIT Press, 1960.
• David Blackwell, Discrete dynamic programming, Annals of Mathematical Statistics 33 (1962). doi:10.1214/aoms/1177704593
• Cyrus Derman, Finite State Markovian Decision Processes, Academic Press, 1970.
• Paul J. Schweitzer and Awi Federgruen, The functional equations of undiscounted Markov renewal programming, Mathematics of Operations Research 3 (1978). doi:10.1287/moor.3.4.308

Motivation

A decision maker who does not want to be judged by expected loss alone needs a functional that weighs the bad outcomes more heavily than the good ones and still behaves well under optimization: monotone, convex, unchanged by shifting the loss, and scaling with it. Kusuoka's theorem (mission II) says every such law-invariant functional on an atomless space is a supremum of mixtures of Average Values-at-Risk. The mixtures themselves are the distortion risk functionals of Denneberg (Distorted probabilities and insurance premiums, Methods of Operations Research 63, 1990) and Acerbi (Spectral measures of risk: A coherent representation of subjective risk aversion, Journal of Banking and Finance 26, 2002, doi:10.1016/S0378-4266(02)00281-9), also called spectral risk measures: integrals of the quantile function against a nondecreasing density $\sigma$. They are the risk functionals used throughout Pflug and Pichler's book (doi:10.1007/978-3-319-08843-3) in the multistage objectives of Chapters 5 and 6, and the three representations proved in Sections 3.2 to 3.4 — as a supremum over densities $Z$, as a maximum over uniform variables, and as an infimum over convex-conjugate constraints — are what make them computable inside an optimization model. This mission formalizes those representations.

Setting

Fix a probability space $(\Omega,\mathcal F,P)$. A random variable $Y:\Omega\to\mathbb R$ is a loss, and the functionals act on $L^\infty$, the bounded measurable ones. The Value-at-Risk at level $u\in(0,1]$ is the lower quantile $\mathsf{V@R}_u(Y)=\inf\{y:P(Y\le y)\ge u\}$ and the Average Value-at-Risk at level $\alpha\in[0,1)$ is $\mathsf{AV@R}_\alpha(Y)=\frac{1}{1-\alpha}\int_\alpha^1\mathsf{V@R}_u(Y)\,du$, extended to $\alpha=1$ by the essential supremum (mission II). A risk functional satisfies the four axioms (M), (C), (T), (H) of Definition 3.2.

A distortion function is a nonnegative, nondecreasing $\sigma:[0,1)\to[0,\infty)$ with $\int_0^1\sigma=1$, and the distortion risk functional with density $\sigma$ is

$$\mathcal R_\sigma(Y)\;=\;\int_0^1\sigma(u)\,\mathsf{V@R}_u(Y)\,du .$$

The Average Value-at-Risk is the case $\sigma_\alpha=(1-\alpha)^{-1}\mathbf 1_{[\alpha,1)}$. A random variable $Z$ is dominated by $\sigma$, $Z\preccurlyeq\sigma$, when $Z\in L^1$, $\mathbb E(Z)=1$ and $\mathsf{AV@R}_\alpha(Z)\le\frac{1}{1-\alpha}\int_\alpha^1\sigma$ for every $\alpha\in[0,1)$. A variable $U$ is uniformly distributed when $P(U\le u)=u$ on $[0,1]$. For $h:\mathbb R\to\mathbb R$ the conjugate is $h^*(s)=\sup_y(s\,y-h(y))\in(-\infty,+\infty]$.

Formalization targets

Goal — Theorem 3.16 (printed pp. 105–106)

$$\mathcal R_\sigma(Y)\;=\;\sup\bigl\{\mathbb E(Y\cdot Z)\ :\ Z\preccurlyeq\sigma\bigr\} \qquad(Y\in L^\infty).$$

The distortion functional is a risk functional (printed p. 99)

$\mathcal R_\sigma$ satisfies (M), (C), (T) and (H).

Representation (3.11) (printed p. 100)

There is a probability measure $\mu$ on $[0,1]$ with $\mathcal R_\sigma(Y)=\int_0^1\mathsf{AV@R}_\alpha(Y)\,\mu(d\alpha)$ for all $Y\in L^\infty$.

Corollary 3.18 (printed pp. 107–108)

For $\alpha<1$, $\mathsf{AV@R}_\alpha(Y)=\sup\{\mathbb E(YZ):\mathbb EZ=1,\ \mathsf{AV@R}_p(Z)\le\frac1{1-\alpha}\ \forall p\in[\alpha,1]\} =\sup\{\mathbb E(YZ):\mathbb EZ=1,\ 0\le Z\le\frac1{1-\alpha}\}$; for $\alpha=1$ the latter with the upper bound dropped.

Corollary 3.19 (printed p. 108)

On an atomless space, $\mathcal R_\sigma(Y)=\max\{\mathbb E(Y\cdot\sigma(U)):U\text{ uniform}\}$, the maximum attained.

Theorem 3.22 (printed p. 111)

$\mathcal R_\sigma(Y)=\inf\{\mathbb E(h(Y)):\int_0^1h^*(\sigma(u))\,du\le 0\}$ over measurable $h:\mathbb R\to\mathbb R$.

Significance

Theorem 3.16 identifies the conjugate of a distortion functional: $\mathcal R_\sigma^*(Z)$ is $0$ when $Z\preccurlyeq\sigma$ and $+\infty$ otherwise, so $\mathcal R_\sigma$ is the support function of an explicit convex set of densities. Everything else follows from that. Corollary 3.18 recovers the classical dual of the Average Value-at-Risk and shows that its constraints can be relaxed to the levels in $[\alpha,1]$; Corollary 3.19 turns the supremum into a maximum over couplings and identifies the maximizer, the co-monotone one, which is how distortion functionals are evaluated on scenario trees; Corollary 3.21, not stated here, combines the theorem with Kusuoka's representation to give the dual of every version independent risk functional. Theorem 3.22 goes the other way and writes $\mathcal R_\sigma$ as an infimum, which is what converts a minimax problem — minimize a supremum over $Z$ — into a plain minimization over the decision and an auxiliary function $h$, the form in which the book solves risk-averse programs in Chapters 5 and 6.

Everything here is proved in the source and in the cited papers; the platform has the finite-scenario Artzner–Delbaen representation of coherent risk measures but nothing on distortion functionals, and Mathlib has no risk-measure material and no quantile-based representation theory. The definitions of this mission sit on those of mission II and are reused as such.

Difficulty

The obvious route to Theorem 3.16 is Fenchel–Moreau duality: $\mathcal R_\sigma$ is convex and lower semicontinuous on $L^\infty$, so it equals its biconjugate, and the task is to compute $\mathcal R_\sigma^*(Z)=\sup_Y\,\mathbb E(YZ)-\mathcal R_\sigma(Y)$. The step where the naive computation stalls is bounding $\mathbb E(YZ)$ by a quantity that depends only on the distributions of $Y$ and $Z$: this is the rearrangement (Chebyshev, Hardy–Littlewood) inequality $\mathbb E(YZ)\le\int_0^1G_Y^{-1}(u)G_Z^{-1}(u)\,du$, with equality for co-monotone couplings. Given it, $\mathcal R_\sigma^*(Z)=\sup_Y\int_0^1G_Y^{-1}(G_Z^{-1}-\sigma)$, and testing with indicator-type $Y$ shows the supremum is $0$ exactly when the upper-tail averages of $Z$ are dominated by those of $\sigma$, which is the constraint $Z\preccurlyeq\sigma$. Neither the rearrangement inequality nor the identification of the conjugate is available in Mathlib.

For the corollaries the difficulties are concrete. Corollary 3.18 needs $Z\ge 0$ to be deduced from the constraints, which the source does by contradiction through the value $p=P(Z<0)$. Corollary 3.19 needs the co-monotone coupling to exist, which is where the atomless hypothesis enters, and needs $\sigma(U)$ to be feasible for (3.15), which uses $G_{\sigma(U)}^{-1}=\sigma$. Theorem 3.22 needs, besides the inequality $\mathbb E(h(Y))\ge\mathcal R_\sigma(Y)$ from Fenchel–Young, an admissible $h$ that nearly attains it; the source builds it from $\sigma$ and $G_Y^{-1}$ via Corollary 3.23.

Formalization scope

Random variables are functions $\Omega\to\mathbb R$ on an arbitrary measurable space with a probability measure, exactly as in mission II, whose MemLinfty, valueAtRisk, averageValueAtRisk, IsRiskFunctional, Atomless and IsKusuokaMeasure are imported and not redefined. The distortion density is a function on $\mathbb R$ constrained on $[0,1)$; its integrability on $(0,1)$ is part of the definition, and both $\int_0^1\sigma$ and $\mathcal R_\sigma$ integrate over the open interval, which excludes the level $0$ where the quantile formula is not meaningful and changes no value.

Every supremum and infimum is over a subtype of functions satisfying the constraints, and the prose of each item records why the family is nonempty and bounded, so that no real supremum takes its junk value. Pointwise constraints on $Z$ are almost sure. The constraint of (3.15) is stated for $\alpha\in[0,1)$: at $\alpha=1$ the source's expression is $0/0$ and means the limit, and the constraint there is implied by the others; (3.19) at $\alpha=1$ is stated separately with the upper bound dropped. The conjugate takes values in the extended reals, and admissibility in Theorem 3.22 requires $h^*(\sigma(u))$ finite almost everywhere, integrable, with integral at most $0$, and $h(Y)$ integrable — the conditions under which the expectations in (3.24) are the integrals the source means.

Two hypotheses are added beyond the printed statements and are flagged in the items: atomless for Corollary 3.19, because otherwise no uniform variable need exist and the maximum would be over an empty set; and integrability of $h(Y)$ in Theorem 3.22, without which Lean's integral of a non-integrable $h(Y)$ is $0$. Theorem 3.16 itself carries no atomless hypothesis, and the mission notes check the identity on two-point spaces by hand. A trivializing formalization — a constraint set that is empty or unbounded, so that the supremum is $0$ — is excluded by the constant density $Z\equiv 1$, feasible for every $\sigma$, and by $Z\ge 0$.

Welcome contributions beyond the milestones: Corollary 3.15 (the representation through the distribution function for $Y\ge 0$), Corollary 3.21 (the dual of a version independent risk functional through its Kusuoka set), Corollary 3.23, and the explicit mixing measure (3.12).

Selected references

• Georg Ch. Pflug and Alois Pichler, Multistage Stochastic Optimization, Springer, 2014, Sections 3.2–3.4. doi:10.1007/978-3-319-08843-3
• Carlo Acerbi, Spectral measures of risk: A coherent representation of subjective risk aversion, Journal of Banking and Finance 26 (2002). doi:10.1016/S0378-4266(02)00281-9
• Shigeo Kusuoka, On law invariant coherent risk measures, Advances in Mathematical Economics 3 (2001). doi:10.1007/978-4-431-67891-5_4
• Alois Pichler, The natural Banach space for version independent risk measures, Insurance: Mathematics and Economics 53 (2013). doi:10.1016/j.insmatheco.2013.07.005
• Georg Ch. Pflug and Werner Römisch, Modeling, Measuring and Managing Risk, World Scientific, 2007. doi:10.1142/6478

Motivation

Perturbative quantum field theory computes correlation functions of a field by expanding around a Gaussian (free) theory. Every term of that expansion is a Feynman diagram, and the rule that turns a diagram into a number is Wick's theorem: the expectation of a product of Gaussian field modes is the sum, over all ways of pairing the modes up, of the product of the two-point functions of the pairs. The same identity is known in probability and statistics as the Isserlis theorem (L. Isserlis, 1918) and is the standard tool for computing moments of Gaussian vectors; in random-matrix theory the counting of pairings it produces is the origin of the Catalan-number asymptotics of Wigner's semicircle law.

The uploaded source is the Wikipedia article Feynman diagram, which states Wick's theorem for the free scalar field and then, in the section Higher Gaussian moments — completing Wick's theorem, verifies the one-variable case by direct Gaussian integration. This mission formalizes that content: the combinatorics of pairings, the one-dimensional Gaussian moment formulas, and the multivariate identity itself.

Setting

Fix $d, n \in \mathbb{N}$ and work on $\mathbb{R}^d$ with coordinates $x_1,\dots,x_d$. Let $\mu$ be a centered Gaussian measure on $\mathbb{R}^d$: a Gaussian probability measure all of whose coordinate means vanish, $\int x_i \, d\mu(x) = 0$ for every $i$. Its covariance (in the physics reading, the propagator) is

$$G_{ij} \;=\; \int x_i x_j \, d\mu(x).$$

A pairing of the labels $\{0,1,\dots,2n-1\}$ is a partition of these $2n$ labels into $n$ unordered pairs; equivalently, a permutation $\sigma$ of the labels with $\sigma\circ\sigma = \mathrm{id}$ and $\sigma(i)\neq i$ for all $i$ (a fixed-point-free involution). The set of pairings is written $P_n$. For a weight $W_{ab}$ indexed by labels, the Wick sum is

$$\mathrm{Wick}(W) \;=\; \sum_{\sigma \in P_n} \ \prod_{\substack{i \,:\, i < \sigma(i)}} W_{i\,\sigma(i)},$$

the inner product ranging over the $n$ pairs of $\sigma$, each counted once through its smaller element.

In the article's field-theory notation the labels are momenta $k_1,\dots,k_{2n}$, the coordinates are the field modes $\phi(k_j)$, and the two-point function carries the momentum-conserving delta function, $\langle \phi(k)\phi(k')\rangle = \delta(k-k')/k^2$. This mission works with the finite-dimensional Gaussian vector rather than the field, so the delta functions are absorbed into the covariance matrix $G$.

Formalization targets

Goal — Wick's theorem (Isserlis' theorem)

For a centered Gaussian measure $\mu$ on $\mathbb{R}^d$ and any labels $k_1,\dots,k_{2n} \in \{1,\dots,d\}$,

$$\int \prod_{j=1}^{2n} x_{k_j} \, d\mu(x) \;=\; \sum_{\sigma\in P_n} \ \prod_{i < \sigma(i)} \left( \int x_{k_i} x_{k_{\sigma(i)}} \, d\mu(x) \right).$$

No hypothesis is imposed on the covariance: it may be singular and the labels $k_j$ may repeat, which is exactly the situation the article's "completing Wick's theorem" section addresses.

Supporting targets (milestones)

1. $\int_{\mathbb{R}} e^{-a x^{2}/2}\,dx = \sqrt{2\pi/a}$ for $a>0$.
2. $\int_{\mathbb{R}} x^{2n} e^{-a x^{2}/2}\,dx = \dfrac{(2n-1)!!}{a^{n}}\sqrt{2\pi/a}$ for $a>0$.
3. $\int x^{2n}\,d\mathcal{N}(0,v) = (2n-1)!!\, v^{n}$ for a real Gaussian law of variance $v \ge 0$.
4. $\#P_n = (2n-1)!!$.
5. Correlation functions of odd order vanish: $\int \prod_{j=1}^{2n+1} x_{k_j}\,d\mu = 0$.
6. The four-point function: $\langle x_{k_1}x_{k_2}x_{k_3}x_{k_4}\rangle$ equals the sum of the three products $G_{k_1k_2}G_{k_3k_4} + G_{k_1k_3}G_{k_2k_4} + G_{k_1k_4}G_{k_2k_3}$.

Targets 1–3 are the article's displayed Gaussian integrals, target 4 is its pairing count, targets 5–6 are the two explicit consequences it records for the field correlators.

Significance

Wick's theorem is the computational content of every Feynman-diagram expansion: once it is available, a perturbative term is a finite sum over diagrams, and the symmetry factors of diagrams are bookkeeping on the pairing set $P_n$. On the probabilistic side it gives all moments of a Gaussian vector in closed form, which is the entry point to Gaussian chaos expansions, Wiener–Itô integrals, and moment methods for random matrices.

Mathlib (revision 0df444a) has real Gaussian measures gaussianReal, the general class IsGaussian of Gaussian measures on a topological vector space, the Gaussian integral $\int e^{-bx^2} = \sqrt{\pi/b}$, and the double factorial Nat.doubleFactorial, but no higher-moment formula for Gaussian measures and no Isserlis/Wick statement. The mission therefore produces new library-level content, not a re-derivation of existing formal results; the result itself has been classical since 1918 (Isserlis) and 1950 (Wick).

Difficulty

The obvious route — expand the characteristic function $\exp(-\tfrac12 t^{\mathsf T} G t)$ and differentiate $2n$ times at $t=0$ — requires differentiating under an integral sign $2n$ times and identifying the resulting combinatorial sum with a sum over pairings; both steps are where the formal work lies. Integrability is not automatic from the statement and has to be established (Gaussian measures have moments of all orders, but the product $\prod_j x_{k_j}$ must be shown integrable before any manipulation). The naive attempt to reduce to the independent case by diagonalizing $G$ meets a second difficulty: the change of variables must be tracked through the pairing sum, and $G$ may be singular, so no invertible whitening transform exists in general. The one-variable case (milestone 3) is not a special case to be waved through either: it is the statement the article singles out, because a naive "each mode pairs with a distinct partner" argument fails when all labels coincide.

Formalization scope

The ambient space is EuclideanSpace ℝ (Fin d); measures are Mathlib Measures and Gaussianity is the Mathlib class IsGaussian, which is defined by every continuous linear functional pushing forward to a real Gaussian law. Centering is stated as an explicit hypothesis on the coordinate means, so the measure is not assumed standard and the covariance is unconstrained (in particular degenerate covariances, and repeated labels $k_i = k_j$, are included). Integrals are Bochner integrals, which return $0$ for non-integrable functions; the statements are nonetheless non-vacuous because Gaussian measures integrate all polynomials.

Pairings are formalized as fixed-point-free involutions of Fin (2 * n) and the pair product ranges over $\{i : i < \sigma(i)\}$, so each pair contributes once. The case $n = 0$ is included: the empty product is $1$, the unique pairing of the empty label set is the identity, and both sides of the goal equal $1$. The double factorial is Mathlib's Nat.doubleFactorial, evaluated at 2 * n - 1 in truncated natural subtraction, so the $n = 0$ value is $0!! = 1$.

Contributions welcome: the Gaussian moment lemmas (milestones 1–3) as standalone Mathlib-style results, the pairing count (milestone 4) as pure combinatorics independent of the analysis, and any reduction of the goal to the independent-coordinate case.

Selected references

• L. Isserlis, On a formula for the product-moment coefficient of any order of a normal frequency distribution in any number of variables, Biometrika 12 (1918), 134–139. DOI: 10.1093/biomet/12.1-2.134
• G. C. Wick, The evaluation of the collision matrix, Physical Review 80 (1950), 268–272. DOI: 10.1103/PhysRev.80.268
• Feynman diagram, Wikipedia. https://en.wikipedia.org/wiki/Feynman_diagram

Motivation

Reinforcement learning formalizes a scenario supervised learning cannot: an agent that actively interacts with an environment, choosing actions that change both the state it observes next and the reward it receives, rather than passively receiving an i.i.d. labeled sample. Every practical treatment of this scenario — from classical dynamic programming to modern deep reinforcement learning — is built on the Markov decision process (MDP), a model in which the effect of an action depends only on the current state, not on the full history that led to it. Two questions define the theory this mission covers: given a fixed way of acting (a policy), what value does it obtain, and how is that value actually computed rather than merely characterized as the solution of a fixed-point equation? Mohri, Rostamizadeh and Talwalkar's chapter 17 answers both for the stationary, infinite-horizon discounted case, and this mission targets its two central results: that a fixed policy's value is not just characterized but uniquely determined by a linear system with an explicit closed-form solution (Theorem 17.10), and that the optimal value function — obtained instead by choosing the best action at every state — can be computed by an iterative algorithm guaranteed to converge regardless of where it starts (Theorem 17.11).

Setting

A (finite) Markov decision process consists of a finite set of states $S$, a finite set of actions $A$, a transition kernel $P[s'\mid s,a]$ giving the distribution over the next state $s'$ after taking action $a$ at state $s$, and an expected reward $\mathbb E[r(s,a)]$ for that transition. A (stationary) policy $\pi:S\to\Delta(A)$ assigns each state a distribution over actions — possibly, but not necessarily, a point mass on a single action. Fixing $\pi$ turns the MDP into an ordinary Markov chain on $S$: at each step the agent is at some state $s$, draws $a\sim\pi(s)$, receives (expected) reward $\mathbb E[r(s,a)]$, and moves to a state drawn from $P[\cdot\mid s,a]$. For a discount factor $\gamma\in[0,1)$, the value of $\pi$ at $s$ is the expected discounted sum of future rewards starting from $s$,

$$V_\pi(s) = \mathbb E_{a_t\sim\pi(s_t)}\Big[\sum_{t=0}^{+\infty}\gamma^t r(s_t,a_t) \;\Big|\; s_0=s\Big],$$

and the state-action value function $Q_\pi(s,a)$ is the analogous quantity for taking $a$ at $s$ and then following $\pi$. Marginalizing the raw kernel and reward over the mixed action $\pi(s)$ gives the induced transition matrix $P_{s,s'}=P[s'\mid s,\pi(s)]=\sum_a \pi(s)(a) P[s'\mid s,a]$ and induced reward vector $R_s=\mathbb E[r(s,\pi(s))]=\sum_a\pi(s)(a)\,\mathbb E[r(s,a)]$ — the objects that turn $\pi$'s value into a genuinely linear-algebraic quantity. A policy $\pi^*$ is optimal if $V_{\pi^*}(s)\ge V_\pi(s)$ for every policy $\pi$ and every state $s$; write $V^*$ for its value function.

Formalization targets

Theorem 17.10 (goal). For a finite MDP and a fixed policy $\pi$, the matrix $I-\gamma P$ (with $P$ the policy-induced transition matrix) is invertible, and $\pi$'s value function is the unique solution of the Bellman equations, given in closed form by

$$V_\pi = (I-\gamma P)^{-1} R.$$

Proposition 17.9 (milestone). The value function itself satisfies the linear system that Theorem 17.10 solves:

$$\forall s\in S,\quad V_\pi(s) = \mathbb E_{a\sim\pi(s)}[r(s,a)] + \gamma\sum_{s'} P[s'\mid s,\pi(s)]\,V_\pi(s').$$

Theorem 17.7 (milestone). A policy $\pi$ is optimal if and only if it places probability only on $Q_\pi$-maximizing actions: for every $(s,a)$ with $\pi(s)(a)>0$, $a\in \operatorname{argmax}_{a'} Q_\pi(s,a')$.

Theorem 17.11 (milestone). The Bellman optimality operator $\Phi$, $[\Phi(V)](s)=\max_{a} \{\mathbb E[r(s,a)]+\gamma\sum_{s'}P[s'\mid s,a]V(s')\}$, is a $\gamma$-contraction for $\lVert\cdot\rVert_\infty$; consequently, for any starting vector $V_0$, the value-iteration sequence $V_{n+1}=\Phi(V_n)$ converges to a fixed point of $\Phi$.

Significance

Theorem 17.10 is what makes policy evaluation on a finite MDP an exact, finite computation rather than an infinite limit: instead of summing an infinite discounted series or solving an implicit fixed-point equation numerically, a single $|S|\times|S|$ matrix inversion gives the policy's value at every state simultaneously. It is also the base case every planning algorithm in the chapter builds on: policy iteration alternates optimizing a policy with exactly this evaluation step. Theorem 17.11 gives the complementary guarantee for the harder problem of finding the optimal value function directly, without fixing a policy first: value iteration converges from any starting point, with a convergence rate ($O(\log(1/\epsilon))$ iterations for $\epsilon$-accuracy) that follows from the same contraction argument. Together, the two results are the mathematical content behind why dynamic-programming planning for finite MDPs is tractable at all — the discount factor $\gamma<1$, not any structural assumption on rewards or transitions, is what buys both the uniqueness in Theorem 17.10 and the convergence in Theorem 17.11. Formalizing them requires reproducing this linear-algebraic and metric content precisely, not just asserting the conclusions: an invertibility claim asserted without the operator-norm argument, or a convergence claim without the contraction property, would state something true by fiat rather than the book's actual result. No faithful prior art exists on the platform for this exact model (see Formalization scope).

Difficulty

The obvious shortcut for Theorem 17.10 is to assert $I-\gamma P$ is invertible without proof — true, but not what the book does, and not informative about why it holds. The genuine content is that $P$, being row-stochastic (every row of $P$ sums to exactly $1$, since $\pi(s)$ and $P[\cdot\mid s,a]$ are both proper distributions), has operator norm $\lVert P\rVert_\infty=1$ exactly, so $\lVert\gamma P\rVert_\infty=\gamma<1$ strictly; this rules out $1$ as an eigenvalue of $\gamma P$, which is exactly what invertibility of $I-\gamma P$ requires. The same $\gamma<1$ fact, applied differently, drives Theorem 17.11: showing $\Phi$ is $\gamma$-Lipschitz requires bounding $\Phi(V)(s)-\Phi(U)(s)$ by comparing the maximizing action for $V$ against the same action's value under $U$ (not $U$'s own maximizer), since the two suprema need not be attained at the same action — a step easy to state incorrectly as a direct comparison of two maxima. Both theorems fail if $\gamma=1$ is allowed: the discounted setting's central asset, a strict contraction, disappears exactly at that boundary.

Formalization scope

States and actions are modeled as finite types (Fintype S, Fintype A); the raw kernel and reward P : S → A → S → ℝ, Er : S → A → ℝ are unconstrained functions, with IsTransitionKernel asserting the required distribution property explicitly rather than assuming it silently. A policy is π : S → A → ℝ with IsPolicy π asserting π s is a distribution over A for every s — deliberately not π : S → A or a PMF-valued function, since Theorem 17.7's own quantifier ("for any pair (s,a) with π(s)(a) > 0") requires treating π(s) as a genuine mixture. PolicyValue is defined as the actual infinite discounted expectation (via an explicit state-occupation-distribution recursion), not as the Bellman fixed point — so that Proposition 17.9 (the value function satisfies the linear system) and Theorem 17.10 (that system has a unique, invertible-matrix solution) are both non-vacuous claims about the same object, rather than one being definitionally true of the other. The trivializing formalization this rules out is asserting IsUnit (1 - γ • P) as a bare hypothesis, or defining V_π as (1-γP)⁻¹R and calling the resulting identity a theorem; both would erase the mission's actual content. Two platform modules model related MDPs (BertsekasSSPModel, a stochastic-shortest-path model with a termination-probability deficit rather than exact row-stochasticity, and FoundationsRL.RLBasics, a finite-horizon episodic model indexed by layer) — neither specializes exactly to this chapter's stationary, always-continuing, infinite-horizon discounted convention, so every definition here is drafted fresh rather than imported. This chunk covers §17.2–17.4.2 (the MDP model, policy value, Bellman's equations, value and policy iteration); §17.4.3 (the linear-programming formulation) and §17.5 (stochastic-approximation learning algorithms — TD(0), Q-learning, SARSA) are out of scope, since they require a stochastic-approximation convergence substrate this mission does not build.

Selected references

• Mohri, M., Rostamizadeh, A., and Talwalkar, A. Foundations of Machine Learning, 2nd ed., chapter 17. MIT Press, 2018.
• Bellman, R. Dynamic Programming. Princeton University Press, 1957.
• Puterman, M. L. Markov Decision Processes: Discrete Stochastic Dynamic Programming. Wiley, 1994.

Comparing random outcomes without comparing numbers

Operations research and applied probability are full of questions of the shape "is X riskier, larger, or later than Y?" when X and Y are random: two queueing disciplines' waiting times, two inventory policies' shortfall, two insurance portfolios' claim size. Comparing their means answers a weaker question than the modeler usually wants, since two distributions can share a mean while one is uniformly the "worse" one to face. Shaked and Shanthikumar's Stochastic Orders (Springer, 2007) is the standard reference for the family of orders built to make "X is stochastically no larger than Y" precise, and this mission formalizes the foundational member of that family: the usual stochastic order, together with its central structural theorem.

The usual stochastic order

Let $X$ be a real-valued random variable on a probability space $(\Omega,\mu)$, and let $Y$ be a real-valued random variable on a (possibly different) probability space $(\Omega',\nu)$. $X$ is smaller than $Y$ in the usual stochastic order, written $X \le_{st} Y$, if

$$P\{X > x\} \le P\{Y > x\} \quad \text{for every } x \in (-\infty,\infty).$$

Equivalently, $P\{X \le x\} \ge P\{Y \le x\}$ for every $x$: $Y$'s cumulative distribution function never exceeds $X$'s. Two features of the definition matter for everything that follows. First, it is a comparison of distributions, not of jointly defined variables — $X$ and $Y$ need not live on the same probability space, and the order depends only on their laws. Second, the inequality must hold for every threshold $x$ simultaneously; it is not a comparison of two summary statistics such as means or medians, and a distribution can fail $X \le_{st} Y$ even while $E[X] \le E[Y]$.

A closely related order compares residual lifetimes rather than raw tail probabilities. $X$ is smaller than $Y$ in the hazard rate order, written $X \le_{hr} Y$, if the survival functions $\bar F(x) = P\{X>x\}$, $\bar G(x) = P\{Y>x\}$ satisfy $\bar F(x)\bar G(y) \ge \bar F(y)\bar G(x)$ for all $x \le y$ — the general, density-free form of the condition that for absolutely continuous distributions reduces to a pointwise comparison of hazard rates $r(t) = f(t)/\bar F(t)$. It is one of several stronger orders (alongside the likelihood ratio order) that the book shows all imply $\le_{st}$, and it is the mission's second definition.

Formalization targets

Goal: the coupling characterization (Theorem 1.A.1)

$$X \le_{st} Y \iff \exists\,(\Omega'',\rho),\ \hat X,\hat Y : \Omega'' \to \mathbb{R} \text{ with } \hat X =_{st} X,\ \hat Y =_{st} Y,\ \text{and } P\{\hat X \le \hat Y\} = 1.$$

This is the order's Strassen-type coupling theorem: $X \le_{st} Y$ holds exactly when copies of $X$ and $Y$ can be realized on one common probability space so that, almost surely, the copy of $X$ never exceeds the copy of $Y$. It is the weakest possible statement with this shape (existence of some coupling on some space, not a canonical one), and it is the archetype every later chapter's own coupling theorem — for the hazard rate order, the convex order (via a martingale coupling), the multivariate orders — restates for its own comparison.

Supporting milestones

• Theorem 1.A.2, the same characterization restated via a common random source: $X \le_{st} Y$ iff there is a random variable $Z$ and functions $\psi_1 \le \psi_2$ with $X =_{st} \psi_1(Z)$ and $Y =_{st} \psi_2(Z)$.
• Theorem 1.A.3(a), closure under monotone transformations: $X \le_{st} Y$ and $g$ increasing gives $g(X) \le_{st} g(Y)$ (reversed for $g$ decreasing).
• Theorem 1.A.3(b), closure under increasing functions of independent families and, as its corollary, under convolution: independent $X_i \le_{st} Y_i$ and $\psi$ increasing gives $\psi(X_1,\dots,X_m) \le_{st} \psi(Y_1,\dots,Y_m)$, in particular $\sum_i X_i \le_{st} \sum_i Y_i$.
• Theorem 1.B.1, the first link of the book's implication chain: $X \le_{hr} Y \implies X \le_{st} Y$.

Significance

The usual stochastic order is the weakest and most widely used of the book's orders, and the coupling theorem is the tool that makes it tractable: comparing tail probabilities for every threshold is an infinite family of inequalities, while exhibiting one coupling on one space settles the comparison in a single stroke. This is the technique behind, among many applications in the book and the wider literature, comparing waiting times across queueing disciplines, bounding the effect of parameter uncertainty on a risk measure, and proving monotonicity results for Markov chains by coupling their sample paths. The closure properties (Theorem 1.A.3) are what let such comparisons survive composition with a monotone payoff function or a sum of independent components — exactly the operations a queueing or inventory model performs on its primitives.

Formalizing this theorem establishes, in one place, both the statement and the exact shape of its existential witness (a genuine ∃ over a probability space and two identically distributed random variables, not merely two abstractly-known-to-exist objects), which is the piece every later chapter's own coupling theorem in this series needs to get right independently. No proof is attempted in this mission; the book's own three-line construction ($\hat X = F^{-1}(U)$, $\hat Y = G^{-1}(U)$ for a uniform $U$) is a natural target for a future proof-bearing pass.

Difficulty

The definitional trap is the "for all $x$" in $\le_{st}$: it is tempting to shortcut it into a comparison of one summary statistic (a mean, an essential supremum, a specific quantile), which is uniformly weaker and would make Theorem 1.A.1 either false or a different, easier theorem. The coupling theorem itself has a second trap in its existential structure: the common probability space $(\Omega'',\rho)$ must be genuinely quantified — fixed neither to $\Omega$ nor $\Omega'$ in advance — and "$\hat X =_{st} X$" is equality of the two random variables' laws, not almost-sure equality of $\hat X$ and $X$ themselves (they typically do not share a domain before the construction). Theorem 1.A.3(b)'s general statement, for an arbitrary increasing $\psi:\mathbb{R}^m \to \mathbb{R}$ on independent families, is easy to understate as only its convolution corollary; the book gives both, and the general statement is the one with the actual content (the sum is recovered by taking $\psi$ to be the coordinate sum, itself increasing).

Formalization scope

Random variables are represented as measurable functions into $\mathbb{R}$ from arbitrary measurable spaces, and "equality in law" is Mathlib's native ProbabilityTheory.IdentDistrib, which is defined for two (possibly different) probability spaces and requires almost-everywhere measurability rather than fixing an ambient space in advance — matching the book's own treatment of $X \le_{st} Y$ as a comparison between distributions. UsualOrder μ ν X Y and HazardRateOrder μ ν X Y both take $X:\Omega\to\mathbb{R}$ on $(\Omega,\mu)$ and $Y:\Omega'\to\mathbb{R}$ on a separate $(\Omega',\nu)$, so no theorem here can be trivialized by silently forcing $X$ and $Y$ onto one shared space. HazardRateOrder is stated via the survival-function cross-product inequality (Eq. 1.B.4), the form that needs no absolute-continuity hypothesis, rather than via a derivative-based hazard rate, so it applies to the same general random variables as UsualOrder without an extra regularity hypothesis the book itself does not require at this level of generality. Independence in Theorem 1.A.3(b) is Mathlib's ProbabilityTheory.iIndepFun; the random variable $Z$ of Theorem 1.A.2 is left valued in an arbitrary measurable space, matching the book's unrestricted statement, rather than narrowed to $\mathbb{R}$-valued for convenience. A trivializing formalization is ruled out explicitly: $\le_{st}$ is not formalized as a comparison of means, medians, or any single moment, and the coupling theorems' ∃ is a genuine existential over a probability space, not a hypothesis smuggling in the witness's existence.

This mission draws on no platform prior art (repeated searches for "stochastic order", "coupling", "Strassen" and "identically distributed" turned up nothing on-topic as of 2026-09-18); it is the first mission in a nine-chapter series covering the whole book, and its own definitions are not imported by later chapters — each restates locally what it needs, per the series' own convention that draft items cannot import another draft's definitions. Reusable beyond this mission: the UsualOrder/HazardRateOrder shape (parametrizing over two separate probability spaces) is the pattern every later chapter's own order definition follows.

Selected references

• M. Shaked and J. G. Shanthikumar, Stochastic Orders, Springer Series in Statistics, Springer, 2007. https://doi.org/10.1007/978-0-387-34675-5

Motivation

Every convergence-rate result in the preceding chapters is an upper bound: some specific estimator (the Lasso, PCA, kernel ridge regression) achieves a given error rate. A natural and much harder question is the complementary one: is that rate actually the best any procedure could achieve, no matter its computational cost? Answering this requires a theory of lower bounds that holds simultaneously for every conceivable estimator — a fundamentally different kind of argument from constructing and analyzing one particular algorithm. Wainwright's High-Dimensional Statistics: A Non-Asymptotic Viewpoint (Cambridge University Press, 2019), Chapter 15, develops this theory, unifying classical techniques (Le Cam, Assouad, Fano) under one reduction: converting continuous estimation into discrete hypothesis testing.

Setting

Let $\mathcal X$ be a sample space and $\mathcal P$ a class of probability distributions on $\mathcal X$. A functional $\theta: \mathcal P \to \Omega$ assigns each distribution a parameter of interest. An estimator is a measurable map $\hat\theta: \mathcal X \to \Omega$. Fix a semi-metric $\rho:\Omega\times\Omega\to[0,\infty)$ — symmetric, triangle-inequality-satisfying, $\rho(\theta,\theta)=0$, but possibly $\rho(\theta,\theta')=0$ for $\theta\ne\theta'$ — and an increasing $\Phi:[0,\infty)\to[0,\infty)$. The minimax risk is

$$\mathfrak M\bigl(\theta(\mathcal P); \Phi\circ\rho\bigr) \;:=\; \inf_{\hat\theta} \sup_{P\in \mathcal P} \mathbb E_P\bigl[\Phi\bigl(\rho(\hat\theta, \theta(P))\bigr)\bigr],$$

the smallest worst-case expected loss achievable by any measurable estimator (Eq. (15.2)).

Given a 2δ-separated set $\{\theta_1,\dots,\theta_M\} \subseteq \theta(\mathcal P)$ (every pair satisfies $\rho(\theta_j,\theta_k)\ge 2\delta$) with representative distributions $P_{\theta_1},\dots,P_{\theta_M}$, Wainwright constructs a testing problem: sample $J$ uniformly from $[M]$, then $Z\sim P_{\theta_J}$; write $Q$ for the resulting joint law of $(J,Z)$. A test function $\psi:\mathcal X\to[M]$ attempts to recover $J$ from $Z$; its error probability is $Q[\psi(Z)\ne J]$.

Formalization targets

Goal (Proposition 15.1, "From estimation to testing")

For any increasing $\Phi$ and any 2δ-separated set with its induced joint testing measure $Q$,

$$\mathfrak M\bigl(\theta(\mathcal P); \Phi\circ\rho\bigr) \;\ge\; \Phi(\delta) \inf_{\psi} Q\bigl[\psi(Z)\ne J\bigr].$$

This mission formalizes Proposition 15.1 alone (see Formalization scope below for why, and what a follow-up mission would add to reach Fano's method proper, Proposition 15.12).

Significance

Proposition 15.1 is the single reduction every subsequent technique in the chapter specializes: Le Cam's two-point method (Lemma 15.9, $M=2$, bounding the testing error via total variation distance), Fano's method (Proposition 15.12, bounding it via mutual information $I(Z;J)$ and Fano's inequality), and Assouad's method (a different, hypercube-based packing). Formalizing it in full generality — general $\Phi$, general semi-metric, general $M$-ary packing set — gives a single reusable lemma that a future mission proving any of these specific bounds can build on directly, rather than re-deriving the reduction each time.

The theorem is already proved in the source; this mission's contribution is a machine-checked formal statement (and, eventually, proof) of the reduction, in a form composing with any future formalization of the chapter's testing-error bounds (total variation, Fano, or otherwise).

Difficulty

The proof combines two ingredients that must each be kept in their sharpest form: Markov's inequality applied to $\Phi(\rho(\hat\theta,\theta))$ (which only needs $\Phi$ increasing, not convex or any specific shape — a premature specialization to $\Phi(t)=t^2$ would silently prove a weaker, less reusable statement), and the reduction of any estimator to a test via nearest- packing-point assignment (Eq. (15.4)), which uses the triangle inequality on $\rho$ in a specific direction (bounding $\rho(\theta_k,\hat\theta)$ from below via $\rho(\theta_j,\theta_k)$ and $\rho(\theta_j,\hat\theta)$) to show that a small estimation error forces the induced test to be correct. Getting the direction and strictness of every inequality right — non-strict separation, but strict distance in the "test is correct" event — is where a naive restatement goes wrong.

Formalization scope

$\mathcal X$, $\Omega$ are arbitrary measurable spaces; the distribution class $\mathcal P$ is realized as an indexed family measure : Idx → Measure 𝒳 rather than a bare set of measures, composing directly with θ : Idx → Ω. The semi-metric $\rho$ is a bare function with explicit nonnegativity/reflexivity/symmetry/triangle-inequality hypotheses, matching the book's own footnote definition, rather than Mathlib's PseudoMetricSpace typeclass (kept self-contained, no extra instance machinery). The minimax risk is valued in ENNReal via the lower Lebesgue integral ∫⁻, not the Bochner integral, specifically to avoid the non-integrable-loss junk value 0 that a Bochner-integral formalization would silently introduce — a trivializing formalization would use ∫ (Bochner) here, letting a non-integrable loss vanish and making the inequality easier to satisfy than the book's actual claim; this mission does not do that. The joint testing measure $Q$ is characterized by its slice-measure equations directly on the product space $[M]\times\mathcal X$, avoiding Mathlib's general conditional/disintegration machinery while remaining exactly equivalent to "$J$ uniform, $Z\mid J=j\sim P_{\theta_j}$".

Disclosed major scope decision. BRIEF.md recommended Proposition 15.12 (the Fano bound itself, Φ(δ)(1 - (I(Z;J)+\log 2)/\log M)) as this mission's goal. That statement requires, in addition to everything above, a formalized notion of mutual information $I(Z;J)$ between a finite-valued and a general (possibly continuous) random variable, and its use of a Fano-type inequality (Eq. (15.31), itself deferred by the book to "Section 15.4" and not fully quoted in the brief). Building a faithful mutual-information formalization general enough for this setting (finite $J$, arbitrary measurable $Z$) — matching Mathlib's or this repo's existing, narrower information-theoretic developments (SourceCoding.*, built for a different, channel-coding purpose per BRIEF.md's own prior-art note) or building one from scratch — is substantially more than this session's remaining budget after building and self-reviewing 08-pca and 07-sparse-linear. This mission instead formalizes Proposition 15.1, the foundational reduction Proposition 15.12 itself specializes (via a particular bound on $\inf_\psi Q[\psi(Z)\ne J]$), so that a follow-up mission can add the mutual-information/Fano step on top of HighDimStat.Minimax.estimation_to_testing without redoing this reduction. This mission's name, fixed from missions/README.md, still names "Fano's Method" as the series slot this chunk occupies; its actual content is the reduction step every method in that family (including Fano's) shares — recorded explicitly here and in STATUS.md, not left implicit.

Selected references

• Wainwright, M. J. High-Dimensional Statistics: A Non-Asymptotic Viewpoint. Cambridge University Press, 2019. Chapter 15. DOI: 10.1017/9781108627771.
• Le Cam, L. "Convergence of estimates under dimensionality restrictions." Annals of Statistics, 1(1), 1973, 38–53.
• Fano, R. M. Transmission of Information: A Statistical Theory of Communications. MIT Press, 1961.
• Yu, B. "Assouad, Fano, and Le Cam." In Festschrift for Lucien Le Cam, Springer, 1997, 423–435.

Ordering how strongly two variables move together

The first eight chapters of Shaked and Shanthikumar's Stochastic Orders (Springer, 2007) compare single random variables — location, variability, convexity. Chapter IX turns to a different question: given two coordinates of a random vector, how strongly do they tend to move together? "Large values of $X_1$ go with large values of $X_2$" is a qualitative property (positive dependence); comparing how much two bivariate distributions exhibit it, holding their individual marginals fixed, is what this chapter's positive dependence orders make precise. This mission formalizes the two central ones — the PQD order and the supermodular order — and the supermodular order's closure properties, the chapter's own capstone result.

The PQD and supermodular orders

Let $X=(X_1,\dots,X_n)$ and $Y=(Y_1,\dots,Y_n)$ be two random vectors, with joint survival function $\bar F(x)=P\{X_1>x_1,\dots,X_n>x_n\}$ and joint distribution function $F(x)=P\{X_1\le x_1,\dots,X_n\le x_n\}$ (similarly $\bar G, G$ for $Y$). $X$ is smaller than $Y$ in the PQD order, written $X\le_{PQD}Y$, if

$$\bar F(x)\le\bar G(x)\ \text{ and }\ F(x)\le G(x) \quad\text{for every }x\in\mathbb{R}^n.$$

(It follows, though it need not be assumed, that $X$ and $Y$ then share the same univariate marginals — the two inequalities together pin the marginals down, unlike a single one alone.) This is Lehmann's positive-quadrant-dependence order in its general multivariate form: $Y$'s coordinates cluster together, in both the upper and lower "quadrants," at least as strongly as $X$'s.

A function $\varphi:\mathbb{R}^n\to\mathbb{R}$ is supermodular if $\varphi(x)+\varphi(y)\le\varphi(x\wedge y)+\varphi(x\vee y)$ for the coordinatewise meet/join — the function-level notion the platform's Topkis/Supermodularity series already formalizes for games and lattices, reused here directly. $X$ is smaller than $Y$ in the supermodular order, written $X\le_{sm}Y$, if $E[\varphi(X)]\le E[\varphi(Y)]$ for every supermodular $\varphi$ for which the two expectations exist. Because the indicator of any upper or lower orthant is itself supermodular, $X\le_{sm}Y\implies X\le_{PQD}Y$ (Eq. 9.A.17): the supermodular order is strictly finer than the PQD order, and for $n=2$ the two coincide.

Formalization targets

Goal: closure properties of the supermodular order (Theorem 9.A.9(a),(c))

$$(X_1,\dots,X_n)\le_{sm}(Y_1,\dots,Y_n) \implies (g_1(X_1),\dots,g_n(X_n))\le_{sm}(g_1(Y_1),\dots,g_n(Y_n))$$

whenever every $g_i$ is increasing, or every $g_i$ is decreasing (part a); and

$$X\le_{sm}Y \implies X_I\le_{sm}Y_I \quad\text{for every } I\subseteq\{1,\dots,n\}$$

(part c, closure under marginalization). These are two of Theorem 9.A.9's five closure properties — the theorem the brief for this mission recommends as its goal — chosen as the ones that need no further definitional machinery beyond the orders themselves (composition and a coordinate projection, both plain pushforwards of the vector's law).

Supporting milestone

• Theorem 9.A.4, the general multivariate closure of the PQD order under independent pairing: if $X\le_{PQD}Y$, $U\le_{PQD}V$, with $X,U$ independent and $Y,V$ independent, then $(\varphi_1(X_1,U_1),\dots,\varphi_n(X_n,U_n))\le_{PQD}(\varphi_1(Y_1,V_1),\dots,\varphi_n(Y_n,V_n))$ for all increasing $\varphi_i$ — the closure property the chapter opens with (as Theorem 9.A.1, the bivariate case), generalized to $n$ dimensions.

Significance

Positive dependence orders formalize a comparison operations research and risk management need constantly but rarely state precisely: two portfolios, insurance lines, or queueing networks with identical individual risk profiles can still differ sharply in how their components co-move, and that co-movement — not the marginals — is often what drives tail risk, correlated failures, or aggregate variability. The supermodular order is the standard tool for comparing this directly: Theorem 9.A.9's closure properties are what let a comparison established for primitive components survive the operations a model actually performs on them — relabeling by a monotone transform (part a), reading off a subset of coordinates (part c), combining independent sub-vectors (part b), conditioning on a covariate (part d), or passing to a distributional limit (part e). Chapter VI's multivariate stochastic order and Chapter VII's multivariate convex order are the two other multivariate orders in the book; the supermodular order is the one built specifically to compare dependence structure rather than location or spread, and it connects directly to the platform's existing Topkis/Supermodularity series, whose function-level predicate this mission reuses rather than restates.

Formalizing Theorem 9.A.9 fixes the exact shape a supermodular-order closure claim takes when built as a pushforward of the vector's law — the natural Mathlib-idiomatic rendering once the order is stated on laws rather than on variables tied to a fixed ambient space — for any future mission in this series or elsewhere that needs to state a closure property of a vector-valued stochastic order. No proof is attempted; both drafted parts of Theorem 9.A.9 have short book proofs (part (a) from the fact that composing a supermodular function with all-increasing or all-decreasing coordinate maps is again supermodular; part (c) the book calls "easy to prove"), making them plausible future proof targets.

Difficulty

The chapter's central definitional trap is that the PQD order's defining condition — even in its general multivariate form — determines the marginals as a consequence of the two joint inequalities holding together, not as a separate hypothesis to add; adding same-marginals as an extra explicit condition on top of Eqs. (9.A.13)-(9.A.14) would not be false, but it would present as a hypothesis something the theorem's own two conditions already force. A second trap is specific to Theorem 9.A.9(a): the $g_i$ must be uniformly increasing or uniformly decreasing across all coordinates, not an arbitrary per-coordinate mix — a mixed-monotonicity version is a materially different, unproven claim, since it is exactly the all-same-direction condition that keeps a composition with a supermodular function supermodular. A third, more basic trap common to every function-class order in this book: "for every supermodular $\varphi$" must be a genuine universal quantifier with the two expectations' existence stated as an explicit hypothesis, not a global assumption or a specific test function standing in for the whole class.

Formalization scope

Unlike the univariate orders formalized elsewhere in this series (which take random variables $X:\Omega\to\mathbb{R}$, $Y:\Omega'\to\mathbb{R}$ on two possibly-different probability spaces), PQDOrder and SupermodularOrder here are stated directly on the two vectors' laws — measures $P,Q$ on $\iota\to\mathbb{R}$ for a finite index type $\iota$ — since the book's own comparisons never reference a joint law of $X$ and $Y$ together, only their separate distributions. This is also what lets Theorem 9.A.9(a)'s coordinatewise composition and (c)'s marginalization be stated as plain pushforwards (Measure.map) of one law, without carrying an ambient sample space through the statement. The index type is left an arbitrary finite type (Fintype ι, not fixed to Fin n) so the same two declarations serve both a full $n$-vector and any marginal sub-vector obtained by restricting to a subset $I$ of coordinates — the shape Theorem 9.A.9(c) itself needs. Supermodularity reuses the platform's own Supermodularity.Monotonicity.SupermodularOn (from the Topkis series, via a kind: reference item) with its relativizing set argument fixed to Set.univ, its plain unrelativized form — checked to match the book's φ(x)+φ(y)≤φ(x∧y)+φ(x∨y) on ι → ℝ's own coordinatewise lattice structure exactly, not a games- or player-scoped variant. Independence in Theorem 9.A.4 is formalized by taking the joint law of the two independent vectors to be the product measure of their marginals — the standard way to construct an independent coupling with prescribed marginals — rather than adding a separate independence hypothesis about pre-existing random variables. Measurable hypotheses are added on every map a Measure.map is taken along, guarding against the pushforward's junk-zero-measure convention for a non-measurable map; none narrows what the book's own theorems claim, since every map used (a monotone/antitone composition, a coordinate projection) is genuinely measurable. A trivializing formalization is ruled out explicitly: the supermodular quantifier ranges over the reused, already-audited Topkis predicate rather than a hand-rolled or restricted one, and the "for all increasing $\varphi_i$" of Theorem 9.A.4 is a genuine universal over jointly-monotone functions of two arguments, not a fixed example.

This mission draws on the platform's own Topkis/Supermodularity series for its supermodular function-level building block (Supermodularity.Monotonicity.SupermodularOn, reused as a reference item — the strongest prior-art connection of the whole nine-chapter series, per the series plan) but on no prior art for the orders themselves (repeated searches for "PQD", "positive dependence", "Fréchet bound" and "supermodular" returned nothing on-topic beyond the Topkis series as of 2026-09-18). Reusable beyond this mission: the laws-only, arbitrary-Fintype- index formalization pattern is available to any later mission in this series needing a multivariate order (Chapter VI's multivariate stochastic order, Chapter VII's multivariate convex order) that wants marginalization or coordinatewise composition stated as plain pushforwards.

Selected references

• M. Shaked and J. G. Shanthikumar, Stochastic Orders, Springer Series in Statistics, Springer, 2007. https://doi.org/10.1007/978-0-387-34675-5
• E. L. Lehmann, "Some concepts of dependence", Annals of Mathematical Statistics, 37(5), 1966, 1137–1153. https://doi.org/10.1214/aoms/1177699260

Motivation

Principal component analysis (PCA) is one of the oldest and most widely used tools in multivariate statistics: given data with covariance matrix $\Sigma$, project onto the top eigenvector(s) of $\Sigma$ to find the directions of maximal variance. In practice one never observes $\Sigma$ itself, only a perturbed version — a sample covariance matrix $\hat\Sigma = \Sigma + P$, with $P$ the (random) estimation error. The natural question, asked since at least Davis and Kahan (1970) and Wedin (1972), is: how close is the top eigenvector of $\hat\Sigma$ to that of $\Sigma$? Wainwright's High-Dimensional Statistics: A Non-Asymptotic Viewpoint (Cambridge University Press, 2019), Chapter 8, gives a self-contained, sharp, non-asymptotic answer, phrased entirely in terms of two deterministic quantities: the eigengap of $\Sigma$, and a single scalar summarizing how the perturbation $P$ couples to the top eigendirection.

Unlike most of the results in this book series, this one is a statement of pure linear algebra: no probability, no concentration inequality, no sample size is needed to state or prove it. Randomness enters only afterward, when $P$ is instantiated as an actual sampling error and bounded using the machinery of earlier chapters (Corollary 8.7, out of this mission's scope).

Setting

Let $\Sigma \in \mathbb R^{d\times d}$ be a symmetric positive semidefinite matrix. Say $\theta \in \mathbb R^d$ is a maximal unit eigenvector of $\Sigma$ if $\|\theta\|_2=1$ and $\theta$ maximizes the Rayleigh quotient $\langle\theta,\Sigma\theta\rangle$ over the whole unit sphere $\mathcal S^{d-1}$ — the variational characterization of the top eigenvector/eigenvalue pair, matching Eq. (8.14) of the book. Write $\gamma_1(\Sigma) := \langle\theta^*,\Sigma\theta^*\rangle$ for the corresponding top eigenvalue. Say $\Sigma$ has eigengap $\nu>0$ at $\theta^*$ if every unit vector $v$ orthogonal to $\theta^*$ satisfies $\langle v,\Sigma v\rangle \le \gamma_1(\Sigma) - \nu$ — the Courant-Fischer variational form of the book's $\nu := \gamma_1(\Sigma)-\gamma_2(\Sigma)$.

For a symmetric perturbation matrix $P \in \mathbb R^{d\times d}$, write $|\!|\!|P|\!|\!|_2 := \sup_{\|v\|_2=1}|\langle v,Pv\rangle|$ for its $\ell_2$-operator norm, and

$$\tilde p \;:=\; P\theta^* - \langle P\theta^*,\theta^*\rangle\,\theta^*$$

for the component of $P\theta^*$ orthogonal to $\theta^*$ — the piece of the perturbation that actually couples the top eigendirection to the rest of the space (Eq. (8.11)). Note $\|\tilde p\|_2$ can be far smaller than $|\!|\!|P|\!|\!|_2$: a perturbation can be large in every direction yet barely move the top eigenvector, if its interaction with $\theta^*$ specifically is small.

Formalization targets

Goal (Theorem 8.5)

Let $\Sigma$ be symmetric positive semidefinite with maximal unit eigenvector $\theta^*$ and eigengap $\nu>0$. For any symmetric $P$ with $|\!|\!|P|\!|\!|_2 < \nu/2$, and any maximal unit eigenvector $\hat\theta$ of $\hat\Sigma := \Sigma+P$ with $\langle\hat\theta,\theta^*\rangle \ge 0$,

$$\|\hat\theta - \theta^*\|_2 \;\le\; \frac{2\|\tilde p\|_2}{\nu - 2|\!|\!|P|\!|\!|_2}.$$

Milestone (Lemma 8.6, the PCA basic inequality)

Under the same eigengap hypothesis, with $\Psi(\Delta;P) := \langle\Delta,P\Delta\rangle

• 2\langle\Delta,P\theta^\rangle$and$\Delta := \hat\theta-\theta^$,

$$\nu\bigl(1-\langle\hat\theta,\theta^*\rangle^2\bigr) \;\le\; |\Psi(\Delta;P)|.$$

Significance

The bound isolates exactly what drives eigenvector instability: not the raw size of the perturbation but its projection onto the top eigendirection, rescaled by the inverse eigengap. This explains, in one inequality, the qualitative phenomenon Example 8.4 illustrates numerically (a tiny perturbation can move the eigenvector far when the eigengap is small) and quantifies exactly how far. It is the deterministic engine behind every consistency result for PCA in the rest of the chapter: Corollary 8.7 (rates for the spiked covariance model) and the sparse-PCA guarantee (Theorem 8.10) both specialize this same bound, after bounding $\|\tilde p\|_2$ and $|\!|\!|P|\!|\!|_2$ using concentration for a specific random design. It is also a sharp instance of the general Davis-Kahan-type perturbation theory for symmetric matrices, phrased with an explicit, non-asymptotic constant rather than an $O(\cdot)$.

The theorem is already proved in the source; this mission's contribution is a machine-checked formal statement (and, eventually, proof) of the bound and its supporting basic inequality, composing with any future formalization of the chapter's probabilistic corollaries.

Difficulty

The proof is genuinely variational, not spectral: it never diagonalizes $\Sigma$ or $\hat\Sigma$, only uses that $\theta^*$ and $\hat\theta$ are optimal for their respective Rayleigh-quotient maximizations. The one place a naive argument fails is in bounding $|\Psi(\Delta;P)|$ itself (the proof of Lemma 8.6): a direct Cauchy-Schwarz bound on $\langle\Delta,P\Delta\rangle$ using $|\!|\!|P|\!|\!|_2$ alone would produce a bound in terms of $|\!|\!|P|\!|\!|_2$ throughout, not the sharper $\|\tilde p\|_2$ the theorem actually delivers; getting the sharper dependence requires decomposing $\Delta$ along $\theta^*$ and its orthogonal complement and tracking the two pieces separately (the $\varrho$, $z$ decomposition on p. 244). The sharpness of the threshold $|\!|\!|P|\!|\!|_2 < \nu/2$ is also not a proof artifact: the book's own $2\times 2$ example ($\Sigma=\mathrm{diag}(2,1)$, $P=\mathrm{diag}(-1/2,1/2)$) shows the perturbed matrix can lose a unique maximal eigenvector exactly at $|\!|\!|P|\!|\!|_2 = \nu/2$.

Formalization scope

$\Sigma$, $P$, $\hat\Sigma=\Sigma+P$ are Matrix (Fin d) (Fin d) ℝ; "symmetric" is M.transpose = M; "positive semidefinite" is the quadratic-form condition $\langle v, \Sigma v\rangle \ge 0$ for all $v$, stated directly rather than via a Mathlib PosSemidef typeclass. "Maximal unit eigenvector" and "eigengap" are both given their variational (Rayleigh-quotient) characterizations rather than defined through Mathlib's enumerated matrix-eigenvalue API — mathematically equivalent to the book's spectral definitions by the Courant-Fischer theorem, and the same characterization the book's own proof works with throughout (Eq. (8.14)). The operator norm $|\!|\!|\cdot|\!|\!|_2$ is likewise realized variationally, valid because this mission only ever applies it to symmetric matrices, exactly as the book does in this chapter. $\tilde p$ is realized as a basis-independent vector in $\mathbb R^d$ (the component of $P\theta^*$ orthogonal to $\theta^*$) rather than the book's basis-dependent $\mathbb R^{d-1}$ representative; its $\ell_2$-norm — the only quantity the theorem's conclusion uses — is identical either way.

Deliberate scope decision, disclosed here rather than silently: the book's own Theorem 8.5 concludes that $\hat\Sigma$ "has a unique maximal eigenvector $\hat\theta$" satisfying the bound — asserting both existence and uniqueness as part of the theorem, on top of the quantitative bound. This mission formalizes only the quantitative bound, for an arbitrary maximal unit eigenvector $\hat\theta$ of $\hat\Sigma$ satisfying the sign condition $\langle\hat\theta,\theta^*\rangle\ge0$ — exactly what the book's own proof establishes (the proof never separately argues existence or uniqueness; both are consequences that could be derived from the bound together with a compactness argument for existence, left to a future extension) and exactly what every downstream use in the chapter (Examples, Corollary 8.7) actually invokes. A trivializing formalization would instead drop the sharp threshold $|\!|\!|P|\!|\!|_2 < \nu/2$ to ≤, or conflate the general operator norm with the symmetric-matrix Rayleigh-quotient characterization for a non-symmetric $P$; this mission does neither. Corollary 8.7 (spiked covariance rates) and Theorem 8.10 (sparse PCA, needing the uniform deviation condition of Eq. (8.26)) are out of this mission's scope; a companion mission formalizing them on top of this one's pca_eigenvector_perturbation_bound is natural future work, together with a proof of existence of a maximal unit eigenvector via compactness of the sphere.

Selected references

• Wainwright, M. J. High-Dimensional Statistics: A Non-Asymptotic Viewpoint. Cambridge University Press, 2019. Chapter 8. DOI: 10.1017/9781108627771.
• Davis, C., Kahan, W. M. "The rotation of eigenvectors by a perturbation. III." SIAM Journal on Numerical Analysis, 7(1), 1970, 1–46.
• Wedin, P.-Å. "Perturbation bounds in connection with singular value decomposition." BIT Numerical Mathematics, 12(1), 1972, 99–111.