Probability

Motivation

A retailer that carries several substitutable products (brands of the same good, sizes of the same item) faces demands that are linked in two ways: a customer who finds the first-choice product out of stock may buy another one, and the first-choice demands themselves are statistically dependent. Inventory planners routinely assume multivariate normal demand and need to know how the dependence between demands affects what the assortment earns. In the single-product newsvendor and in classical risk-pooling arguments, positive correlation is known to hurt; with substitution it is not obvious, because substitution is itself a form of pooling whose value depends on whether two products run out together.

Netessine and Rudi (SSRN 303779, 2002; journal version in Operations Research 51(2), 2003) settle the question for the centralized model with an arbitrary number of products: under multivariate normal demand, the retailer's expected profit decreases in every demand correlation. The authors note that this effect had been observed numerically by Rajaram and Tang (working paper, UCLA, 1999) and Ernst and Kouvelis (1999) but not proved. They prove the analogous statement for each firm in the competitive model with stocking quantities held fixed.

Setting

There are $n$ products. Product $i$ is stocked in quantity $Q_i\ge 0$ at unit cost $c_i$, sold at price $r_i$, and salvaged at value $s_i$, with $r_i>c_i>s_i>0$. Write $u_i=r_i-c_i$ (underage cost) and $o_i=c_i-s_i$ (overage cost). The substitution fraction $a_{ij}\in[0,1]$ is the share of the unmet demand for product $i$ that switches to product $j$; $a_{ii}=0$ and $\sum_j a_{ij}<1$.

The first-choice demand is a random vector $D=(D_1,\dots,D_n)$. The effective demand of product $i$, after substitution from the other products, is

$$D^s_i = D_i + \sum_{j\ne i} a_{ji}\,(D_j-Q_j)^+ ,\qquad x^+=\max(0,x).$$

For one realization, product $k$ earns

$$\pi_k(D) = u_kD^s_k - u_k(D^s_k-Q_k)^+ - o_k(Q_k-D^s_k)^+ ,$$

and a single centralized company owning all products earns $\pi(D)=\sum_k\pi_k(D)$. The expected profits are $\pi(Q)=\mathbb E\,\pi(D)$ (the paper's equation (1)) and $\pi_k(Q)=\mathbb E\,\pi_k(D)$ (equation (9), the profit of firm $k$ when each product is owned by a separate firm).

In this mission $D\sim N(m,\Sigma)$ is multivariate normal with mean $m\in\mathbb R^n$ and positive semidefinite covariance matrix $\Sigma=(\sigma_{ij})$. The correlation of products $i\ne j$ is $\rho_{ij}=\sigma_{ij}/\sqrt{\sigma_{ii}\sigma_{jj}}$. Raising $\rho_{ij}$ means passing from $\Sigma^1$ to $\Sigma^2$ where $\sigma^2_{ij}=\sigma^2_{ji}\ge\sigma^1_{ij}=\sigma^1_{ji}$ and every other entry, in particular every variance, is unchanged.

A function $f:\mathbb R^n\to\mathbb R$ is supermodular if $f(x)+f(y)\le f(x\vee y)+f(x\wedge y)$ for the componentwise maximum and minimum, and submodular if $-f$ is supermodular. A random vector $D^1$ is smaller than $D^2$ in the supermodular order, $D^1\le_{sm}D^2$, if $\mathbb Ef(D^1)\le\mathbb Ef(D^2)$ for every supermodular $f$ whose expectations exist.

Formalization targets

Goal: Proposition 2 (p. 7)

For $\Sigma^2$ obtained from $\Sigma^1$ by raising one covariance $\sigma_{ij}$, $i\ne j$,

$$\pi^{\Sigma^2}(Q)\le\pi^{\Sigma^1}(Q)\quad\text{for every } Q\ge 0,$$

and, for every maximizer $Q^1$ of $\pi^{\Sigma^1}$ and every maximizer $Q^2$ of $\pi^{\Sigma^2}$ over $Q\ge0$,

$$\pi^{\Sigma^2}(Q^2)\le\pi^{\Sigma^1}(Q^1).$$

These are the two readings of "stocking quantities held fixed or adjusted optimally".

Milestones

1. The effective demand $h(D)=D_i+\sum_{j\ne i}a_{ji}(D_j-Q_j)^+$ is a valuation (both super- and submodular).
2. $D\mapsto\min(h(D),Q_i)$ is submodular.
3. The realized profit $\pi(D)=\sum_i[(u_i+o_i)\min(D^s_i,Q_i)-o_iQ_i]$, and each $\pi_k(D)$, is submodular in $D$.
4. (Müller and Scarsini, quoted on p. 7) If $\sigma^1_{ii}=\sigma^2_{ii}$ and $\sigma^1_{ij}\le\sigma^2_{ij}$ for all $i,j$, then

$$N(m,\Sigma^1)\le_{sm}N(m,\Sigma^2).$$

5. (Proposition 5, p. 10) For every firm $k$ and every fixed $Q\ge0$, $\pi_k^{\Sigma^2}(Q)\le\pi_k^{\Sigma^1}(Q)$.

Significance

The result. Proposition 2 says that, under normal demand, positive dependence between substitutable products is costly to the retailer whatever the substitution pattern, the number of products, and the cost parameters, and that re-optimizing the stocking levels cannot undo this. It turns an empirical regularity of numerical studies into a structural statement, and it gives a reason to prefer assortments of products whose demands are less correlated. Proposition 5 gives the same conclusion for each competing firm when stocking quantities are fixed; the paper does not claim it when firms re-optimize, since equilibrium quantities move in ways that need not favour every firm.

Formalizing it. The result is proved in the paper, partly by citation: the comparison of normal laws in the supermodular order is a theorem of Müller and Scarsini that is quoted, not proved, there. No machine-checked version of that comparison is known on the platform: the existing Slepian and Sudakov–Fernique comparisons concern expected suprema and tail probabilities, not expectations of arbitrary supermodular functions. A formal proof of milestone 4 would therefore be reusable well beyond inventory theory, wherever supermodular functionals of Gaussian vectors are compared (portfolio risk, dependence orderings, stop-loss premiums).

Difficulty

The deterministic steps (milestones 1–3) are elementary lattice arguments. The difficulty is milestone 4. The natural route is Gaussian interpolation: connect the two laws through $N(m,(1-t)\Sigma^1+t\Sigma^2)$ and differentiate $\mathbb Ef$ in $t$. That argument needs $f$ twice differentiable with controlled growth, while the statement must hold for every supermodular $f$ that is merely integrable under the two endpoint laws, such as the piecewise-linear profits of this model, and for singular covariance matrices. Passing from smooth to general supermodular functions, without assuming integrability along the interpolation path, is where the work lies. The paper's own route to the "adjusted optimally" clause, a total derivative $d\pi/d\rho_{ij}$ along a differentiable optimal path, presumes a smoothness of the optimizer that is not available, since the expected profit need not be concave in $Q$; the goal is stated as a comparison of optimal values, which avoids it.

Formalization scope

All declarations live in DemandSubstitution.Correlation. Products are indexed by Fin n, 0-based. The model data and the standing assumptions $r_i>c_i>s_i>0$, $a_{ij}\in[0,1]$, $a_{ii}=0$, $\sum_ja_{ij}<1$ form a structure Model n; $u_i$ and $o_i$ are definitions. The paper's §2 assumption that demand is continuous with positive support is dropped, as the paper drops it: a normal law charges negative demands. $N(m,\Sigma)$ is Mathlib's multivariateGaussian on EuclideanSpace ℝ (Fin n), pushed forward to Fin n → ℝ; every statement assumes both covariance matrices positive semidefinite, because Mathlib's construction returns a Dirac mass otherwise. Correlations are encoded through covariances with fixed variances, which is equivalent to raising $\rho_{ij}$. "Decreasing" is weak. Expectations are Bochner integrals; the profit integrands are Lipschitz and integrable under every normal law, and this integrability is to be proved, not assumed. Supermodularity and the supermodular order reuse the platform definitions SupermodularOn and SupermodularOrder.

A trivializing formalization is ruled out: the Müller–Scarsini comparison is a milestone theorem, never a hypothesis of the goal or of Proposition 5, and the "adjusted optimally" clause quantifies over every pair of maximizers, not over one chosen optimizer.

Contributions welcome: integrability of Lipschitz functions under multivariateGaussian; Gaussian interpolation for supermodular functions; the supermodular comparison itself, which would be reusable on its own.

Selected references

• S. Netessine and N. Rudi, Centralized and Competitive Inventory Models with Demand Substitution, Simon School Working Paper OP 02-01, 2002, SSRN 303779. https://doi.org/10.2139/ssrn.303779 (published in Operations Research 51(2):329–335, 2003, https://doi.org/10.1287/opre.51.2.329.12788)
• A. Müller and M. Scarsini, Some remarks on the supermodular order, Journal of Multivariate Analysis 73(1):107–119, 2000. https://doi.org/10.1006/jmva.1999.1867
• D. M. Topkis, Supermodularity and Complementarity, Princeton University Press, 1998. https://doi.org/10.1515/9781400822539
• K. Rajaram and C. S. Tang, The impact of product substitution on retail merchandising, working paper, University of California, Los Angeles, 1999 (reference [16] of the source; no stable link).
• D. M. Topkis, Minimizing a submodular function on a lattice, Operations Research 26(2):305–321, 1978. https://doi.org/10.1287/opre.26.2.305
• R. Ernst and P. Kouvelis, The effects of selling packaged goods on inventory decisions, Management Science 45(8):1142–1155, 1999. https://doi.org/10.1287/mnsc.45.8.1142

Motivation

Two-stage stochastic linear programs model a decision $x$ taken before uncertainty is revealed, followed by a corrective recourse action once it is. The classical model needs the full probability distribution of the uncertain data, which in practice is rarely known; often only estimates of the mean and the second moments are available, and the computed plan can be sensitive to the distribution assumed. The minimax (distributionally robust) approach replaces the single distribution by the class of all distributions with given first and second moments and optimizes against the worst member of that class. Bertsimas, Doan, Natarajan and Teo (Math. Oper. Res. 35(3), 2010) combine this with risk aversion, through a piecewise-linear disutility of the second-stage cost, and determine which instances can be solved as semidefinite programs.

When the randomness sits in the right-hand side of the second stage, the paper shows that the problem is NP-hard in general (Theorem 3.1). This mission concerns the tractable special case that follows (§3.2): if the extreme points of the dual feasible region of the recourse problem are known, the whole minimax problem is one semidefinite program.

Setting

A first-stage decision $x\in\mathbb R^n$ is feasible if $x\in X=\{x : Ax=b,\ x\ge 0\}$ and costs $c'x$. Given a recourse matrix $W\in\mathbb R^{r\times d}$, a technology matrix $T\in\mathbb R^{r\times n}$, a constant cost vector $q\in\mathbb R^d$ and a realized right-hand side $h\in\mathbb R^r$, the second-stage cost is written as the optimal value of the dual of the recourse linear program,

$$\mathcal Q(h,x)=\max_{p}\ (h-Tx)'p\quad\text{s.t.}\quad W'p\le q .$$

Risk aversion enters through the disutility $\mathbb U(t)=\max_{k=1,\dots,K}(\alpha_k t+\beta_k)$ with $\alpha_k\ge 0$. The right-hand side $\tilde h$ is random, and its distribution is only known to lie in the moment class

$$\mathcal P=\{P : \mathbb E_P[\tilde h]=\mu,\ \mathbb E_P[\tilde h\tilde h']=Q\}.$$

The worst-case expected disutility is $\hat Z(x)=\sup_{P\in\mathcal P}\mathbb E_P[\mathbb U(\mathcal Q(\tilde h,x))]$, and the problem is

$$Z=\min_{x\in X}\big(c'x+\hat Z(x)\big).$$

The moment problem defining $\hat Z(x)$ has the dual (12): minimize $Q\cdot Y+\mu'y+y_0$ over symmetric $Y$, $y\in\mathbb R^r$, $y_0\in\mathbb R$ with $h'Yh+y'h+y_0\ge\mathbb U(\mathcal Q(h,x))$ for every $h$, where $Q\cdot Y=\sum_{i,j}Q_{ij}Y_{ij}$.

The standing assumptions are: $X$ nonempty and bounded (Assumption 1); complete recourse, $\{Ww : w\ge 0\}=\mathbb R^r$ (Assumption 2); a nonempty dual polyhedron $\{p : W'p\le q\}$ (Assumption 3); $Q\succ\mu\mu'$ (Assumption 4); and the extreme points $p_1,\dots,p_N$ of $\{p : W'p\le q\}$ are known (Assumption 5).

Formalization targets

Goal: Theorem 3.2 (p. 589)

The problem is equivalent to the semidefinite program

$$\hat Z_{\rm SDP}=\min_{x,Y,y,y_0}\ c'x+Q\cdot Y+\mu'y+y_0\ \ \text{s.t.}\ \begin{pmatrix}Y & \tfrac12(y-\alpha_k p_i)\\ \tfrac12(y-\alpha_k p_i)' & y_0+\alpha_k p_i'Tx-\beta_k\end{pmatrix}\succeq 0\ \ \forall k,i,\quad Ax=b,\ x\ge 0, \tag{14}$$

in the sense that $\hat Z(x)$ equals the value of the per-$x$ program (15) for every $x\in X$, and $Z=\hat Z_{\rm SDP}$, both being genuine finite infima.

Milestones (§3.1 and the proof of Theorem 3.2, pp. 588–590)

1. Strong duality for (11)/(12): $\hat Z(x)=\hat Z_D(x)$.
2. The constraints of (12) are equivalent to $h'Yh+(y-\alpha_kp)'h+y_0+\alpha_kp'Tx-\beta_k\ge 0$ for all $h$, all $p$ with $W'p\le q$ and all $k$.
3. Every feasible $Y$ of (12) is positive semidefinite.
4. $\mathcal Q(h,x)=\max_{i=1,\dots,N}(h-Tx)'p_i$.
5. A quadratic $h'Yh+2v'h+c$ is nonnegative on $\mathbb R^r$ iff the bordered matrix $\begin{pmatrix}Y&v\\v'&c\end{pmatrix}$ is positive semidefinite (already proved on the platform).
6. The feasible sets of (12) and (15) coincide, so $\hat Z_D(x)$ is the value of (15).

Significance

The result. Theorem 3.2 turns a two-level problem — a minimum over $x$ of a supremum over an infinite-dimensional set of distributions — into one convex program of polynomial size in $n$, $r$, $K$ and $N$, jointly in the first-stage decision and the dual multipliers. The first-stage decision enters the linear matrix inequalities only through the corner entry $y_0+\alpha_kp_i'Tx-\beta_k$, which is affine in $x$; this is what keeps the joint program convex. Together with the NP-hardness of the general case (Theorem 3.1), it locates the source of intractability in the enumeration of dual vertices. The same program (15) is the starting point of the construction of an extremal distribution (Theorem 3.3), which the paper uses to stress-test first-stage solutions.

Formalizing it. The result is proved in the paper, with strong duality for the moment problem cited from the literature. Nothing in this chain is formalized: neither the moment-problem duality under the interior condition $Q\succ\mu\mu'$, nor the vertex representation of a bounded linear program over a pointed polyhedron given by its extreme points. The lemma that a nonnegative quadratic corresponds to a positive semidefinite bordered matrix is already proved on the platform and is reused.

Difficulty

The reduction from (12) to (15) is elementary once $\mathcal Q(h,x)=\max_i(h-Tx)'p_i$ is available, but that identity is a genuine statement about polyhedra: the dual region $\{p : W'p\le q\}$ is unbounded, so compactness arguments such as Krein–Milman do not apply, and one needs both that the polyhedron is pointed (here a consequence of complete recourse, which forces $W$ to have full row rank) and that a bounded linear objective attains its maximum at a vertex. The other substantial step is strong duality between the moment problem over measures and the semi-infinite dual (12), which requires the moment vector $(\mu,Q)$ to lie in the interior of the moment cone; without $Q\succ\mu\mu'$ a duality gap or an empty class $\mathcal P$ can occur. Finally, the equality $Z=\hat Z_{\rm SDP}$ requires exchanging an infimum over $x$ with the infimum over the dual variables, which needs each per-$x$ program to be feasible and bounded.

Formalization scope

Vectors are Fin n → ℝ, matrices Matrix (Fin r) (Fin d) ℝ, and pieces and extreme points are indexed by Fin K and Fin N (zero-based). The paper's letter $p$ is both the dimension of $q$ and the dual variable; in Lean the dimension is $d$ and the dual vector is π. The class $\mathcal P$ consists of probability measures on $\mathbb R^r$ with square-integrable coordinates, mean $\mu$ and uncentred second-moment matrix $Q$, built on the published MomentDRO.Conf.Setting; the paper's densities are replaced by measures. $\mathcal Q(h,x)$ is the real supremum over the dual polyhedron. Every printed min/max is a real infimum/supremum of a set of objective values, and the statements assert with IsLUB/IsGLB that these are genuine finite bounds; attainment of the minima is not claimed. $Y$ is required symmetric throughout. The bordered matrix is the published WorstCaseVaR.KnownMoments.bordered and $\succeq 0$ is Matrix.PosSemidef.

Readings of the source:

• As printed, Assumption 2 (complete recourse) and Assumption 3 ("for all $q$") contradict each other when $r\ge 1$. The section-wise reading is used: complete recourse together with dual feasibility at the single constant $q$ of §3.
• "Assumption i" in Theorem 3.2 and in its proof is read as Assumption 5.
• "The strong duality condition" is read as Assumption 4.
• Assumption 5 is the hypothesis that the list $p_1,\dots,p_N$ is exactly the set of extreme points of the dual region; "explicitly known" has no further mathematical content.
• (11) integrates over $\mathbb R^p$ where $\mathbb R^r$ is meant.

A formalization in which $\hat Z_{\rm SDP}$ is defined through $\hat Z_D(x)$, or in which positive semidefiniteness is replaced by a hand-written quadratic condition, would make the goal a restatement of strong duality; (14) is stated as a single program in $(x,Y,y,y_0)$ with Matrix.PosSemidef constraints. The hypotheses of the goal are satisfiable ($r=1$, $W=[1,-1]$, $q=(1,1)$, $\mu=0$, $Q=1$, $X=\{1\}$).

Contributions welcome: the vertex representation of bounded linear programs over pointed polyhedra, and duality for moment problems with first and second moments, are both reusable well beyond this mission.

Selected references

• D. Bertsimas, X. V. Doan, K. Natarajan, C.-P. Teo, Models for Minimax Stochastic Linear Optimization Problems with Risk Aversion, Mathematics of Operations Research 35(3):580–602, 2010. https://doi.org/10.1287/moor.1100.0445
• K. Isii, On sharpness of Tchebycheff-type inequalities, Annals of the Institute of Statistical Mathematics 14:185–197, 1962. https://doi.org/10.1007/BF02868641
• A. Shapiro, On duality theory of conic linear problems, in Semi-Infinite Programming, Kluwer, 2001, pp. 135–165. https://doi.org/10.1007/978-1-4757-3403-4_7
• L. El Ghaoui, M. Oks, F. Oustry, Worst-case value-at-risk and robust portfolio optimization: a conic programming approach, Operations Research 51(4):543–556, 2003. https://doi.org/10.1287/opre.51.4.543.16101

Motivation

Linear decision rules are a standard device for making multistage optimization under uncertainty computationally tractable: a recourse decision is restricted to be an affine function of the uncertain data observed so far. Their weakness is that a recourse variable that must stay nonnegative forces the rule to be very conservative. Goh and Sim, Distributionally Robust Optimization and Its Tractable Approximations (Operations Research, 2010), repair this with deflected linear decision rules, in which the violated part of an affine rule is corrected and paid for. The cost of the correction is the expected positive part $\mathbb E_{\mathbb P}((r^0 + r'\tilde\zeta)^+)$ of an affine function of the uncertainty.

The distribution $\mathbb P$ is not known. The paper assumes only partial distributional information (a support set, a set containing the mean, a covariance, directional deviations) and bounds the worst case of this expectation over every distribution consistent with that information. Such worst-case bounds go back to Scarf's 1958 mean–variance newsvendor bound and form the core of distributionally robust optimization. This mission formalizes the first bound of §5, which uses support and mean information only.

Setting

Let $\tilde\zeta$ be a random vector in $\mathbb R^{n}$ (in the paper, the segregated uncertainty vector of dimension $n = N_E$). Two sets describe what is known about it:

• a support set $\mathcal V \subseteq \mathbb R^n$ with $\mathbb P(\tilde\zeta \in \mathcal V) = 1$;
• a mean set $\widehat{\mathcal V} \subseteq \mathbb R^n$ containing the mean $\hat\zeta = \mathbb E_{\mathbb P}(\tilde\zeta)$. A known mean is the case where $\widehat{\mathcal V}$ is a single point.

The family of distributions consistent with this information is

$$\mathbb F_1 = \{\mathbb P : \hat\zeta = \mathbb E_{\mathbb P}(\tilde\zeta) \in \widehat{\mathcal V},\ \mathbb P(\tilde\zeta \in \mathcal V) = 1\}.$$

For data $r^0 \in \mathbb R$ and $r \in \mathbb R^n$ the worst-case value is $\sup_{\mathbb P \in \mathbb F_1} \mathbb E_{\mathbb P}((r^0 + r'\tilde\zeta)^+)$, with $x^+ = \max\{x, 0\}$. The paper's candidate for this value is the deterministic quantity of (21),

$$\pi^1(r^0, r) = \inf_{s \in \mathbb R^{n}} \Phi(s), \qquad \Phi(s) = \sup_{\hat\zeta \in \widehat{\mathcal V}} \{s'\hat\zeta\} + \sup_{\zeta \in \mathcal V} \max\{r^0 + r'\zeta - s'\zeta,\ -s'\zeta\}.$$

Both quantities live in the extended reals: following the paper (p. 909), an infeasible maximization has value $-\infty$ and an infeasible minimization $+\infty$.

Formalization targets

Goal: Theorem 1 (p. 909)

For closed convex $\mathcal V, \widehat{\mathcal V}$ with $\widehat{\mathcal V} \cap \operatorname{int} \mathcal V \neq \emptyset$, and all $r^0, r$,

$$\sup_{\mathbb P \in \mathbb F_1} \mathbb E_{\mathbb P}\big((r^0 + r'\tilde\zeta)^+\big) = \pi^1(r^0, r).$$

Milestones

1. Upper-bound part of Theorem 1. $\sup_{\mathbb P \in \mathbb F_1} \mathbb E_{\mathbb P}((r^0 + r'\tilde\zeta)^+) \le \pi^1(r^0, r)$, for arbitrary sets $\mathcal V, \widehat{\mathcal V}$.
2. Remark after Theorem 1, first case (p. 910). If $\widehat{\mathcal V} \subseteq \mathcal V$ is nonempty and $r^0 + r'\zeta \ge 0$ on $\mathcal V$, then the worst-case value equals $r^0 + \sup_{\hat\zeta \in \widehat{\mathcal V}} r'\hat\zeta$, and $\Phi(r) = \pi^1(r^0, r)$.
3. Remark after Theorem 1, second case (p. 910). If $\widehat{\mathcal V} \subseteq \mathcal V$ is nonempty and $r^0 + r'\zeta \le 0$ on $\mathcal V$, then the worst-case value is $0$, and $\Phi(0) = \pi^1(r^0, r) = 0$.
4. Tightness part of Theorem 1. Under the goal's hypotheses, $\pi^1(r^0, r) \le \sup_{\mathbb P \in \mathbb F_1} \mathbb E_{\mathbb P}((r^0 + r'\tilde\zeta)^+)$.

Significance

The result itself. Theorem 1 converts an optimization over an infinite-dimensional set of probability measures into a finite-dimensional convex program over the given sets. When $\mathcal V$ and $\widehat{\mathcal V}$ are representable by linear and second-order cone constraints, the epigraph of $\pi^1$ is a pair of robust constraints over those sets and is again tractable (Remark, p. 910). The bound is tight, so no better bound can be built from support and mean information alone. In the paper, $\pi^1$ is one of three building blocks; their infimal convolution (Theorem 4) prices the deflection in the deflected linear decision rules of §6, and the experiments there compare these rules against plain linear decision rules.

Formalizing it. The proof of Theorem 1 is in the paper's online appendix, not in the main text. A search of the Prove2Me library (October 2026) found no formal statement of this or any other mean-and-support bound. A formal proof has to settle the extended-real bookkeeping (empty families, unbounded suprema) and the exact hypotheses under which tightness holds. As printed, the theorem omits a condition: with no constraint qualification it fails when $\mathbb F_1$ is empty, which can happen for disjoint closed convex sets. This mission states the corrected version.

Difficulty

The upper-bound part is weak duality and needs no structure on the sets. The difficulty is tightness. It amounts to a strong-duality statement between the worst case over distributions and the infimum over the multiplier $s$, and strong duality for convex problems over unbounded sets needs a constraint qualification. Two obvious approaches fail. The first is to search for a single worst-case distribution: for unbounded $\mathcal V$ the supremum over $\mathbb F_1$ is generally not attained, and it may be $+\infty$. The second is to restrict to compact $\mathcal V$, which excludes the paper's main use, a segregation of an unbounded support such as $(-\infty, -1] \cup [-1, 1] \cup [1, \infty)$. Mathlib has separating-hyperplane theorems but no Fenchel duality theorem, so this duality has to be built.

Formalization scope

• A distribution is a measure $P$ on Fin n → ℝ with IsProbabilityMeasure P, and the random vector is the identity map. Every coordinate is required to be $P$-integrable, so the mean is a genuine expectation; the mean is MomentDRO.Conf.meanVec from the published definition MomentDRO_Conf_Setting.
• $\mathbb P(\tilde\zeta \in \mathcal V) = 1$ is the almost-everywhere statement ∀ᵐ ζ ∂P, ζ ∈ V, which needs no measurability of $\mathcal V$.
• The positive part is max (r0 + r ⬝ᵥ ζ) 0.
• The worst-case value, both suprema and the infimum in (21) are taken in EReal. An empty supremum is $-\infty$, an empty infimum $+\infty$, and an unbounded supremum $+\infty$. The paper's "inf" is ⨅; attainment is claimed only where the Remark claims it, via the objective $\Phi$ (pi1Obj) at $s = r$ and $s = 0$.
• The paper's tractable conic representable sets are modelled as closed and convex sets. The goal and its tightness part also assume $\widehat{\mathcal V} \cap \operatorname{int}\mathcal V \ne \emptyset$. This condition is not on the page; it is consistent with the paper's full-dimensionality assumptions (pp. 905–906), and it implies $\mathbb F_1 \neq \emptyset$.
• The Remark's two cases assume $\widehat{\mathcal V} \subseteq \mathcal V$ and $\widehat{\mathcal V} \ne \emptyset$. The paper's construction of $\mathcal V, \widehat{\mathcal V}$ from one set $\mathcal H$ (§4.5) gives the inclusion, and without it the first case is false.
• The segregation map, the matrices $F$, $g$ and the set $\mathcal H$ of §4 are context only and are not formalized: the theorem is stated for an arbitrary random vector with the listed properties.
• A real-valued supremum, a family without integrability (for which the mean would be a junk value), or a $\pi^1$ defined through the worst case would each trivialize the statement; the definitions avoid all three, and $\pi^1$ is defined directly from the sets.

Needed infrastructure: finite mixtures of Dirac masses as members of $\mathbb F_1$, integration of almost-surely bounded functions, separation of convex sets in $\mathbb R^n$, and a Fenchel-type duality for a concave function over a closed convex set under an interior-point condition. The duality and the concave-envelope facts are reusable well beyond this mission, including for Theorems 2–4 of the same paper. Proofs of any milestone, and reusable lemmas on extended-real suprema of integrals, are welcome.

Selected references

• J. Goh and M. Sim, Distributionally Robust Optimization and Its Tractable Approximations, Operations Research 58(4, Part 1 of 2):902–917, 2010. https://doi.org/10.1287/opre.1090.0795
• H. Scarf, A Min–Max Solution of an Inventory Problem, in K. Arrow, ed., Studies in the Mathematical Theory of Inventory and Production, Stanford University Press, 1958, 201–209.
• W. Chen and M. Sim, Goal-Driven Optimization, Operations Research 57(2):342–357, 2009 (cited by the paper for the infimal convolution of bounds).
• R. T. Rockafellar, Convex Analysis, Princeton University Press, 1970 (Theorem 31.1, Fenchel duality). https://doi.org/10.1515/9781400873173

Motivation

The sample average approximation (SAA) method replaces an expectation in a stochastic program by the average over a random sample and solves the resulting deterministic problem. It is the standard way to solve two-stage stochastic linear programs with large or continuous scenario sets, and it underlies the stochastic counterpart and sample-path optimization methods of the 1990s. A basic question is how large the sample has to be.

For smooth problems the SAA optimal solution $\hat x_N$ converges to the true one at the rate $N^{-1/2}$ of the central limit theorem. A. Shapiro and T. Homem-de-Mello (SIAM J. Optim. 11(1), 2000) observed that for problems whose optimum is a sharp minimum, as is typical of two-stage linear programs with finitely many scenarios, the picture is different. In §2 of the paper they prove that the SAA problem then returns the exact true solution with probability one for all sufficiently large $N$. In §3, formalized here, they show that the probability of the complementary event, the SAA solution not being exactly $\bar x$, decays exponentially in $N$. One consequence they draw is that a modest sample can find the optimum of a problem with an astronomically large number of scenarios. This mission is the second of three on the paper; the first formalizes the almost sure statement (Theorem 2.1), the third the polyhedral case with multiple optima (Theorem 3.2).

Setting

Let $P$ be a probability measure on a sample space $(\Omega,\mathcal F)$, $\Theta\subseteq\mathbb R^m$, and $h:\mathbb R^m\times\Omega\to\mathbb R$. The true problem (1.1) is

$$\min_{x\in\Theta}\ f(x):=\mathbb E_P\,h(x,\omega).$$

Given an i.i.d. sample $\omega^1,\omega^2,\dots$ with law $P$, the SAA problem (1.2) is

$$\min_{x\in\Theta}\ \hat f_N(x):=\frac1N\sum_{j=1}^N h(x,\omega^j).$$

$A$ and $A_N$ denote the sets of optimal solutions of (1.1) and (1.2). For a function $g$, $g'(x,d)=\lim_{t\downarrow0}(g(x+td)-g(x))/t$ is its directional derivative, and $h'_\omega(x,d)$ is that of $h(\cdot,\omega)$. $T_\Theta(\bar x)$ is the tangent cone of $\Theta$ at $\bar x$, and $S^{m-1}=\{d:\|d\|=1\}$.

The standing hypotheses are those of Theorem 2.1 of the paper: (i) $h(\cdot,\omega)$ is convex for every $\omega$; (ii) $f$ is finite valued; (iii) $\Theta$ is closed and convex; (iv) Assumption (A): $A=\{\bar x\}$ and $\bar x$ is a sharp minimum, i.e. there is $c>0$ with

$$f(x)\ge f(\bar x)+c\|x-\bar x\|\qquad\text{for all }x\in\Theta. \tag{2.2}$$

Section 3 adds Assumption (B): with $\eta(d,\omega):=h'_\omega(\bar x,d)$ (3.1), there is $\kappa>0$ such that $\sup_{d\in S^{m-1}}|\eta(d,\omega)|\le\kappa$ for $P$-almost every $\omega$. It holds whenever $\Omega$ is finite.

The event of interest is

$$\mathcal E_N:=\{A_N\neq\emptyset\text{ and }A_N=\{\bar x\}\},\tag{3.7}$$

that the SAA problem has the unique optimal solution $\hat x_N=\bar x$.

For a direction $d$, $M_d(t)=\mathbb E_P e^{t\eta(d,\omega)}$ and $I_d(\alpha)=\sup_{t\in\mathbb R}[t\alpha-\log M_d(t)]$ (3.5) is the one-dimensional rate function of $\eta(d,\omega)$. On the Banach space $Z=C(S^{m-1})$, $I(z)=\sup_{z^*\in Z^*}\{z^*(z)-\log M(z^*)\}$ (3.4) is the rate function of the random element $\eta(\cdot,\omega)$.

Formalization targets

Goal: Theorem 3.1 (p. 11)

Under (i)–(iv) and Assumption (B) there is a constant $\beta>0$ with

$$\limsup_{N\to\infty}\frac1N\log P(\mathcal E_N^c)\le-\beta.\tag{3.8}$$

The goal asserts the existence of an exponential rate and does not fix its value.

Stronger: the explicit rate (p. 12)

With $c$ from (2.2) and $\kappa$ from (B),

$$\limsup_{N\to\infty}\frac1N\log P(\mathcal E_N^c)\le-\frac{c^2}{2\kappa^2},$$

which is (3.8) combined with (3.15).

Milestones

In the order of the proof:

• the reduction $P(\mathcal E_N^c)\le P(\zeta_N\in F)$ with $\zeta_N=\hat f'_N(\bar x,\cdot)$ and $F=\{z:\inf_{T_\Theta(\bar x)\cap S^{m-1}}z\le0\}$ (p. 11);
• for a random variable $X$ with $|X|\le\kappa$ and $\Lambda(t)=\log\mathbb E e^{tX}$, the bounds $|\Lambda'|\le\kappa$, $|\Lambda''|\le\kappa^2$ (3.11) and $|\Lambda'(t)-\Lambda'(s)|\le\kappa^2|t-s|$ (3.12);
• the quadratic lower bound $I_d(\alpha)\ge|\alpha-\bar\alpha_d|^2/(2\kappa^2)$ (3.13);
• $I_d(0)\ge c^2/(2\kappa^2)$ on $T_\Theta(\bar x)\cap S^{m-1}$ (3.14);
• the comparison $I(z)\ge\sup_{d}I_d(z(d))$ (3.6);
• $\beta=\inf_{z\in F}I(z)\ge c^2/(2\kappa^2)$ (3.15).

Significance

Theorem 3.1 says that at a sharp minimum the SAA method is exact except on an event of probability $O(e^{-\beta N})$. The sample size needed for a given confidence therefore grows only logarithmically in the confidence level, and it does not depend on the number of scenarios. This is the basis of the paper's numerical observation that problems with very many scenarios are solved exactly from samples of a few hundred. The explicit constant $c^2/(2\kappa^2)$ relates the rate to two interpretable quantities: the sharpness $c$ of the optimum and the bound $\kappa$ on the subgradients of $h(\cdot,\omega)$ at $\bar x$.

The result is proved in the paper. As far as is known it has no machine-checked proof. A complete formalization produces three things. First, a Lean proof that the SAA problem recovers a sharp optimum exactly with exponentially small failure probability. Second, the elementary facts about cumulant generating functions of bounded variables, (3.11)–(3.13), which are reusable in any Chernoff-type argument. Third, the link between the rate function of a $C(K)$-valued random element and its one-dimensional projections, (3.6).

Difficulty

The event $\mathcal E_N^c$ is controlled by the positivity of a whole function, $\zeta_N(d)=\hat f_N'(\bar x,d)$, on the infinite set $T_\Theta(\bar x)\cap S^{m-1}$. A Chernoff bound for each fixed direction $d$ is easy, but there are uncountably many directions, so a union bound over $d$ is not available. The paper handles this with Cramér's theorem in the separable Banach space $C(S^{m-1})$, and that infinite-dimensional large deviations upper bound is not in Mathlib. A proof must either develop that bound or find a finite-dimensional substitute that keeps an exponential rate.

A second difficulty is conversion. Cramér's theorem gives a rate in terms of $\inf_{z\in F}I(z)$, and that infimum must be shown positive. Positivity comes from the sharpness constant $c$ through the strong convexity of the one-dimensional rate functions, which needs the bounds (3.11)–(3.13) on the cumulant generating function.

Formalization scope

• The decision space is EuclideanSpace ℝ (Fin m). The scenario space $(\Omega,\mathcal F,P)$ carries the law of one scenario. The sample is one i.i.d. sequence of $\Omega$-valued maps on a separate probability space $(S,Q)$, and $\hat f_N$ uses its first $N$ terms. Optimal sets are written without infima: $A_N=\{\bar x\}$ encodes "nonempty and equal to $\{\bar x\}$".
• Directional derivatives are the limUnder of the difference quotient as $t\downarrow0$. Every statement assumes (i), under which the limit exists. The tangent cone is Mathlib's tangentConeAt, which for convex $\Theta$ is the closure of the cone of feasible directions.
• Disclosed pins.
  • "f finite valued" is $P$-integrability of every $h(x,\cdot)$, which includes measurability.
  • Assumption (B) is stated pointwise in $d\in S^{m-1}$ under the almost sure quantifier.
  • $\mathcal E_N^c$ is not shown measurable, so its probability is an outer measure, and an upper bound on it is the stronger claim.
  • The rate (3.8) is stated in the equivalent form: there is $\beta>0$ such that for every $\varepsilon>0$, eventually $Q(\mathcal E_N^c)\le e^{-(\beta-\varepsilon)N}$. A logarithm of a real probability is never taken.
  • Rate functions $I_d$ and $I$ take values in $(-\infty,+\infty]$. They use the published extended-real logMGF and legendre, under which $I_d(\alpha)=+\infty$ for $|\alpha|>\kappa$ is a genuine value, not a default.
  • The one-dimensional bounds (3.11)–(3.13) are stated for an arbitrary bounded measurable random variable.
• Trivializing readings are ruled out. The goal's $\beta$ is fixed before $\varepsilon$ and $N$, so it cannot depend on them. The event is uniqueness and equality, $A_N=\{\bar x\}$, not $\bar x\in A_N$. The standing hypotheses are jointly satisfiable: a sorry-free check exhibits them for $h(x,\omega)=\|x\|$, $\Theta=\mathbb R^m$ and a product-measure sample.
• Infrastructure. A complete development needs: directional derivatives of finite convex functions and their continuity in the direction; interchange of directional derivative and expectation (Proposition 2.2 of the paper, formalized in the first mission of the series); derivatives of the cumulant generating function (Mathlib's ProbabilityTheory.cgf and its analyticity); and an exponential upper bound for empirical means in $C(S^{m-1})$ or a substitute. The cumulant generating function lemmas and any Banach-space Cramér bound are reusable beyond this mission. Contributions are welcome for any milestone, and for alternative proofs of the goal that avoid Cramér's theorem in function space.

Selected references

• A. Shapiro and T. Homem-de-Mello, On the rate of convergence of optimal solutions of Monte Carlo approximations of stochastic programs, SIAM J. Optim. 11(1):70–86, 2000. https://doi.org/10.1137/S1052623498349541 (formalized from the authors' preprint, Stochastic Programming E-Print Series)
• J.-D. Deuschel and D. W. Stroock, Large Deviations, Academic Press, 1989 (Cramér's theorem in Banach spaces, cited as [8]).
• J.-B. Hiriart-Urruty and C. Lemaréchal, Convex Analysis and Minimization Algorithms II, Springer, 1993, Theorem X.4.2.2 (strong convexity of conjugates, cited as [12]). https://doi.org/10.1007/978-3-662-06409-2
• A. J. Kleywegt, A. Shapiro and T. Homem-de-Mello, The sample average approximation method for stochastic discrete optimization, SIAM J. Optim. 12(2):479–502, 2002. https://doi.org/10.1137/S1052623499363220

Motivation

A telephone call center must decide how many agents $N$ to staff against an offered load $R=\lambda/\mu$ (arrival rate over service rate). Staff too few and callers wait and hang up; staff too many and agents sit idle. For large centers the industry rule of thumb is square-root safety staffing, $N\approx R+\beta\sqrt R$, and the question is what service level such a rule buys.

For the model without abandonment, the $M/M/N$ queue, the answer is the theorem of Halfin and Whitt (Oper. Res. 1981): the probability of waiting has a nondegenerate limit in $(0,1)$ exactly when $\sqrt N(1-\rho_N)\to\beta>0$. Real callers are impatient, and Garnett, Mandelbaum and Reiman (M&SOM 2002) prove the analogue for the Erlang-A queue, where abandonment changes both the limits and the admissible range of $\beta$. Their Theorem 4 is the theoretical basis of the "quality-and-efficiency-driven" (QED) staffing regime used in call-center practice.

Timeline. Palm (1943) introduced exponential patience; Riordan (1962) computed the conditional mean wait of the model; Jagerman (1974) gave the asymptotics of the Erlang-B formula; Halfin and Whitt (1981) proved the square-root law for $M/M/N$ ($\beta>0$ only); Fleming, Stolyar and Simon (1994) found the limiting stationary law of the scaled Erlang-A queue; Garnett, Mandelbaum and Reiman (2002) proved the three-way equivalence below, for $\beta$ of either sign.

Setting

The Erlang-A queue ($M/M/N+M$) has Poisson arrivals of rate $\lambda>0$, $N\ge1$ identical agents with exponential service of rate $\mu>0$, an unlimited waiting room served in order of arrival, and independent exponential patience times of rate $\theta>0$: a caller whose wait exceeds their patience abandons. The number in system is a birth–death process with birth rate $\lambda$ and death rate $\min(k,N)\mu+(k-N)^+\theta$; its stationary law is

$$\pi_k\propto\frac{(\lambda/\mu)^k}{k!}\ (k\le N),\qquad \pi_k\propto\prod_{j=N+1}^{k}\frac{\lambda}{N\mu+(j-N)\theta}\cdot\frac{(\lambda/\mu)^N}{N!}\ (k>N).$$

The three performance measures are:

• the delay probability $P_N\{W>0\}=\sum_{k\ge N}\pi_k$, the chance an arriving caller must wait;
• the abandonment probability $P_N\{Ab\}=\sum_{n\ge0}\pi_{N+n}\,\frac{(n+1)\theta}{N\mu+(n+1)\theta}=\frac{\theta}{\lambda}\sum_{k>N}(k-N)\pi_k$;
• the Erlang-B blocking probability $P_N\{Bl\}=E(\lambda/\mu,N)$ of the loss system $M/M/N/N$.

Write $\varphi,\Phi$ for the standard normal density and distribution function, $h(x)=\varphi(x)/[1-\Phi(x)]$ for its hazard rate, and

$$w(x,y)=\Bigl[1+\frac{h(-xy)}{y\,h(x)}\Bigr]^{-1}.$$

In the asymptotic statements the $N$-th system has $N$ agents, arrival rate $\lambda_N$, the fixed service rate $\mu$, patience rate $\theta_N$, and traffic intensity $\rho_N=\lambda_N/(N\mu)$.

Formalization targets

Goal: Theorem 4 (p. 219)

With $\theta_N\equiv\theta\in(0,\infty)$, the following are equivalent:

$$\lim_{N\to\infty}\sqrt N(1-\rho_N)=\beta\in\mathbb R,\qquad \lim_{N\to\infty}P_N\{W>0\}=\alpha\in(0,1),\qquad \lim_{N\to\infty}\sqrt N\,P_N\{Ab\}=\Delta\in(0,\infty),$$

in which case $\alpha=w(-\beta,\sqrt{\mu/\theta})$ and $\Delta=[\sqrt{\theta/\mu}\,h(\beta\sqrt{\mu/\theta})-\beta]\,\alpha$.

Milestones

1. The weak limit of the scaled stationary queue $(Q_N(\infty)-N)/\sqrt N$ (proof of Theorem 2*, Part 2, p. 225).
2. Lemma 1: $P_N\{W>0\}\to w(-\beta,\sqrt{\mu/\theta})$.
3. Three exact finite-$N$ identities for $P_N\{Ab\}$ (Riordan, Palm, and the final form in $P_N\{W>0\}$ and $P_N\{Bl\}$), p. 226.
4. The Erlang-B asymptotic $\sqrt N\,P_N\{Bl\}\to h(-\beta)$.
5. Lemma 2: $\sqrt N\,P_N\{Ab\}\to[\sqrt{\theta/\mu}\,h(\beta\sqrt{\mu/\theta})-\beta]\,w(-\beta,\sqrt{\mu/\theta})$.
6. The Remark on the extremes: $\beta=-\infty$ iff $\alpha=1$ iff $\Delta=\infty$, and $\beta=+\infty$ iff $\alpha=0$ iff $\Delta=0$.

Lemmas 1 and 2 and the stationary limit are stated with $\theta_N\to\theta\in(0,\infty)$, the paper's setting for them.

Significance

The result. Theorem 4 states that a single parameter, the scaled staffing gap $\beta$, determines both service measures in the limit, and that no other staffing sequence gives a nondegenerate delay probability or abandonment of order $1/\sqrt N$. Because $\beta$ may be negative, a call center with impatient callers can be staffed below the offered load and still keep abandonment at order $1/\sqrt N$, which is impossible without abandonment. The explicit $\alpha(\beta)$ and $\Delta(\beta)$ turn a target delay probability or abandonment level into a value of $\beta$, which is the staffing rule of Table 4 of the paper.

Formalizing it. The theorem is proved in the paper, partly by citation (Fleming et al. for the stationary limit, Jagerman for Erlang-B, Riordan and Palm for exact formulas); none of it has been machine-checked. A formalization supplies the cited inputs, the exact identities, and the converse directions, which the paper treats briefly. The $\theta=0$ counterpart, the Halfin–Whitt theorem, is posed separately on the platform (QueueingFundamentals.BirthDeath.halfin_whitt, HalfinWhitt81.Stationary.theorem_1) and is not a target here.

Difficulty

Every quantity is an explicit series, but its limit is a local-limit problem. The $N$ terms below $N$ form a Poisson-type sum whose normalized tail must be shown to converge to a Gaussian integral uniformly enough to pass through the normalization. The terms above $N$ form a product whose logarithm must be expanded to second order in $(k-N)/\sqrt N$. The two pieces enter with comparable weight, so neither can be dropped. Termwise convergence alone does not suffice: the series need tail bounds uniform in $N$. The converse directions require monotonicity of $\alpha$ and $\Delta$ in $\beta$ and the extreme regimes, not only the limits at finite $\beta$.

Formalization scope

• Model. All objects are explicit functions of $(N,\lambda,\mu,\theta)$ defined in ErlangA.Staffing.Model; the stationary law is the explicit normalized sequence, not an arbitrary steady state. The sequences are $\lambda,\theta:\mathbb N\to\mathbb R$ with $\lambda_N,\theta_N>0$ for $N\ge1$, and $\mu>0$ is fixed.
• Delay. $P_N\{W>0\}=\sum_{k\ge N}\pi_k$ (arrivals finding all agents busy). The paper's proof writes $P\{Q_N(\infty)>N\}$, which differs by $\pi_N\to0$.
• Abandonment. $P_N\{Ab\}$ is the arrival form above. It equals the balance form $(\theta/\lambda)\sum_{k>N}(k-N)\pi_k$ of equation (2).
• Limits. Every limit is Filter.Tendsto … atTop of the $N$-th quantity; $\pm\infty$ limits use atTop/atBot. The ranges $\beta\in\mathbb R$, $\alpha\in(0,1)$, $\Delta\in(0,\infty)$ are part of the goal: without them the equivalence is false (the delay probability also converges, to $1$, when $\sqrt N(1-\rho_N)\to-\infty$).
• Goal shape. A conjunction of List.TFAE of the three existential conditions with the "in which case" clause. Neither part may be dropped, and the equivalence may not be weakened to one direction or to $\beta>0$.
• Misprints. Two printed formulas are corrected: the Riordan term $(\lambda_N\theta)$ is read as $(\lambda_N/\theta)$, and a parenthesis is restored in the density $f$ on p. 225.
• Infrastructure. Needed: tail bounds for Poisson-type series, a local central limit estimate for the Poisson distribution near its mean, Gaussian tail asymptotics $h(t)\sim t$, and the incomplete gamma function. The Erlang-B asymptotic and the Gaussian facts are reusable for every many-server queue. Solutions to any milestone, and to auxiliary lemmas such as summability of the stationary weights, are welcome.

Selected references

• O. Garnett, A. Mandelbaum, M. Reiman, Designing a call center with impatient customers, Manufacturing & Service Operations Management 4(3), 208–227, 2002. https://doi.org/10.1287/msom.4.3.208.7753
• S. Halfin, W. Whitt, Heavy-traffic limits for queues with many exponential servers, Operations Research 29(3), 567–587, 1981. https://doi.org/10.1287/opre.29.3.567
• D. L. Jagerman, Some properties of the Erlang loss function, Bell System Technical Journal 53, 1974 (as cited in the paper's reference list).
• P. J. Fleming, A. Stolyar, B. Simon, Heavy traffic limit for a mobile phone system loss model, Proc. 2nd Int. Conf. Telecommunication Systems, Modeling and Analysis, Nashville, 158–176, 1994.

Motivation

Many decision problems in operations research have a finite set of alternatives and an objective that is an expectation: choose a configuration, a route or a stocking level $x$ from a finite set $\mathcal S$ to minimize the expected cost $\mathbb E\,G(x, W)$, where $W$ is random. When the expectation cannot be computed in closed form, the sample average approximation (SAA) method draws an i.i.d. sample $W^1, \dots, W^N$, replaces the expectation by the sample average, and solves the resulting deterministic problem. Its optimal value $\hat v_N$ is used in practice as an estimate of the true optimal value $v^*$, and statistical bounds on the quality of a candidate solution are built from replications of $\hat v_N$ (Mak, Morton and Wood, 1999). Confidence statements of that kind need the distribution of $\hat v_N$, at least asymptotically.

Kleywegt and Shapiro (preprint, two-author version; published with Homem-de-Mello in SIAM J. Optim. 12 (2002)) gave the asymptotic distribution of $\hat v_N$ for discrete problems in their Proposition 2.3. The analogous result for problems with continuous decision variables is due to Shapiro (1991); the discrete statement needs neither compactness nor Lipschitz continuity in $x$, and requires finite variances only at the optimal solutions.

Setting

Let $\mathcal S$ be a nonempty finite set, $\mathcal W$ a measurable space, and $G : \mathcal S \times \mathcal W \to \mathbb R$ a function such that for every $x \in \mathcal S$ the map $G(x, \cdot)$ is measurable and $\mathbb E\,|G(x, W)| < \infty$, where $W$ is a random element of $\mathcal W$. The true problem (1.1) is

$$v^* = \min_{x \in \mathcal S} g(x), \qquad g(x) = \mathbb E\, G(x, W),$$

and its set of optimal solutions is $\mathcal S^* = \{x \in \mathcal S : g(x) = v^*\}$. Given an i.i.d. sample $W^1, W^2, \dots$ of $W$, the sample average function and the SAA problem (2.1) are

$$\hat g_N(x) = \frac1N \sum_{n=1}^N G(x, W^n), \qquad \hat v_N = \min_{x \in \mathcal S} \hat g_N(x),$$

with optimal set $\hat{\mathcal S}_N$. Because $\mathcal S$ is finite, all minima are attained. For $x \in \mathcal S$ the variance $\sigma^2(x) = \operatorname{Var}\{G(x, W)\}$ of (2.7) may or may not exist; when the variances exist on a subset $T \subseteq \mathcal S$, the autocovariance of $G$ on $T$ is the matrix $\big(\operatorname{Cov}(G(x, W), G(x', W))\big)_{x, x' \in T}$. Convergence in distribution is written $\Rightarrow$.

Formalization targets

Goal: Proposition 2.3, (2.8)

If the variances $\sigma^2(x)$ exist for every $x \in \mathcal S^*$, then

$$\sqrt N\,(\hat v_N - v^*) \;\Rightarrow\; \min_{x \in \mathcal S^*} Y(x),$$

where $(Y(x))_{x \in \mathcal S^*}$ is a Gaussian vector with zero mean and $\operatorname{Cov}(Y(x), Y(x')) = \operatorname{Cov}(G(x, W), G(x', W))$. The statement is general in the number of optimal solutions: when $\mathcal S^*$ has more than one point, the limit is the minimum of a correlated Gaussian vector, which is not normal and has negative mean.

Milestones (in the order of the paper's argument)

1. Almost surely, $\hat v_N - \min_{x \in \mathcal S^*} \hat g_N(x) = 0$ for all $N$ large enough.
2. (2.6): $\sqrt N\,[\hat v_N - \min_{x \in \mathcal S^*} \hat g_N(x)] \to 0$ almost surely and in probability.
3. $\sqrt N\,[\min_{x \in \mathcal S^*} \hat g_N(x) - v^*] = \min_{x \in \mathcal S^*} \{\sqrt N\,[\hat g_N(x) - g(x)]\}$.
4. The joint central limit theorem: for $T \subseteq \mathcal S$ with finite variances,

$$\Big(\sqrt N\,[\hat g_N(x) - g(x)]\Big)_{x \in T} \;\Rightarrow\; (Y(x))_{x \in T}.$$

Companions

• (2.9): if $\mathcal S^* = \{x^*\}$, then $\sqrt N\,(\hat v_N - v^*) \Rightarrow N(0, \sigma^2(x^*))$.
• Negative bias (pp. 4–5): for every $N \ge 1$, $\mathbb E[\hat v_N] \le v^*$.

Significance

The result. Proposition 2.3 identifies the first-order fluctuations of the SAA optimal value. It shows that $\hat v_N$ is $\sqrt N$-consistent, that it is asymptotically normal when the true optimum is unique, and that its asymptotic bias is of order $N^{-1/2}$ and equal to $\mathbb E[\min_{x \in \mathcal S^*} Y(x)] < 0$ when several solutions are optimal and their costs are not perfectly correlated. Confidence intervals and stopping rules for SAA, and the statistical lower bounds used to assess candidate solutions, rest on this description.

Formalizing it. The result is classical and proved on paper; it has not been machine-checked. A formal proof requires a multivariate central limit theorem, which this Mathlib does not contain (only the one-dimensional CLT is available), together with the transfer of convergence in distribution through the continuous map $y \mapsto \min_x y(x)$ and through perturbations that vanish in probability. These components are general probability results whose use extends well beyond stochastic programming. The negative-bias inequality and the singleton case (2.9) are included as companion targets.

Difficulty

The argument on the page is short, and each step names a standard theorem. The formal difficulty is that two of those theorems are not available in the needed form. The central limit theorem is used for a random vector indexed by $\mathcal S^*$, with an arbitrary, possibly singular, covariance matrix; a coordinatewise application of the one-dimensional CLT gives only the marginal limits, which do not determine the law of the minimum. The step from (2.6) to (2.8) combines convergence in distribution with a term that tends to zero in probability, which requires a Slutsky-type lemma on the real line. The tempting simplification of assuming a unique optimal solution removes exactly the case in which the limit is not normal, and is therefore not the goal.

Formalization scope

The feasible set is a Finset of an arbitrary type X, assumed nonempty; minima over it and over $\mathcal S^*$ are Finset.inf', so they are attained. $\mathcal S^*$ is the Finset of minimizers of $g$. The sample is a sequence W : ℕ → Ω → 𝒲 of measurable maps on a probability space, mutually independent (iIndepFun) and identically distributed with W 0; the paper's $W^n$ is W (n - 1), and $g(x)$ is the expectation of $G(x, W^1)$. Measurability of $G(x, \cdot)$ is Borel measurability. "The variance exists" is square integrability (MemLp … 2), assumed in the goal only on $\mathcal S^*$, exactly as in the Proposition. Convergence in distribution is Mathlib's TendstoInDistribution; the limit is any random vector $Y$, on any probability space, whose law is multivariateGaussian 0 Σ on the Euclidean space $\mathbb R^{\mathcal S^*}$, with $\Sigma$ the autocovariance matrix of $G$ on $\mathcal S^*$. That matrix is positive semidefinite under the variance hypothesis, so the Gaussian law is genuine (Mathlib's multivariateGaussian degenerates to a point mass only for a matrix that is not positive semidefinite). The milestone CLT is stated jointly, for every subfamily $T \subseteq \mathcal S$ with finite variances; the page states it for $T = \mathcal S$ in terms of marginals plus covariances.

A trivializing formalization is ruled out: the limit is indexed by $\mathcal S^*$ and not by $\mathcal S$, its covariance is the autocovariance of $G$ and not a free matrix, square integrability is assumed wherever a covariance is used, and the negative-bias statement asserts integrability of $\hat v_N$ so that the Bochner integral is not its junk value $0$.

A complete development needs the multivariate central limit theorem for i.i.d. square-integrable vectors (for instance through characteristic functions or the Cramér–Wold device), the continuous mapping theorem, a Slutsky lemma, and the strong law of large numbers applied to finitely many functions. The multivariate CLT and the Slutsky lemma are reusable beyond this mission. Contributions to any of these, and proofs of any item, are welcome.

Selected references

• A. J. Kleywegt and A. Shapiro, The sample average approximation method for stochastic discrete optimization, preprint (two-author version), https://edoc.hu-berlin.de/server/api/core/bitstreams/2dcf2d5b-5fb9-44d1-b36a-5e378c803717/content. Published with T. Homem-de-Mello in SIAM J. Optim. 12(2):479–502, 2002. https://doi.org/10.1137/S1052623499363220
• A. Shapiro, Asymptotic analysis of stochastic programs, Annals of Operations Research 30:169–186, 1991. https://doi.org/10.1007/BF02204815
• W. K. Mak, D. P. Morton and R. K. Wood, Monte Carlo bounding techniques for determining solution quality in stochastic programs, Operations Research Letters 24:47–56, 1999. https://doi.org/10.1016/S0167-6377(98)00054-6
• A. Shapiro, D. Dentcheva and A. Ruszczyński, Lectures on Stochastic Programming: Modeling and Theory, SIAM, 2009. https://doi.org/10.1137/1.9780898718751