Operations Research

Motivation

Fast-fashion retailers such as Zara replenish their stores several times a week from a central warehouse, and the warehouse holds a limited stock of each garment that must be split across hundreds of stores. Caro and Gallien built the store-level sales model behind the shipment-allocation system they developed with Zara and tested in a live pilot (Caro and Gallien, Inventory Management of a Fast-Fashion Retail Network, working paper, August 2, 2007; published in Operations Research 58(2), 2010). One display policy shapes the model. A garment (a reference) comes in several sizes, a few of which are designated major sizes. As soon as one major size runs out in a store, the whole reference is removed from the shop floor, so the remaining stock of the other sizes stops selling.

The allocation optimization only works if expected store sales have the right shape. They must grow with inventory, show decreasing marginal returns, and show complementarity across sizes. Proposition 1 of the paper establishes these properties for the model. The paper also cites the analogous result of Lu and Song (2003) for assemble-to-order systems as the template for its proof.

Setting

Fix a finite set of sizes $\mathcal S$, split into the major sizes $\mathcal S^+$ and the minor sizes $\mathcal S^- = \mathcal S \setminus \mathcal S^+$. Time $t \ge 0$ is measured from the last replenishment, and the next replenishment comes at time $T > 0$.

Demand is random. On a probability space, each size $s$ has a counting process $N_s(t)$, the number of sale opportunities for $s$ in $[0, t]$. It is a Poisson process with rate $\lambda_s > 0$:

• $N_s(0) = 0$;
• its paths are non-decreasing and right-continuous with jumps of size one;
• the increment $N_s(t) - N_s(u)$ is Poisson with mean $\lambda_s (t-u)$;
• increments over disjoint intervals are independent.

The processes of different sizes are mutually independent. In Lean this family is IsPoissonFamily lam N P.

An inventory vector $q = (q_s) \in \mathbb N^{\mathcal S}$ gives the stock of each size right after replenishment. The virtual stockout time of size $s$ is

$$\tau_s(q_s) = \inf\{t \ge 0 : N_s(t) = q_s\},$$

the time at which size $s$ would run out if it stayed on display. For a set of sizes $\mathcal A$, write $\tau_{\mathcal A} = \min_{s \in \mathcal A} \tau_s(q_s)$ and $a \wedge b = \min(a, b)$. Under the major-size policy, the random sales over one period are

$$G(q) = \sum_{s \in \mathcal S^+} N_s(\tau_{\mathcal S^+} \wedge T) + \sum_{s \in \mathcal S^-} N_s(\tau_{\mathcal S^+ \cup \{s\}} \wedge T). \tag{1}$$

The expected sales function is $g(q) = \mathbb E[G(q)]$ (expectedSales). The proof also uses $h^{\mathcal A}(q) = \mathbb E[\tau_{\mathcal A} \wedge T]$ (hA) and the marginal difference $\Delta_s f(q) = f(q + e_s) - f(q)$ (delta), where $e_s$ is the $s$-th unit vector.

Formalization targets

Goal: Proposition 1

For every $q \in \mathbb N^{\mathcal S}$ and all sizes $s \ne s'$:

$$g \text{ is non-decreasing}, \qquad \Delta_s g(q + e_s) \le \Delta_s g(q), \qquad \Delta_s g(q) \le \Delta_s g(q + e_{s'}),$$

and $g$ is supermodular on the lattice $\mathbb N^{\mathcal S}$:

$$g(q) + g(q') \le g(q \vee q') + g(q \wedge q').$$

The statement is uniform in the rates, the horizon and the choice of major sizes, so it covers every instance of the model.

Milestones (Appendix §5.1 and §3.1.3, in the order the proof uses them)

1. On every sample path, $q \mapsto G(q)$ is non-decreasing.
2. Identity (2): $g(q) = \lambda_{\mathcal S^+} h^{\mathcal S^+}(q) + \sum_{s \in \mathcal S^-} \lambda_s h^{\mathcal S^+ \cup \{s\}}(q)$, where $\lambda_{\mathcal S^+} = \sum_{s \in \mathcal S^+} \lambda_s$.
3. For $s \in \mathcal A$, $\Delta_s h^{\mathcal A}(q)$ equals an integral over $[0, T]$ of Poisson probabilities, and that integral equals $\frac{1}{\lambda_s}\mathbb P(\tau_s(q_s+1) \le \tau_{\mathcal A \setminus \{s\}} \wedge T)$.
4. $\Delta_s h^{\mathcal A}$ is non-increasing in $q_s$ and non-decreasing in $q_{s'}$ for $s' \ne s$.
5. On every sample path, $q \mapsto \tau_{\mathcal A} \wedge T$ is supermodular.
6. $h^{\mathcal A}$ is supermodular.

Significance

The result. Proposition 1 justifies treating each store's expected sales as a concave, complementary function of the size profile. The paper's allocation MIP is built on that structure. It approximates $g$ by a lower envelope of tangents and embeds the approximation in a mixed integer program that allocates warehouse stock across stores. The proposition explains why sending a unit of a major size can raise the value of the minor sizes already in a store. It also explains why the marginal value of a size falls as more of it is shipped. Without these properties a greedy or envelope-based allocation would have no structural footing.

Formalizing it. The result is proved on paper, with one error in the printed proof. The last line of the Appendix display gives $\Delta_s h^{\mathcal A}(q)$ as a probability, but the quantity is a time. The correct value carries a factor $1/\lambda_s$ and has $q_s + 1$ in place of $q_s$. This mission states the corrected identity. As far as is known, no machine-checked proof of any part of the result exists. A complete development would give a verified link from a continuous-time stochastic model to the discrete convexity properties that inventory optimization relies on. It would also produce reusable facts about Poisson processes and their hitting times along the way.

Difficulty

The monotonicity of $g$ is pathwise and elementary. The other properties are not pathwise properties of $G$. Under the major-size policy, a realisation of $G$ need not have decreasing increments in $q_s$. The concavity and the complementarity appear only after taking expectations and rewriting $g$ through the stopping times, using identity (2). That identity is an optional-sampling statement for the compensated Poisson processes at the bounded random times $\tau_{\mathcal A} \wedge T$. These times depend on several independent processes at once, so the relevant filtration is the joint one. Mathlib has optional sampling only in discrete time.

The marginal-difference formula needs three ingredients:

• the tail-integral representation of $\mathbb E[\tau_{\mathcal A} \wedge T]$;
• the factorisation of $\mathbb P(\tau_{\mathcal A} > t)$ through independence across sizes;
• the Erlang law of the hitting times $\tau_s(k)$.

The first idea, proving concavity of $g$ directly from (1) path by path, fails for the reason above.

Formalization scope

All objects live in the namespace FastFashion.Structure:

• Sizes: a Fintype with decidable equality. The major sizes are a Finset Sp and the minor sizes are its complement. No nonemptiness assumption is placed on the sizes or on $\mathcal S^+$; $\mathcal S^+ = \emptyset$ is allowed, as in the paper.
• Inventories: vectors $q : \mathcal S \to \mathbb N$ with the componentwise order and lattice operations. $e_s$ is Pi.single s 1.
• Stockout times: $\tau_s(q_s)$ takes values in WithTop ℝ. It is $+\infty$ on paths that never reach $q_s$, so it is never Lean's junk value $0$, and $\tau_\emptyset = +\infty$. The truncation $\tau_{\mathcal A} \wedge T$ (stopMin) is a real number in $[0, T]$.
• Expectations and probabilities: expectations are Bochner integrals. They need no integrability hypothesis, because $0 \le G(q) \le \sum_s N_s(T)$ and $0 \le \tau_{\mathcal A} \wedge T \le T$. Probabilities are real numbers.
• Supermodularity: Topkis's lattice definition, the published platform definition Supermodularity.Monotonicity.SupermodularOn, applied on the whole lattice.
• Poisson family: the processes are independent as whole processes. The σ-algebras generated by all times are independent, not only the one-dimensional marginals.

Two formalizations would trivialize the mission, and neither is used. $g$ is defined as the expectation of $G$ from (1), never by formula (2), and $h^{\mathcal A}$ is never defined by a product-of-probabilities integral. Proposition 1 states both readings of "supermodular": the per-coordinate marginal-difference inequalities and the lattice inequality.

Deviations from the printed text:

• The last line of the Appendix display is replaced by $\frac{1}{\lambda_s}\mathbb P(\tau_s(q_s+1) \le \tau_{\mathcal A\setminus\{s\}} \wedge T)$.
• "non-decreasing in $x_s$" in Proposition 1 is read as $q_s$.

Infrastructure a complete development needs:

• hitting times of counting processes and their Erlang laws;
• a continuous-time optional sampling theorem for the compensated Poisson process, or a direct computation from independent increments;
• the tail formula $\mathbb E[X] = \int_0^\infty \mathbb P(X > t)\,dt$ for bounded non-negative $X$;
• lattice facts about minima of single-variable increasing functions.

The Poisson-process and hitting-time lemmas are reusable well beyond this mission and are welcome as separate contributions. The companion mission on the paper's tangent-envelope approximation uses the same model layer.

Selected references

• F. Caro, J. Gallien, Inventory Management of a Fast-Fashion Retail Network, working paper, August 2, 2007; published in Operations Research 58(2):257–273, 2010. https://doi.org/10.1287/opre.1090.0698
• D. M. Topkis, Supermodularity and Complementarity, Princeton University Press, 1998. https://doi.org/10.1515/9781400822539
• Y. Lu, J.-S. Song, Order-Based Cost Optimization in Assemble-to-Order Systems, Operations Research 53(1):151–169, 2005 (working paper 2003). https://doi.org/10.1287/opre.1040.0146
• I. Karatzas, S. E. Shreve, Brownian Motion and Stochastic Calculus, 2nd ed., Springer, 1991. https://doi.org/10.1007/978-1-4612-0949-2

Why stability of two-stage programs matters

A linear two-stage stochastic program with fixed recourse chooses a here-and-now decision $x$ before a random vector $\xi$ is observed, and pays a recourse cost afterwards. It is the basic model of stochastic programming, used in capacity planning, energy and supply-chain models. In practice the distribution $\mu$ of $\xi$ is never known exactly: it is estimated from data, replaced by a discrete scenario approximation, or perturbed for robustness. The question is how much the optimal value and the optimal decisions can move when $\mu$ is replaced by a nearby $\nu$, and in which distance between probability measures "nearby" should be measured.

Rachev and Römisch (Math. Oper. Res. 27(4), 2002, doi:10.1287/moor.27.4.792.304; the mission cites the authors' preprint from edoc.hu-berlin.de) answer this with the method of probability metrics: they derive from the structure of the integrand a canonical distance on measures for each class of models. For linear two-stage programs the canonical distance is the Fortet–Mourier metric of order 2.

Timeline. Robinson and Wets (1987) studied qualitative stability of two-stage programs with respect to weak convergence of measures. Walkup and Wets (1969) had earlier established the piecewise-bilinear structure of the second-stage value function used here. Römisch and Schultz (1991) proved a quantitative stability result for two-stage programs with complete recourse in the Wasserstein metric $W_2$. The present paper (published 2002) replaces complete recourse by relatively complete recourse plus dual feasibility, and $W_2$ by the metric $\zeta_2$, which it bounds by a multiple of $W_2$ (p. 13).

Setting

Let $X\subseteq\mathbb R^m$ be a nonempty polyhedron and $\Xi\subseteq\mathbb R^s$ a polyhedron. Fix $c\in\mathbb R^m$, an $(r,\overline m)$-matrix $W$, and data $q(\xi)\in\mathbb R^{\overline m}$, $h(\xi)\in\mathbb R^r$ and an $(r,m)$-matrix $T(\xi)$ depending affine linearly on $\xi$. The second-stage value is

$$\Phi(u,t)=\inf\{uy:\ Wy=t,\ y\ge0\},$$

with $\operatorname{pos}W=\{Wy:y\ge0\}$ and $D=\{u:\{z:W'z\le u\}\ne\emptyset\}$. The integrand is $f_0(\xi,x)=cx+\Phi(q(\xi),h(\xi)-T(\xi)x)$ when $h(\xi)-T(\xi)x\in\operatorname{pos}W$ and $q(\xi)\in D$, and $+\infty$ otherwise. For a Borel probability measure $\nu$ on $\Xi$ the problem is

$$\min\Big\{\int_\Xi f_0(\xi,x)\,\nu(d\xi):x\in X\Big\},$$

with optimal value $v(\nu)$ and solution set $S(\nu)$. Two assumptions are made: (A1) $h(\xi)-T(\xi)x\in\operatorname{pos}W$ and $q(\xi)\in D$ for all $(\xi,x)\in\Xi\times X$; (A2) $\mu$ has a finite second moment.

The Fortet–Mourier metric $\zeta_2$ on $\mathcal P_2(\Xi)$, the probability measures on $\Xi$ with finite second moment, is

$$\zeta_2(\mu,\nu)=\sup_{f\in\mathcal F_2(\Xi)}\Big|\int_\Xi f\,d\mu-\int_\Xi f\,d\nu\Big|,\quad \mathcal F_2(\Xi)=\{f:|f(\xi)-f(\tilde\xi)|\le\max\{1,\|\xi\|,\|\tilde\xi\|\}\|\xi-\tilde\xi\|\}.$$

For an open bounded $\mathcal U\supseteq S(\mu)$, the growth function is $\psi(\tau)=\inf\{\int_\Xi f_0(\xi,x)\,\mu(d\xi)-v(\mu): d(x,S(\mu))\ge\tau,\ x\in X\cap\operatorname{cl}\mathcal U\}$, and $\Psi(\eta)=\eta+\psi^{-1}(2\eta)$ with $\psi^{-1}(t)=\sup\{\tau\ge0:\psi(\tau)\le t\}$.

Formalization targets

Goal: Theorem 3.3 (p. 13)

Under (A1), (A2), $S(\mu)\ne\emptyset$ and $\mathcal U$ an open bounded neighbourhood of $S(\mu)$, there are $L>0$ and $\delta>0$ such that for every $\nu\in\mathcal P_2(\Xi)$ with $\zeta_2(\mu,\nu)<\delta$:

$$|v(\mu)-v(\nu)|\le L\,\zeta_2(\mu,\nu),\qquad \emptyset\ne S(\nu)\subseteq S(\mu)+\Psi(L\,\zeta_2(\mu,\nu))\,\mathbb B.$$

The constants are existential, so the goal does not depend on any particular estimate of them.

Milestones

1. Lemma 3.1 (p. 11): $\Phi$ is finite and continuous on the polyhedral cone $D\times\operatorname{pos}W$, piecewise of the form $C_ju\cdot t$ on finitely many polyhedral cones with disjoint interiors, convex in $t$ and concave in $u$.
2. Proposition 3.2 (p. 11): $f_0$ is a normal convex integrand, and for $\|x\|\le r$ it satisfies $|f_0(\xi,x)-f_0(\tilde\xi,x)|\le Lr\max\{1,\|\xi\|,\|\tilde\xi\|\}\|\xi-\tilde\xi\|$, a Lipschitz bound in $x$ with factor $\hat L\max\{1,\|\xi\|^2\}$, and $|f_0(\xi,x)|\le Kr\max\{1,\|\xi\|^2\}$.
3. $\mathcal P_{\mathcal F_{\mathcal U}}\supseteq\mathcal P_2(\Xi)$ (p. 13): every measure with finite second moment satisfies the integrability conditions of the general theory.
4. Corollary 2.8 (p. 10): the general Lipschitz stability theorem in the metric $\zeta_g$ for convex models whose integrands lie, up to a factor $L$, in the class $\mathcal F_g$.

Significance

The theorem says that, under the two natural well-posedness conditions of linear recourse models, optimal values are Lipschitz continuous and solution sets move by at most $\Psi(L\zeta_2)$, where $\zeta_2$ only sees test functions with the local Lipschitz growth of the integrand. Since $\zeta_2\le(1+\int\|\xi\|^2d\mu+\int\|\xi\|^2d\nu)^{1/2}W_2$, the result contains the earlier $W_2$ stability theorem.

The results are proved in the paper, except Lemma 3.1, which it cites from Walkup and Wets (1969). None of them is formalized in Mathlib or on Prove2Me. A formalization would provide a checked model of the two-stage recourse function, the Fortet–Mourier metrics, and an end-to-end quantitative stability theorem for stochastic programs, none of which exist in Mathlib.

Difficulty

The obvious argument bounds $|v(\mu)-v(\nu)|$ by $\sup_x|\int f_0(\cdot,x)\,d(\mu-\nu)|$ and then by $L\zeta_2$. Two steps fail. First, the supremum ranges over all of $X$, which may be unbounded, while $f_0(\cdot,x)/L$ lies in $\mathcal F_2$ only for $x$ in a bounded set; the argument must localize to $X\cap\operatorname{cl}\mathcal U$ and then show, using convexity, that for $\nu$ close to $\mu$ the localized and global problems coincide. Second, the Lipschitz estimate in $\xi$ for $f_0$ requires the piecewise-bilinear structure of $\Phi$ (Lemma 3.1) and a chaining argument across the polyhedral pieces of $\Xi$; the integrand is not differentiable.

Formalization scope

• $\mathbb R^n$ is EuclideanSpace ℝ (Fin n). The paper never fixes its norm and all its constants are existential, so the Euclidean norm is a faithful instance. Second-stage vectors $u,t,y,z$ are coordinate vectors Fin k → ℝ; index sets are 0-based.
• A polyhedron is a finite intersection of closed half-spaces. The standing assumptions are: $X$ a nonempty polyhedron, $\Xi$ a polyhedron, $q,h,T$ affine maps.
• A measure in $\mathcal P(\Xi)$ is a Borel probability measure on $\mathbb R^s$ with $\nu(\Xi^c)=0$; $\int_\Xi f_0(\xi,x)\,\nu(d\xi)$ is the extended-real integral of the published DupacovaWets.Consistency.expect over $\nu$ restricted to $\Xi$.
• $\Phi$, $f_0$, $v$ and $\psi$ are extended-real infima ($+\infty$ over an empty set); "min" in $\psi$ is read as an infimum. $\zeta_2$, $\zeta_g$, $\psi^{-1}$ and $\Psi$ take values in $[0,\infty]$, and the distance to a set is the extended distance ($+\infty$ to $\emptyset$).
• To rule out trivializations: every value bound states that $v(\mu)$, $v(\nu)$ are finite before comparing them, so the extended-real arithmetic $\infty-\infty$ cannot make it vacuous; in $\zeta_2$ a test function that is not integrable for both measures contributes $+\infty$, never a junk $0$.
• Disclosed readings: Proposition 3.2 leaves its radius $r$ unquantified; it is stated with the constants first and then every $r\ge1$. Lemma 3.1's "$\operatorname{pos}W\times D$" is read as $D\times\operatorname{pos}W$, in the order of $\Phi$'s arguments. Corollary 2.8 adds $\mu,\nu\in\mathcal P_{\mathcal F_{\mathcal U}}(\Xi)$, an assumption of Theorem 2.2 that its proof invokes; Theorem 3.3 needs no such addition, because milestone 3 supplies it.
• Needed infrastructure: LP duality and the Walkup–Wets decomposition of $\Phi$; normal integrands and measurability of their infima; integrability under moment conditions; the general stability theorems (Theorems 2.2–2.3 of the paper) for $d=0$. The definitions of $\zeta_p$, $\mathcal F_g$ and normal integrands are reusable well beyond this mission. Contributions toward any milestone, and proofs of the auxiliary facts (for example, that every $f\in\mathcal F_2(\Xi)$ is integrable for $\nu\in\mathcal P_2(\Xi)$), are welcome.

Selected references

• S. T. Rachev, W. Römisch, Quantitative stability in stochastic programming: The method of probability metrics, Mathematics of Operations Research 27(4) (2002) 792–818. https://doi.org/10.1287/moor.27.4.792.304 (cited here from the authors' preprint, edoc.hu-berlin.de)
• D. W. Walkup, R. J-B Wets, Lifting projections of convex polyhedra, Pacific Journal of Mathematics 28 (1969) 465–475 (the paper's [55]).
• S. M. Robinson, R. J-B Wets, Stability in two-stage stochastic programming, SIAM Journal on Control and Optimization 25 (1987) 1409–1416 (the paper's [37]).
• W. Römisch, R. Schultz, Stability analysis for stochastic programs, Annals of Operations Research 30 (1991) 241–266 (the paper's [40]).
• R. T. Rockafellar, R. J-B Wets, Variational Analysis, Springer (the paper's [38]; normal integrands, Chapter 14).
• S. T. Rachev, Probability Metrics and the Stability of Stochastic Models, Wiley, Chichester 1991 (the paper's [32]).

Motivation

A Markov decision process is specified by numbers: transition probabilities and rewards. In practice these numbers are estimated from data, elicited from experts, or produced by aggregating the states of a larger model, and in each case what is actually known is a range for each parameter rather than its value. Givan, Leach and Dean (Artificial Intelligence, 2000) introduced bounded-parameter MDPs (BMDPs) to reason about a whole family of MDPs at once: every transition probability is only known to lie in a closed interval. Their motivation was state-space aggregation, where grouping states with similar but not identical dynamics yields exactly such interval bounds, and the same object later reappeared as the interval or "box" uncertainty set of robust MDPs (Nilim and El Ghaoui, 2005; Iyengar, 2005) and of optimistic exploration in reinforcement learning (Strehl and Littman, 2008).

Over such a family, the value of a policy is no longer a number but an interval, and "optimal" needs a definition. This mission formalizes the paper's answer: two total orders on intervals, the corresponding optimal policies, and the Bellman-like equations their interval value functions satisfy.

Setting

An exact MDP $M=\langle Q,A,F,R\rangle$ has a finite set $Q$ of states, a finite nonempty set $A$ of actions, transition probabilities $F_{pq}(\alpha)\ge 0$ with $\sum_{q}F_{pq}(\alpha)=1$, a reward $R(q)\in\mathbb R$ for each state and a discount rate $0\le\gamma<1$. A policy is a map $\pi:Q\to A$; $\Pi$ is the set of all policies. The value function $V_{M,\pi}(p)$ is the expected discounted sum of rewards $\sum_{t\ge 0}\gamma^t\,\mathbb E[R(X_t)\mid X_0=p]$ along the Markov chain with transitions $F_{pq}(\pi(p))$. The operators

$$VI_{M,\pi}(v)(p)=R(p)+\gamma\sum_{q}F_{pq}(\pi(p))\,v(q),\qquad VI_{M,\alpha}(v)(p)=R(p)+\gamma\sum_{q}F_{pq}(\alpha)\,v(q)$$

act on value functions $v:Q\to\mathbb R$. $V_1\le_{\mathrm{dom}}V_2$ is the statewise order.

A BMDP $M_\updownarrow$ gives for each $p,q,\alpha$ an interval $[F_\downarrow{}_{pq}(\alpha),F_\uparrow{}_{pq}(\alpha)]\subseteq[0,1]$ with $\sum_q F_\downarrow{}_{pq}(\alpha)\le 1\le\sum_q F_\uparrow{}_{pq}(\alpha)$, a reward $R$ and a discount rate $\gamma$. Its members $M\in M_\updownarrow$ are the exact MDPs with that $R$ and $\gamma$ and with $F^M_{pq}(\alpha)\in[F_\downarrow{}_{pq}(\alpha),F_\uparrow{}_{pq}(\alpha)]$. The interval value of a policy is

$$V_\updownarrow{}_\pi(q)=\Big[\min_{M\in M_\updownarrow}V_{M,\pi}(q),\ \max_{M\in M_\updownarrow}V_{M,\pi}(q)\Big]=[V_\downarrow{}_\pi(q),V_\uparrow{}_\pi(q)].$$

Intervals are compared by the optimistic order, which compares upper bounds first and breaks ties by lower bounds, and the pessimistic order, which compares lower bounds first:

$$[l_1,u_1]\le_{\mathrm{opt}}[l_2,u_2]\iff u_1<u_2\ \vee\ (u_1=u_2\wedge l_1\le l_2),\qquad [l_1,u_1]\le_{\mathrm{pes}}[l_2,u_2]\iff l_1<l_2\ \vee\ (l_1=l_2\wedge u_1\le u_2).$$

A policy $\pi_{\mathrm{opt}}$ is optimistically optimal if $V_\updownarrow{}_{\pi_{\mathrm{opt}}}(q)\ge_{\mathrm{opt}}V_\updownarrow{}_\pi(q)$ for all $\pi$ and $q$; pessimistically optimal is defined with $\ge_{\mathrm{pes}}$. The optimal interval value functions are $V_\updownarrow{}_{\mathrm{opt}}=V_\updownarrow{}_{\pi_{\mathrm{opt}}}$ and $V_\updownarrow{}_{\mathrm{pes}}=V_\updownarrow{}_{\pi_{\mathrm{pes}}}$.

Formalization targets

Goal: Theorem 9, equations (25) and (26)

At every state $p$,

$$V_\updownarrow{}_{\mathrm{opt}}(p)=\max_{\alpha\in A,\ \le_{\mathrm{opt}}}\Big[\min_{M\in M_\updownarrow}VI_{M,\alpha}(V_\downarrow{}_{\mathrm{opt}})(p),\ \max_{M\in M_\updownarrow}VI_{M,\alpha}(V_\uparrow{}_{\mathrm{opt}})(p)\Big],$$

and the same with $\mathrm{pes}$ in place of $\mathrm{opt}$. The maximum over actions is for the total order on intervals: attained by some action and an upper bound for every action.

Milestones

In the order the proof uses them: Theorem 3 (second sentence: $V_{M,\pi}$ is the unique fixed point of $VI_{M,\pi}$); Theorem 6 (comparison: $u\le_{\mathrm{dom}}VI_{M,\pi}(u)$ implies $u\le_{\mathrm{dom}}V_{M,\pi}$, and its three variants); Lemma 1 (every member's values and backups are bracketed by finitely many order-maximizing MDPs); Lemma 2 (statewise composition of member MDPs dominates, or is dominated by, both components); Theorem 7 ($\pi$-maximizing and $\pi$-minimizing MDPs exist among the order-maximizing ones); Corollary 1 ($V_\downarrow{}_\pi$ and $V_\uparrow{}_\pi$ are a minimum and a maximum for $\le_{\mathrm{dom}}$); Lemma 3 (statewise composition of policies improves the upper, resp. lower, bounds); Theorem 8 (optimistically and pessimistically optimal policies exist).

Significance

Theorem 9 is the interval analogue of the Bellman optimality equation. It is what makes the optimal interval value functions computable: the paper's interval value iteration algorithms $IVI_\updownarrow{}_{\mathrm{opt}}$ and $IVI_\updownarrow{}_{\mathrm{pes}}$ (Section 5) iterate exactly the right-hand sides of (25) and (26). The lower half of (26) is the max-min equation of robust dynamic programming with rectangular box uncertainty, and the upper half of (25) is the max-max equation behind optimistic planning; Theorem 9 states both, together with the tie-breaking second component, in one model. Theorems 7 and 8 are of independent use: a single member MDP is worst (or best) for a policy at all states at once, and the partial orders on interval value functions nevertheless have maxima over policies.

None of these results has a machine-checked proof. The platform has robust Bellman equations in other models (costs, general rectangular sets, a sup over nature) but nothing on BMDPs, interval value functions, or the orders $\le_{\mathrm{opt}}$ and $\le_{\mathrm{pes}}$. The mission produces the paper's optimality theory end to end, from the comparison principle for exact MDPs to the two Bellman-like equations, on a model shared with the companion mission on interval policy evaluation.

Difficulty

The orders $\le_{\mathrm{opt}}$ and $\le_{\mathrm{pes}}$ are total on intervals but only partial on interval value functions, and they are lexicographic. The obvious construction of an optimal policy, switching statewise to whichever of two policies has the better interval, fails: the composed policy need not dominate either component in $\le_{\mathrm{opt}}$, because lower bounds can get worse at states where the upper bounds do not change (the paper says so explicitly, eq. (20)). Existence therefore needs a two-stage argument. A second obstacle is that $V_\downarrow{}_\pi$ and $V_\uparrow{}_\pi$ are defined statewise, as minima over an uncountable family; that one member attains them at every state simultaneously is a theorem (Theorem 7), and the Bellman-like equations rest on it. Finally, the inner minima in (25) act on $V_\downarrow{}_{\mathrm{opt}}$, the lower bound of an optimistically chosen policy, which is not the componentwise best lower bound over policies; the equation couples the two components through the tie-break.

Formalization scope

$Q$ and $A$ are finite types, $A$ is nonempty; $Q$ may be empty. An exact MDP is a structure with transition function F p α q $=F_{pq}(\alpha)$, stochastic rows, reward $R:Q\to\mathbb R$ and discount $0\le\gamma<1$. $V_{M,\pi}$ is defined as the convergent series $\sum_t\gamma^tP_\pi^tR$, not as a solution of (2). Rewards of the BMDP are tight, as the paper assumes from footnote 3 on; reward intervals are not formalized. Members of $M_\updownarrow$ form a subtype of exact MDPs with the BMDP's $R$ and $\gamma$. The minima and maxima of Definition 3 and of (25)–(26) are the real infimum and supremum over that subtype, which is nonempty and on which values are bounded, so no junk value enters. Intervals are pairs (lower, upper); the orders (17) are spelled out. Policies are deterministic stationary maps $Q\to A$, the paper's $\Pi$. An ordering of $Q$ is a bijection Fin |Q| ≃ Q, and the index $r$ of Definition 1 is the largest index at which expression (9) does not exceed 1 (the page's phrase admits ties; only the largest index keeps the remaining mass inside its interval). Definition 8 is read as $V_\updownarrow{}_{\mathrm{opt}}=V_\updownarrow{}_{\pi_{\mathrm{opt}}}$ for an optimistically optimal $\pi_{\mathrm{opt}}$, so the goal is quantified over all optimistically (pessimistically) optimal policies; Theorem 8 is the milestone that makes this quantification non-vacuous. Lemma 3 is stated as in the Appendix (p. 36): part (d) concerns $\pi_4=\pi_1\oplus_{\mathrm{pes}}\pi_2$, where p. 17 misprints $\pi_3$.

A formalization in which members carry their own discount, rows need not sum to one, the minima range over all functions $Q\to A\to Q\to\mathbb R$, or $V_\updownarrow{}_{\mathrm{opt}}$ is the componentwise supremum over policies would make several statements false or trivial; each is ruled out above.

Contributions welcome beyond proofs of the milestones: a reusable theory of discounted policy evaluation for finite MDPs (the series representation, the fixed-point characterization, monotonicity of $VI_{M,\pi}$), which other missions on finite MDPs can import.

Selected references

• R. Givan, S. Leach, T. Dean, Bounded-parameter Markov decision processes, Artificial Intelligence 122 (2000). https://doi.org/10.1016/S0004-3702(00)00047-3
• A. Nilim, L. El Ghaoui, Robust control of Markov decision processes with uncertain transition matrices, Operations Research 53(5), 2005. https://doi.org/10.1287/opre.1050.0216
• G. Iyengar, Robust dynamic programming, Mathematics of Operations Research 30(2), 2005. https://doi.org/10.1287/moor.1040.0129
• A. Strehl, M. Littman, An analysis of model-based interval estimation for Markov decision processes, Journal of Computer and System Sciences 74(8), 2008. https://doi.org/10.1016/j.jcss.2007.08.009
• M. Puterman, Markov Decision Processes: Discrete Stochastic Dynamic Programming, Wiley, 1994. https://doi.org/10.1002/9780470316887

Why heavy-tailed input matters

Measurements of Ethernet and Internet traffic in the 1990s showed that cumulative traffic is self-similar and long-range dependent: correlations of the input rate decay so slowly that they are not summable, and fluctuations at large time scales do not average out the way Poisson-type models predict (Leland et al. 1994). A widely accepted explanation is that the lengths of individual transmissions (file sizes, session durations) are heavy tailed, with infinite variance. Queueing and capacity-planning calculations then depend on which stochastic process approximates the cumulative input over long horizons.

Two approximations had been proposed: fractional Brownian motion, a Gaussian self-similar process with dependent increments, and α-stable Lévy motion, a heavy-tailed process with independent increments. Mikosch, Resnick, Rootzén and Stegeman (Ann. Appl. Probab. 12 (2002) 23–68) showed that both arise from the same models and that the answer depends on how fast the number of sources grows relative to the time scale. This mission formalizes their Gaussian answer for the superposition of ON/OFF sources: Theorem 4.

• 1995–1997: Willinger, Taqqu, Sherman and Wilson derive fractional Brownian motion from superposed ON/OFF sources with heavy-tailed periods, as an iterated limit: first the number of sources $M\to\infty$, then the time scale $T\to\infty$, and for finite-dimensional distributions only.
• 1998: Heath, Resnick and Samorodnitsky (Math. Oper. Res. 23) give the exact decay of the covariance of a single stationary ON/OFF source.
• 2002: Mikosch, Resnick, Rootzén and Stegeman let $M$ and $T$ grow together and show that fast growth gives fractional Brownian motion as a functional limit in $\mathbb C[0,\infty)$, while slow growth gives stable Lévy motion.

The superposed ON/OFF model

A single ON/OFF source alternates between ON-periods, during which it sends work at rate $1$, and silent OFF-periods. The ON-periods $X_1,X_2,\dots$ are iid with law $F_{\mathrm{on}}$ and the OFF-periods $Y_{\mathrm{off}},Y_1,Y_2,\dots$ iid with law $F_{\mathrm{off}}$, both on $[0,\infty)$, all independent, with means $\mu_{\mathrm{on}}$, $\mu_{\mathrm{off}}$ and $\mu=\mu_{\mathrm{on}}+\mu_{\mathrm{off}}$. The tails are regularly varying: for $x>0$,

$$\bar F_{\mathrm{on}}(x)=x^{-\alpha}L_{\mathrm{on}}(x),\qquad \bar F_{\mathrm{off}}(x)=x^{-\alpha_{\mathrm{off}}}L_{\mathrm{off}}(x),\qquad 1<\alpha<\alpha_{\mathrm{off}}<2,$$

with $L_{\mathrm{on}},L_{\mathrm{off}}$ slowly varying ($L(cx)/L(x)\to1$ for every $c>0$). So both periods have finite mean and infinite variance, and the ON-periods have the heavier tail.

To make the source stationary, a delay $T_0=B(X^{(0)}_{\mathrm{on}}+Y_{\mathrm{off}})+(1-B)Y^{(0)}_{\mathrm{off}}$ is used, where $B$ is Bernoulli with $P(B=1)=\mu_{\mathrm{on}}/\mu$ and $X^{(0)}_{\mathrm{on}},Y^{(0)}_{\mathrm{off}}$ have the integrated-tail laws $F^{(0)}(x)=\mu_F^{-1}\int_0^x\bar F(s)\,ds$. With $T_n=T_0+\sum_{i=1}^n(X_i+Y_i)$ the source's activity is

$$W_t=B\,\mathbf 1_{[0,X^{(0)}_{\mathrm{on}})}(t)+\sum_{n\ge0}\mathbf 1_{[T_n,T_n+X_{n+1})}(t),\qquad t\ge0,$$

a stationary process with $EW_t=\mu_{\mathrm{on}}/\mu$.

The $T$-th model superposes $M=M(T)$ independent copies $W^{(1)},\dots,W^{(M)}$, where $M$ is integer valued, non-decreasing and $M(T)\to\infty$. The cumulative input is $A(t)=\int_0^t\sum_{m=1}^MW^{(m)}_s\,ds$. With the quantile $b(t)=(1/\bar F_{\mathrm{on}})^{\leftarrow}(t)$, the fast growth condition is

$$\text{Condition 2:}\qquad \lim_{T\to\infty}\frac{b(MT)}{T}=\infty,$$

equivalently $MT\bar F_{\mathrm{on}}(T)\to\infty$. The normalisation and limit constants are

$$d_T=[T^{3-\alpha}L_{\mathrm{on}}(T)M]^{1/2},\qquad \sigma_0^2=\frac{2\mu_{\mathrm{off}}^2\Gamma(2-\alpha)/(\alpha-1)}{\mu^3\Gamma(4-\alpha)},\qquad H=\frac{3-\alpha}{2}\in(\tfrac12,1).$$

A standard fractional Brownian motion $B_H$ is a mean-zero Gaussian process on $[0,\infty)$ with continuous paths and $\mathrm{Cov}(B_H(t),B_H(s))=\tfrac12(t^{2H}+s^{2H}-|t-s|^{2H})$.

Formalization targets

Goal: Theorem 4

If Condition 2 holds, then

$$\frac{A(T\cdot)-TM\mu^{-1}\mu_{\mathrm{on}}(\cdot)}{d_T}\ \xrightarrow{d}\ \sigma_0B_H(\cdot)\qquad\text{weakly in }\mathbb C[0,\infty).$$

Intermediate targets

In the order of the paper's argument: the mean $EW_t=\mu_{\mathrm{on}}/\mu$ (p. 27); the covariance decay (2.5) $\gamma_W(h)\sim\frac{\mu_{\mathrm{off}}^2}{(\alpha-1)\mu^3}h^{-(\alpha-1)}L_{\mathrm{on}}(h)$; Condition 2 $\Leftrightarrow T=o(d_T)$ (p. 61); the variance asymptotic (7.1) $\mathrm{Var}(G_T)\sim\sigma_0^2T^{3-\alpha}L_{\mathrm{on}}(T)$ for $G_T=\int_0^T(W_u-EW_u)\,du$; Lemma 13, the one-dimensional limit $d_T^{-1}\sum_mG^{(m)}_{Tt}\to N(0,\sigma_0^2t^{3-\alpha})$; the covariance limit (7.4); and the second-moment bound $E|d_T^{-1}\sum_mG^{(m)}_{Tu}|^2\le c\,u^{1+\varepsilon}$ behind tightness (p. 63). (2.5) and (7.1) are cited by the paper from Heath–Resnick–Samorodnitsky and Willinger–Taqqu–Sherman–Wilson; they are stated here as the paper states them.

Significance

Theorem 4 justifies fractional Brownian motion as the model of aggregate traffic from many heavy-tailed sources in a single, simultaneous limit, replacing the iterated limit of earlier work, and it upgrades finite-dimensional convergence to weak convergence of paths. Functionals of the path, such as the supremum of $A(Tt)-ct$ that governs the buffer content of a fluid queue, therefore inherit Gaussian approximations. Together with Theorem 2 it shows that a single growth condition on $M$ decides between Gaussian and stable approximations.

The result is proved in the paper; it is not formalized anywhere. A formalization adds a machine-checked statement of the theorem and reusable developments of regularly varying functions, the stationary alternating renewal process, fractional Brownian motion, and convergence in $\mathbb C[0,\infty)$.

Difficulty

The obvious route is a central limit theorem for $\sum_{m}G^{(m)}_{Tt}$, a sum of $M$ iid bounded terms; but the number of terms and their distribution change with $T$, so one needs a triangular-array limit theorem whose variance condition rests on the exact asymptotic (7.1). That asymptotic in turn needs the covariance decay (2.5) of a single stationary source, which requires renewal theory for the alternating process with infinite-variance periods. Convergence of finite-dimensional distributions alone does not give the theorem: tightness in $\mathbb C[0,K]$ needs a moment bound on increments that is uniform in $T$ for small time lags, which requires Potter-type bounds for the regularly varying function $x\mapsto EG_x^2$.

Formalization scope

• Time and models. $T\to\infty$ along atTop on $\mathbb R$. Each model $T$ has its own probability space carrying $M(T)$ independent sources; the theorems hold for every such family. Lean's X n, Y n are the paper's $X_{n+1}$, $Y_{n+1}$ and Lean's sources $0,\dots,M-1$ are the paper's $1,\dots,M$. Time in the normalised process runs over $\mathbb R_{\ge0}$.
• Standing hypotheses. $F_{\mathrm{on}},F_{\mathrm{off}}$ are probability measures on $[0,\infty)$ (non-negativity is implicit in "lengths"); (2.1) is stated as "$x\mapsto x^{\alpha}\bar F(x)$ is slowly varying"; $1<\alpha<\alpha_{\mathrm{off}}<2$; $M$ non-decreasing with $M\to\infty$ (§3.2); Condition 2 where §7 assumes it. No other hypothesis is added.
• Junk values ruled out. $b$ is an infimum over $\{x>0:t\bar F_{\mathrm{on}}(x)\le1\}$, never a division by zero; the series defining $W_t$ has non-negative terms and fails to converge only on a null event; equalities of laws come with measurability. The limit is pinned completely: Gaussian, mean zero, the covariance above with $\sigma_H=1$ and $H=(3-\alpha)/2$, continuous paths, multiplied by $\sigma_0$. Fidi convergence in place of the functional limit, a free scale in the limit, or a limit identified only by its marginals would each be a different, weaker theorem.
• Weak convergence is stated in coupling form: for every sequence $T_n\to\infty$ there are one probability space, copies $Y_n$ with continuous paths of the laws of the normalised processes at $T_n$, and a standard fractional Brownian motion $Y'$ such that $Y_n\to\sigma_0Y'$ uniformly on every $[0,K]$ almost surely. On the Polish space $\mathbb C[0,\infty)$ this is equivalent to weak convergence.
• Infrastructure needed: Karamata's theorem and Potter bounds; the stationary alternating renewal process and its covariance asymptotics; a triangular-array central limit theorem; existence of fractional Brownian motion; a moment criterion for tightness in $\mathbb C[0,K]$. The regular-variation and path-space material is reusable well beyond this mission; contributions of any of these pieces, and of any milestone, are welcome.

Selected references

• T. Mikosch, S. Resnick, H. Rootzén and A. Stegeman, Is network traffic approximated by stable Lévy motion or fractional Brownian motion?, Ann. Appl. Probab. 12(1) (2002), 23–68. https://doi.org/10.1214/aoap/1015961155
• D. Heath, S. Resnick and G. Samorodnitsky, Heavy tails and long range dependence in on/off processes and associated fluid models, Math. Oper. Res. 23 (1998), 145–165. https://doi.org/10.1287/moor.23.1.145
• W. Willinger, M. S. Taqqu, R. Sherman and D. V. Wilson, Self-similarity through high-variability: statistical analysis of Ethernet LAN traffic at the source level, IEEE/ACM Trans. Networking 5 (1997), 71–86. https://doi.org/10.1109/90.554723
• W. E. Leland, M. S. Taqqu, W. Willinger and D. V. Wilson, On the self-similar nature of Ethernet traffic (extended version), IEEE/ACM Trans. Networking 2 (1994), 1–15. https://doi.org/10.1109/90.282603
• N. H. Bingham, C. M. Goldie and J. L. Teugels, Regular Variation, Cambridge University Press, 1987. https://doi.org/10.1017/CBO9780511721434
• P. Billingsley, Convergence of Probability Measures, Wiley, 1968 (Theorem 12.3, moment criterion for tightness).

Motivation

Telephone call centers are large service systems: hundreds of agents answer callers who wait in a queue and hang up when their patience runs out. Staffing them is a trade-off between the cost of agents and the quality of service, and in large centers both the efficiency (agents busy almost all the time) and the service level (most callers answered almost immediately) can be high at once. Garnett, Mandelbaum and Reiman (M&SOM 4(3), 2002) study the simplest model that captures customer abandonment, the Erlang-A or $M/M/N+M$ queue, and derive approximations for it in the regime where the number of agents is large.

Their Theorem 2 is the process-level justification of these approximations: the queue-length process, centred at the number of agents and scaled by its square root, converges to a one-dimensional diffusion. The steady-state performance formulas of the paper (delay probability, abandonment probability, waiting times) are computed from this diffusion and its stationary law.

Timeline.

• Halfin and Whitt (1981) proved the corresponding limit for the $M/M/N$ queue without abandonment, in the regime $\sqrt N(1-\rho_N) \to \beta$ with $\beta > 0$ (Oper. Res. 29(3)). The limit is a diffusion with drift $-\mu(\beta + x)$ below $0$ and constant drift $-\mu\beta$ above.
• Fleming, Stolyar and Simon (1994) conjectured the limit with abandonment, with a slightly different centering, and proved the weak limit of the stationary distributions.
• Garnett, Mandelbaum and Reiman (2002) state the process limit for $0 < \theta < \infty$ as Theorem 2, with the extensions $\theta = 0$ and $\theta = \infty$ as Theorem 2* in Appendix C.

Setting

Fix a service rate $\mu > 0$. For each $N \ge 1$ the $N$-th system has $N$ statistically identical agents, Poisson arrivals of rate $\lambda_N > 0$, exponential service times of rate $\mu$, an unlimited waiting room served first come, first served, and exponential patience of rate $\theta_N > 0$ for every caller. A caller whose wait in queue reaches their patience abandons and does not return. The number of callers in the system, $Q_N(t)$, is a birth–death process on $\{0, 1, 2, \dots\}$ with birth rate $\lambda_N$ and death rate $\min(k, N)\mu + (k-N)^+\theta_N$ in state $k$. Its initial value $Q_N(0)$ is an arbitrary random variable.

The traffic intensity is $\rho_N = \lambda_N/(N\mu)$. The centred and scaled process is

$$q_N(t) = \frac{Q_N(t) - N}{\sqrt N}, \qquad t \ge 0 .$$

When $q_N \ge 0$ it counts waiting callers and when $q_N \le 0$ it counts idle agents, both in units of $\sqrt N$.

For $\beta \in \mathbb R$ and $\theta > 0$, the drift is the piecewise-linear function

$$f(x) = \begin{cases} -\mu(\beta + x), & x \le 0,\\ -(\mu\beta + \theta x), & x > 0,\end{cases}$$

and the limit diffusion $q$ solves $dq(t) = f(q(t))\,dt + \sqrt{2\mu}\,db(t)$, where $b$ is a standard Brownian motion independent of $q(0)$. The noise is additive, so this means: $q$ has continuous paths and $q(t) = q(0) + \int_0^t f(q(s))\,ds + \sqrt{2\mu}\,b(t)$ for all $t \ge 0$. Below $0$ it behaves like an Ornstein–Uhlenbeck process with restraining force $\mu$, above $0$ like one with restraining force $\theta$.

Formalization targets

Goal: Theorem 2 (p. 216)

Assume

$$\lim_{N\to\infty}\sqrt N(1-\rho_N) = \beta \in (-\infty,\infty), \qquad \lim_{N\to\infty}\theta_N = \theta \in (0,\infty),$$

and $q_N(0) \Rightarrow \nu$ for a probability law $\nu$ on $\mathbb R$. Then

$$q_N \Rightarrow q \quad\text{in } D[0,\infty),$$

where $q$ is the solution of $dq = f(q)\,dt + \sqrt{2\mu}\,db$ with $q(0) \sim \nu$, and this solution is unique in law.

Milestone: the limit equation is well posed (Theorem 2, p. 216)

On every probability space with a standard Brownian motion $b$ and an independent initial value $X_0$, the equation has a solution adapted to $\sigma(X_0, b(s): s \le t)$, and any two solutions with the same $b$ and initial value are indistinguishable.

Milestone: infinitesimal moments (Appendix C, p. 224)

The infinitesimal expectation $\mu_N(x)$ and variance $\sigma_N^2(x)$ of $q_N$ at the point $x$, as displayed on p. 224, satisfy

$$\lim_{N\to\infty}\mu_N(x) = f(x), \qquad \lim_{N\to\infty}\sigma_N^2(x) = 2\mu \qquad (x \in \mathbb R).$$

Significance

The result. Theorem 2 is the basis for the paper's QED (quality- and efficiency-driven) approximations $Q_N \approx N + q\sqrt N$, which give explicit formulas for the probability of delay, the probability of abandonment and the waiting-time distribution in terms of $\beta$, $\mu$ and $\theta$. Unlike the $M/M/N$ case, $\beta$ may be of either sign: abandonment keeps the system stable even when the offered load exceeds the number of agents. The resulting square-root staffing rule $N = R + \beta\sqrt R$, with $R = \lambda/\mu$, is used in call-center practice.

Formalizing it. The result is proved in the paper by Stone's criteria for birth–death processes, with uniqueness of the limit equation delegated to the literature. No many-server process limit exists on Prove2Me or in Mathlib. A formal proof would give the first machine-checked heavy-traffic process limit for a many-server queue, and the infrastructure (Poisson time-change representation of a birth–death process, path-space weak convergence to a continuous diffusion, well-posedness of a Lipschitz additive-noise SDE) would serve every other diffusion limit in queueing theory.

Difficulty

The infinitesimal-moment computation is elementary; it identifies the limit but does not prove convergence. The hard step is the passage from generator convergence to weak convergence of processes on $D[0,\infty)$: tightness of $\{q_N\}$ in the Skorokhod space, which needs control of the jumps (of size $1/\sqrt N$) and of the unbounded state space, and identification of every limit point as a solution of the equation, which needs uniqueness of the limit. A pointwise limit of the drifts at each fixed $x$ is not enough: the convergence has to be controlled locally uniformly along the paths. A finite-dimensional limit at fixed times is weaker than the theorem.

Formalization scope

• Queue. $Q_N$ is realized by the time-change representation of Mandelbaum, Massey and Reiman (1998), cited by the paper: $Q_N(t) = Q_N(0) + A(\lambda_N t) - S(\mu\int_0^t \min(Q_N(s), N)\,ds) - R(\theta_N \int_0^t (Q_N(s)-N)^+\,ds)$, with $A, S, R$ independent unit-rate Poisson processes independent of $Q_N(0)$. This has the birth–death law above.
• Time. Paths are functions on $\mathbb R$; every condition is imposed for $t \ge 0$ only. Brownian time is $\mathbb R_{\ge 0}$.
• Indexing. $N$ is both the index and the number of agents; every hypothesis on the $N$-th system is imposed for $N \ge 1$. $\lambda_N$, $\theta_N$ are real sequences and $\mu$ is fixed.
• Weak convergence. All $Q_N$ live on one probability space, which is no loss since only their laws matter. $q_N(0) \Rightarrow \nu$ is convergence of $E\,g(q_N(0))$ for every bounded continuous $g$. $q_N \Rightarrow q$ is stated in coupling form (Skorokhod representation with almost-sure uniform convergence on compact time intervals), which is equivalent to $J_1$ weak convergence when the limit is continuous.
• Limit. The SDE is pathwise (Lebesgue integral, no Itô integral). Brownian motion is Mathlib's IsBrownianReal with every path continuous and every coordinate measurable. $q(0)$ is independent of $b$, which is implicit in the paper. The goal asserts existence of a solution with $q_N \Rightarrow q$ and uniqueness of its law; dropping the uniqueness clause, fixing $\nu$ to a point mass, restricting to $\beta > 0$, $\mu = 1$, $\theta = \mu$, or the stationary case is a different theorem.
• Infinitesimal moments. Defined exactly as displayed, with $\lfloor\cdot\rfloor$ the integer part. The display is printed for $\theta = 0$; the milestone is its $0 < \theta < \infty$ instance with the limit $f$ of Theorem 2, pointwise in $x$.
• Not covered. The cases $\theta = 0$ and $\theta = \infty$ of Theorem 2* and the interchange of limits (Part 2).

Contributions welcome: a Poisson time-change existence theorem, tightness criteria in $D[0,\infty)$, well-posedness of Lipschitz SDEs with additive noise, and the martingale-problem approach to identifying limits.

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–588, 1981. https://doi.org/10.1287/opre.29.3.567
• P. J. Fleming, A. Stolyar, B. Simon, Heavy Traffic Limit for a Mobile Phone System Loss Model, Proc. 2nd Int. Conf. on Telecommunication Systems Modeling and Analysis, 1994.
• A. Mandelbaum, W. A. Massey, M. I. Reiman, Strong Approximations for Markovian Service Networks, Queueing Systems 30:149–201, 1998. https://doi.org/10.1023/A:1019112920622
• S. N. Ethier, T. G. Kurtz, Markov Processes: Characterization and Convergence, Wiley, 1986. https://doi.org/10.1002/9780470316658

Motivation

Airlines sell seats on a single flight in several fare classes. Discount classes book first, so the carrier must decide how many seats to protect for later, higher-paying passengers. The classical answer, from Littlewood's rule for two classes to the nested optimality conditions of Brumelle and McGill (1993), computes optimal protection levels from a forecast of the demand distribution of each class. Those forecasts must be built from censored booking data, which is where most of the practical difficulty lies.

van Ryzin and McGill (2000) observed that the optimality conditions themselves can drive an adaptive rule. After each flight one only records whether a fill event occurred (whether demand reached each protection level) and moves the levels by a stochastic-approximation step. No forecast and no optimization is needed. Their Theorem 1 asserts that this rule converges almost surely to the optimal protection levels, with an explicit mean-square rate. It is one of the standard references for model-free revenue management. The convergence proof combines Robbins–Monro theory, the Robbins–Siegmund supermartingale lemma and an induction over fare classes.

Setting

There are $k+1$ fare classes with fares $f_1 > f_2 > \cdots > f_{k+1} > 0$ and discount ratios $r_{i+1} = f_{i+1}/f_1$. On a flight the class demands $X_1, \dots, X_{k+1}$ are independent, nonnegative and continuously distributed. Successive flights $X^1, X^2, \dots$ are independent with the same law.

A protection vector is $\theta = (\theta_1, \dots, \theta_k)$. The fill events are

$$A_i(\theta, X) = \{X_1 > \theta_1,\ X_1 + X_2 > \theta_2,\ \dots,\ X_1 + \cdots + X_i > \theta_i\}, \qquad i = 1, \dots, k.$$

A vector $\theta^*$ is characterised by the optimality condition (3), $P(A_i(\theta^*, X)) = r_{i+1}$ for $i = 1, \dots, k$. The adjustments are $H_i(\theta, X) = r_{i+1} - \mathbf 1(A_i(\theta, X))$, and their means are $h_i(\theta) = r_{i+1} - P(A_i(\theta, X))$. From an arbitrary $\theta^1$ the algorithm updates

$$\theta^{n+1} = \theta^n - \gamma_n H(\theta^n, X^n), \qquad \gamma_n = \frac{A}{n + B},\ A > 0,\ B \ge 0.$$

The interim protection levels $p_i(\theta) = \max\{\theta_j : 1 \le j \le i\}$ are what the airline actually uses to control bookings. The updates themselves monitor $A_i(\theta^n, X^n)$.

Formalization targets

Goal: Theorem 1 (p. 765)

Assume:

• A1: every $X_i$ has bounded support.
• A2 (windowed): for each $R > 0$ there is $\delta_R > 0$ with $(\theta_i - \theta^*_i)\,h_i(\theta_i, \theta^*_{i-1}, \dots, \theta^*_1) \ge \delta_R|\theta_i - \theta^*_i|^2$ whenever $|\theta_i - \theta^*_i| \le R$.
• A3: the distribution functions of the partial sums $X_1 + \cdots + X_i$ are Lipschitz.

Then for $i = 1, \dots, k$,

$$\theta^n_i \to \theta^*_i \ \text{a.s.}, \qquad p_i(\theta^*) = \theta^*_i,$$

and there are $\beta > 0$ and $C$ with

$$E|\theta^n_i - \theta^*_i|^2 \le C\,\gamma_n^{\beta/2^{i-1}}.$$

The exponent is left existential, as in the paper, so any improvement of the rate constant still proves the goal.

Milestones

The milestones follow the paper's own proof, in attack order:

• the base case $i = 1$ (p. 765);
• the Lipschitz bound on $h_{i+1}$ (p. 766);
• Lemma 3, boundedness of the iterates (p. 764);
• the almost-supermartingale inequality (12) (p. 766);
• Lemma 1 (p. 764);
• the a.s. summability of $\gamma_n|\theta^n_i - \theta^*_i|$ (p. 766);
• Lemma 2, Robbins–Siegmund (p. 764);
• the induction step for (9) (p. 766);
• Lemma 4, the rate recursion (p. 764).

Significance

Theorem 1 says that the nested optimality conditions of Brumelle and McGill can be reached by observing binary fill events alone. It also says that the limit is automatically ordered, $\theta^*_1 \le \cdots \le \theta^*_k$, so the interim levels coincide with the limit. The algorithm is distribution-free and needs neither demand forecasts nor uncensoring, which is why it is used in practice and as a baseline for later data-driven and learning-based revenue management.

The theorem is proved in the paper. As far as is known it has no machine-checked proof. Formalizing it requires a vector Robbins–Monro argument with coupled coordinates, the Robbins–Siegmund lemma, and an explicit polynomial rate. Two of the milestones are general probability facts reusable far beyond this mission: Lemma 1, and Robbins–Siegmund, which no platform item states yet. The formalization also corrects the record on two printed slips, described under Formalization scope.

Difficulty

Coordinate $i+1$ is not a Robbins–Monro process by itself. Its fill event $A_{i+1}(\theta, X)$ depends on all of $\theta_1, \dots, \theta_{i+1}$, so the drift of $\theta^n_{i+1}$ is perturbed by the errors of the lower coordinates. Applying a scalar convergence theorem coordinate by coordinate fails at this point. The perturbation has to be bounded via the Lipschitz assumption A3 and shown to be summable along the path, which uses the mean-square rate already established for class $i$. The rate and the almost-sure convergence must therefore be carried through the induction together. The exponent halves at each class. A second difficulty is that A2 is only usable on a bounded set, so boundedness of the iterates (Lemma 3) has to come first and must be uniform.

Formalization scope

Classes are indexed from $1$ in Lean ($\mathbb N$-indexed vectors whose index $0$ and indices above $k+1$ are never read). A flight's demand vector has the product law $\bigotimes_i \nu_i$, and the sample path has the countable product of copies of it. Lean's iterate … n is the paper's $\theta^{n+1}$, so the rate pairs it with $\gamma_{n+1}$. Continuity of a demand law is "every singleton is null". Expectations are Bochner integrals, and integrability is part of each rate conclusion. The fill event reuses the published NestedSeatAlloc.ProbCond.nestEvent.

Pinned and corrected statements:

• A2 as printed is unsatisfiable. It asks for one $\delta$ valid for all $\theta_i$. Since $|h_i| \le 1$, the inequality fails once $|\theta_i - \theta^*_i| > 1/\delta$, so the printed Theorem 1 is vacuously true. The mission uses the windowed A2 above. That is what the proof applies to bounded iterates, and what the paper's sufficient condition (a density bounded below near $\theta^*$) yields. A formalization with the printed A2, or with any other hypothesis no instance satisfies, is ruled out.
• The base case. The page claims the rate "for all $0 < \beta < 1$". That is false in general. The milestone states "for some $\beta > 0$", which is all Theorem 1 uses.
• Integrability. Lemma 2 assumes each $Z_n$ integrable, so that the conditional expectation is genuine. Lemma 3's "bounded (a.s.)" is read as a bound uniform in $n$.
• Standing assumptions. Nonnegative demands and positive fares are added where the argument uses them.

Not posed: §4's modifications (projection onto the capacity, randomized integer levels, booking lead times) and the simulations of §5.

Contributions welcome: proofs of the general lemmas (Lemma 1, Robbins–Siegmund) stated for arbitrary probability spaces, measurability facts for the recursion, and the elementary Lemma 4.

Selected references

• G. van Ryzin and J. McGill, Revenue Management Without Forecasting or Optimization: An Adaptive Algorithm for Determining Airline Seat Protection Levels, Management Science 46(6), 760–775, 2000. https://doi.org/10.1287/mnsc.46.6.760.11936
• S. L. Brumelle and J. I. McGill, Airline Seat Allocation with Multiple Nested Fare Classes, Operations Research 41(1), 127–137, 1993. https://doi.org/10.1287/opre.41.1.127
• H. Robbins and D. Siegmund, A convergence theorem for non negative almost supermartingales and some applications, in Optimizing Methods in Statistics, Academic Press, 233–257, 1971. https://doi.org/10.1016/B978-0-12-604550-5.50015-8
• A. Benveniste, M. Métivier and P. Priouret, Adaptive Algorithms and Stochastic Approximations, Springer, 1990. https://doi.org/10.1007/978-3-642-75894-2
• E. Lukacs, Stochastic Convergence, 2nd ed., Academic Press, 1975.