Probability

Motivation

Measurements of Ethernet and Internet traffic in the 1990s showed that the amount of data offered to a link is bursty on every time scale and that the lengths of transmissions (file sizes, connection durations) have heavy, regularly varying tails. Two families of approximations for the cumulative input were proposed: fractional Brownian motion, which has dependent Gaussian increments, and α-stable Lévy motion, which has independent, heavy-tailed increments. The two lead to very different predictions for buffer overflow and link dimensioning.

T. Mikosch, S. Resnick, H. Rootzén and A. Stegeman (Ann. Appl. Probab. 12 (2002) 23–68) showed that both answers are correct in different regimes, and that the regime is decided by how fast the connection rate grows relative to the time scale. This mission formalizes their result for the infinite source Poisson model under slow growth, where the limit is a totally skewed α-stable Lévy motion. Companion missions treat the ON/OFF model and the fast-growth regime with its fractional Brownian limit.

Setting

A transmission length has law $F_{\mathrm{on}}$ on $[0,\infty)$ with tail $\bar F_{\mathrm{on}}(x)=F_{\mathrm{on}}((x,\infty))$. Condition (2.8) asks

$$\bar F_{\mathrm{on}}(x)=x^{-\alpha}L(x),\qquad x>0,\quad 1<\alpha<2,$$

with $L$ slowly varying ($L(cx)/L(x)\to1$ for every $c>0$). The mean $\mu_{\mathrm{on}}$ is finite and the variance infinite. The quantile function (2.9) is $b(t)=(1/\bar F_{\mathrm{on}})^{\leftarrow}(t)=\inf\{x:1/\bar F_{\mathrm{on}}(x)\ge t\}$.

In the $T$-th model, connections start at the points $(\Gamma_k)_{k\in\mathbb Z}$ of a homogeneous Poisson process on $\mathbb R$ with rate $\lambda=\lambda(T)$, labelled so that $\Gamma_0<0<\Gamma_1$. Connection $k$ transmits at unit rate for a time $X_k$; the $X_k$ are iid with law $F_{\mathrm{on}}$ and independent of the $\Gamma_k$. The number of active connections and the cumulative input are

$$N(t)=\sum_{k\in\mathbb Z}\mathbf 1[\Gamma_k\le t<\Gamma_k+X_k],\qquad A(t)=\int_0^tN(s)\,ds .$$

The rate $\lambda(T)$ is a positive non-decreasing function of $T$. Slow Growth Condition 1 is

$$\lim_{T\to\infty}\frac{b(\lambda T)}{T}=0,$$

where $b(\lambda T)$ is $b$ evaluated at $\lambda(T)\,T$.

The stable law $S_\alpha(\sigma,\beta,\mu)$ is the law with characteristic function $\exp\{-\sigma^\alpha|\theta|^\alpha(1-i\beta\,\mathrm{sign}(\theta)\tan(\pi\alpha/2))+i\mu\theta\}$ for $\alpha\ne1$. α-stable Lévy motion $X_{\alpha,\sigma,\beta}$ has independent stationary increments with $X(t)-X(s)\sim S_\alpha(\sigma(t-s)^{1/\alpha},\beta,0)$. Write

$$C_\alpha=\frac{1-\alpha}{\Gamma(2-\alpha)\cos(\pi\alpha/2)}\quad(5.1).$$

Formalization targets

Goal: Theorem 1 (corrected scale)

Under (2.8) and Condition 1,

$$\frac{A(T\cdot)-T\lambda\mu_{\mathrm{on}}(\cdot)}{b(\lambda T)}\xrightarrow{\ fidi\ }X_{\alpha,C_\alpha^{-1/\alpha},1}(\cdot),$$

convergence of the finite-dimensional distributions to a totally skewed ($\beta=1$) α-stable Lévy motion.

On the scale. The paper prints the limit $X_{\alpha,1,1}$. Its own proof establishes $\lambda T\,P(j_1>b(\lambda T)x)\to x^{-\alpha}$ (pp. 37–38), which makes $\alpha x^{-\alpha-1}dx$ the Lévy measure of the limit. The totally skewed stable law with that Lévy measure has $\sigma^\alpha=\Gamma(1-\alpha)\cos(\pi\alpha/2)=C_\alpha^{-1}$ (Samorodnitsky–Taqqu, Property 1.2.15; also the paper's own criterion on p. 47 with $c=1$, and the $\sigma$ the paper prints in Theorem 2). Since $C_\alpha$ decreases from $2/\pi$ to $0$ on $(1,2)$ ($C_{1.5}\approx0.399$), the printed limit has the wrong scale and the theorem as printed is false. The goal states the corrected scale $C_\alpha^{-1/\alpha}$; milestones quote the page as printed.

Milestones, in attack order

(2.14), the covariance of $N$; Lemma 1 part 1 (Condition 1 $\iff\lambda T\bar F_{\mathrm{on}}(T)\to0\iff\mathrm{Cov}(N_T(0),N_T(T))\to0$); Lemma 2, slow part ($\lambda T^2\bar F_{\mathrm{on}}(T)/b(\lambda T)\to0$); (4.15), negligibility of the boundary pieces $A_2,A_3,A_4$; the tail limit of $j_1$; (4.17), $A_{12}=O_P([\lambda T]^{1/2})=o_P(b(\lambda T))$; (4.20)–(4.21), $A_{13}=o(b(\lambda T))$; (4.19), corrected, the stable limit of $A_{11}$; the corrected one-dimensional limit of $A(T)$; and $EA_{22}=o(b(\lambda T))$ of §4.5, which reduces the fidi convergence to the one-dimensional one.

Significance

The theorem shows that under slow connection growth the cumulative input is asymptotically a process with independent, infinite-variance increments. Long-range dependence present in each model (the covariance (2.14) decays like $h^{-(\alpha-1)}L(h)$) disappears on the time scale $T$, and the heavy tail of the transmission lengths dominates. Together with the fast-growth result (fractional Brownian motion) it explains why both approximations are found in measured traffic, and it identifies the critical quantity $b(\lambda T)/T$.

The result is proved in the paper, modulo the scale. No part of it is machine-checked. The mission produces a formal statement of the infinite source Poisson model with its marked point process structure, a formal statement of stable laws and stable Lévy motion through their characteristic functions, and machine-checked versions of the paper's estimates. Correcting the printed scale is part of the output.

Difficulty

The pieces $A_2,A_3,A_4,A_{12},A_{13},A_{22}$ are controlled by first moments, regular variation and Lemma 2. The central step is the stable limit of $A_{11}$: a sum of a Poisson number of iid, heavy-tailed summands whose law changes with $T$. A classical central limit theorem does not apply because the variance is infinite, and the summands form a triangular array. The paper cites a point-process limit (Resnick's Exercise 4.4.2.8) for this step; a formal proof needs a convergence theorem for row-wise iid triangular arrays to stable laws, with the constant $C_\alpha$ computed exactly.

Formalization scope

Each $T$ has its own probability space carrying the $T$-th model, and every statement quantifies over all such families. The Poisson points are given by their iid exponential spacings $-\Gamma_0,\Gamma_1,\Gamma_{k+1}-\Gamma_k$ $(k\ne0)$, indexed by $\mathbb Z$. (2.8) is stated as "$x\mapsto x^\alpha\bar F_{\mathrm{on}}(x)$ is slowly varying", with lengths non-negative. $b(t)$ is $\inf\{x>0:t\bar F_{\mathrm{on}}(x)\le1\}$, which equals the page's generalized inverse for $t>0$. $N(t)$ is a cardinality and the region sums are sums of non-negative terms; their junk values (an infinite set, a non-summable family) occur only on null events. $A$ is an interval integral. The growth hypothesis is: $\lambda>0$, non-decreasing, Condition 1. No other hypothesis is added.

The fidi limit is pinned by its characteristic function on $\mathbb R^k$ for sorted times $0\le t_1\le\dots\le t_k$: the product of the characteristic functions of $S_\alpha(\sigma(t_j-t_{j-1})^{1/\alpha},1,0)$ evaluated at $\theta_j+\dots+\theta_k$. The existence of the limit law is part of the claim. One-dimensional limits are stated the same way on $\mathbb R$. Convergence in probability to $0$ is $P_T(|Y_T|>\eta)\to0$ for every $\eta>0$.

Trivializing formalizations are excluded:

• $b$ is never the junk value $0$, since its defining set is nonempty and bounded below.
• The model hypothesis is satisfiable for every rate and length law; a sanity file builds it from product measures, together with a Pareto law satisfying (2.8) and a rate satisfying Condition 1.
• The limit object is fixed by its characteristic function, not merely asserted to exist.
• Only fidi convergence is claimed, as on the page. Convergence in the Skorokhod $J_1$ topology fails for this model (Remark after Theorem 1).

A complete development needs regular variation (Karamata's theorem, Potter bounds, generalized inverses), the Poisson random measure of the marked points, and stable limits for triangular arrays. All three are reusable beyond this mission. Proofs of the milestones are welcome separately, as are general Mathlib-style lemmas on regular variation.

Selected references

• T. Mikosch, S. Resnick, H. Rootzén, 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
• G. Samorodnitsky, M. S. Taqqu, Stable Non-Gaussian Random Processes, Chapman & Hall, 1994. https://doi.org/10.1201/9780203738818
• N. H. Bingham, C. M. Goldie, J. L. Teugels, Regular Variation, Cambridge University Press, 1987. https://doi.org/10.1017/CBO9780511721434
• S. I. Resnick, Extreme Values, Regular Variation, and Point Processes, Springer, 1987. https://doi.org/10.1007/978-0-387-75953-1

Motivation

The ferromagnetic Ising model is the basic model of statistical mechanics for interacting binary variables. On a finite graph it is a probability distribution over spin configurations, and its normalizing constant, the partition function, encodes the thermodynamics: the free entropy per vertex determines the magnetization, the correlations, and the location of phase transitions. On lattices this quantity is understood only in special cases. On sparse random graphs, such as random regular graphs and Erdős–Rényi graphs with bounded average degree, physicists predicted it with the cavity method (Mézard and Parisi, The Bethe lattice spin glass revisited, 2001): locally these graphs look like trees, and on trees the model can be solved by a recursion. The same recursion underlies belief propagation, the message-passing algorithm used in coding theory and machine learning.

Timeline:

• Dembo and Montanari (2008/2010), Ising models on locally tree-like graphs, Ann. Appl. Probab. 20 (2010) 565–592: for every uniformly sparse graph sequence converging locally to a Galton–Watson tree whose offspring law has finite mean, the free entropy density converges at every temperature $\beta \ge 0$ and every field $B$, and for $B > 0$ the limit is the Bethe prediction. Earlier rigorous results covered special cases, among them Erdős–Rényi graphs in a restricted range of parameters (reference [7] of the paper) and random regular graphs (reference [11]).
• Dembo, Montanari and Sun (2013), Factor models on locally tree-like graphs, extended the approach to a class of factor models on locally tree-like graphs.

Setting

Let $G_n = (V_n = [n], E_n)$ be a sequence of finite graphs. The Ising measure at inverse temperature $\beta \ge 0$ and magnetic field $B \in \mathbb R$ is

$$\mu(\underline x) = \frac{1}{Z_n(\beta,B)} \exp\Big\{\beta \sum_{(i,j)\in E_n} x_i x_j + B \sum_{i \in V_n} x_i\Big\}, \qquad x_i \in \{+1,-1\},$$

each edge counted once, and the free entropy density is $\phi_n(\beta,B) = \frac1n \log Z_n(\beta,B)$.

A degree distribution is a law $P = \{P_k\}$ on $\{0,1,2,\dots\}$ with $0 < \bar P = \sum_k k P_k < \infty$; its size-biased version is $\rho_k = k P_k/\bar P$, and $\bar\rho = \sum_k (k-1)\rho_k$. The Galton–Watson tree $T(P,\rho,\infty)$ has a root with $k$ children with probability $P_k$, and every other vertex independently has $k-1$ children with probability $\rho_k$; $T(P,\rho,t)$ is its first $t$ generations.

The ball $B_i(t)$ is the subgraph induced by the vertices within distance $t$ of $i$, rooted at $i$. The sequence converges locally to $T(P,\rho,\infty)$ if, for every $t$ and every rooted tree $T$, the fraction of vertices $i$ with $B_i(t) \simeq T$ (root-preserving isomorphism) tends to $\mathbb P\{T(P,\rho,t) \simeq T\}$. It is uniformly sparse if $\lim_{l\to\infty}\limsup_n \frac1n\sum_i |\partial i|\,\mathbb I(|\partial i| \ge l) = 0$.

The cavity recursion acts on laws of a real random variable $h$:

$$h \overset{d}{=} B + \sum_{i=1}^{K-1} \xi(\beta, h_i), \qquad \xi(\beta,h) = \operatorname{atanh}[\tanh\beta \tanh h],$$

with $K \sim \rho$ and $h_i$ i.i.d. copies of $h$ independent of $K$. For $B > 0$ it has a unique fixed point $h^*$ supported on $[0,\infty)$. The Bethe functional of a law of $h$ is

$$\varphi_h(\beta,B) = \frac{\bar P}{2}\log\cosh\beta - \frac{\bar P}{2}\mathbb E\log[1+\tanh\beta\tanh h_1\tanh h_2] + \mathbb E\log\Big\{e^B\prod_{i=1}^L[1+\tanh\beta\tanh h_i] + e^{-B}\prod_{i=1}^L[1-\tanh\beta\tanh h_i]\Big\},$$

with $L \sim P$ independent of the $h_i$.

Formalization targets

Goal: Theorem 2.4

Assume $\{G_n\}$ is uniformly sparse, converges locally to $T(P,\rho,\infty)$, and $\sum_k k\rho_k < \infty$. Then for every $\beta \ge 0$:

$$\lim_{n\to\infty}\phi_n(\beta,B) = \varphi_{h^*}(\beta,|B|) \quad (B \ne 0),$$

and $\varphi_{h^*}(\beta,B)$ has a limit $L_0$ as $B \downarrow 0$, with $\phi_n(\beta,0) \to L_0$. This is the paper's statement that the limit exists for every $B$, equals (2.9) for $B > 0$, is even in $B$, and is continuous at $B = 0$.

Milestones

1. Tools on Ising models: Griffiths' inequality (Theorem 3.1), the GHS inequality (Theorem 3.2), the subtree-marginal identity (Lemma 4.1).
2. Trees: the boundary effect on the root magnetization, $\mathbb E\{m^{\ell,+}-m^{\ell,0}\} \le M/\ell$ (Lemma 4.3); the correlation factorization (Lemma 4.4); exponential decay of root-to-generation correlations (Corollary 4.5, (4.13)).
3. The recursion: monotone convergence to the unique nonnegative fixed point (Lemma 2.3) and its Lipschitz dependence on $\beta$ in Monge–Kantorovich distance (Lemma 4.6).
4. The Bethe functional: the second-order perturbation bound (Lemma 6.1, Remark 6.2) and its consequence, stationarity of $\varphi_h$ at fixed points (Corollary 6.3, with corrected constant).
5. From trees to graphs: edge averages converge to averages over the edge-rooted tree $\bar T(\rho,t)$ (Lemma 6.4), and $\partial_\beta\phi_n = \frac1n\sum_{(i,j)\in E_n}\langle x_ix_j\rangle_n$ (6.10).

Significance

The result. Theorem 2.4 establishes the replica-symmetric (Bethe) formula for the ferromagnetic Ising model on the standard sparse random graph ensembles at all temperatures, including below the critical temperature $\tanh\beta_c = 1/\bar\rho$ where the Gibbs measure on the limiting tree is not unique. It identifies the free energy of random regular graphs, Erdős–Rényi graphs and configuration models with arbitrary degree distributions of finite variance, and the derivative identities that follow give the limiting magnetization and energy per edge. The tree estimates of Section 4 also yield the convergence of belief propagation (the companion mission of this series).

Formalizing it. The theorem is proved; nothing in this mission has been machine-checked before. A formal proof requires infrastructure that Lean's Mathlib does not have: Ising measures with boundary conditions, correlation inequalities, Galton–Watson trees as laws of offspring functions, local weak convergence of graph sequences, distributional fixed-point equations, and the Monge–Kantorovich distance on the real line. While stating the milestones, Corollary 6.3 turned out to be false as printed (the constant $\bar P\bar\rho$ should be $\mathbb E[L^2]$; a counterexample is $P = \delta_1$); the mission states the corrected version, which is what the paper's proof gives.

Difficulty

The obvious route through uniqueness of the Gibbs measure on the tree fails below the critical temperature: different boundary conditions give different limits, so local quantities are not determined by the local graph structure alone. The paper replaces uniqueness by insensitivity to positive boundary conditions when $B > 0$ (Lemma 4.3), which rests on Griffiths and GHS monotonicity and on averaging over the random tree; the deterministic statement is false. Passing from derivatives to the free entropy itself requires that the Bethe functional be stationary in the cavity law (Corollary 6.3), so that the $\beta$-dependence of $h^*$ can be ignored. Finally, the case $B = 0$ at low temperature is reached only as a limit: there is no fixed-point characterization at $B = 0$.

Formalization scope

• Spins are Booleans ($+1$ = true); edges are unordered and counted once. Graphs $G_n$ live on Fin n.
• A field $B_i = +\infty$ (plus boundary) is encoded by pinning vertex $i$ to $+1$, never by a large real field.
• Rooted trees are offspring functions List ℕ → ℕ in Ulam–Harris form; random trees are laws of offspring functions under the product σ-algebra. Galton–Watson trees are infinite product measures. Fields on trees are nonrandom functions of the vertex.
• Laws of cavity fields are probability measures on $\mathbb R$; the recursion is a map on ProbabilityMeasure ℝ, and convergence in distribution is weak convergence. Expectations over $L \sim P$ are series $\sum_l P_l\,\mathbb E[\cdot]$ with product laws.
• The Monge–Kantorovich distance takes values in $[0,\infty]$. Uniform sparsity is computed in $[0,\infty]$ so that an unbounded sequence does not get a junk $\limsup$ of $0$.
• The fixed point $h^*$ in the goal is the one supported on $[0,\infty)$, selected by choice; Lemma 2.3 shows it exists and is unique, so an arbitrary fixed point cannot be substituted. Theorem 2.4 is stated with $\beta \ge 0$ and asserts convergence for $B > 0$, $B < 0$ and $B = 0$ separately; the $B = 0$ clause asserts the existence of the one-sided limit rather than assuming it.

A complete development needs: Griffiths and GHS inequalities for finite Ising models (reusable for any ferromagnetic model); Ulam–Harris random trees with conditional independence; local weak convergence of graph sequences (reusable for random graph theory); distributional recursions with stochastic ordering; and the coupling characterization of the Wasserstein-1 distance on $\mathbb R$. Contributions of any milestone are welcome, as are proofs of the classical inequalities of Section 3 independent of this paper.

Selected references

• A. Dembo, A. Montanari, Ising models on locally tree-like graphs, Ann. Appl. Probab. 20(2) (2010) 565–592. https://arxiv.org/abs/0804.4726v3 , https://doi.org/10.1214/09-AAP627
• M. Mézard, G. Parisi, The Bethe lattice spin glass revisited, Eur. Phys. J. B 20 (2001) 217–233. https://doi.org/10.1007/PL00011099
• R. B. Griffiths, C. A. Hurst, S. Sherman, Concavity of magnetization of an Ising ferromagnet in a positive external field, J. Math. Phys. 11 (1970) 790–795. https://doi.org/10.1063/1.1665211
• T. M. Liggett, Interacting Particle Systems, Springer, 1985 (Theorem IV.1.21, Griffiths' inequality). https://doi.org/10.1007/978-1-4613-8542-4
• A. Dembo, A. Montanari, N. Sun, Factor models on locally tree-like graphs, Ann. Probab. 41(6) (2013) 4162–4213. https://doi.org/10.1214/12-AOP828

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

Motivation

Facility location models decide where to open depots, warehouses or service centres and which customers each one serves. Classical models assume that an open facility always works. In practice facilities are disrupted by weather, strikes, power loss or supplier failure, and a network designed for normal operation can become very expensive when its nearest facilities go down. Reliable facility location models place facilities so that the sum of the fixed cost and the expected transportation cost, taken over random facility failures, is minimal.

Cui, Ouyang and Shen (UCTC-FR-2010-02, 2010; published as Oper. Res. 58(4):998–1011) study the case in which each site fails independently with its own probability. The most direct model of this situation is a scenario-based stochastic program: list every pattern of failures, give it its probability, and let each customer be served optimally in each pattern. That program has $2^J$ scenarios for $J$ candidate sites and cannot be written down for realistic $J$. The paper proposes instead a compact formulation in which each customer receives an ordered list of backup facilities ("levels") and the probability that each level is the one that serves her is carried by a small set of recursive variables. This mission is about the claim that the compact model loses nothing.

The level-assignment technique comes from Snyder and Daskin (Transp. Sci. 39(3):400–416, 2005), who treated the case of one common failure probability for all sites. The uniform-probability model of Snyder and Shen's textbook is on the platform as SupplyChainTheory_disruptions; it is a different model (one $q$, no cap on the number of levels) and is not reused here.

Setting

There are $I$ customers $i = 0,\dots,I-1$ with demand rates $\lambda_i \ge 0$ and $J$ candidate sites $j = 0,\dots,J-1$ with fixed costs $f_j$ and failure probabilities $0 \le q_j < 1$. Failures are independent. Shipping one unit from site $j$ to customer $i$ costs $d_{ij}$, and each unit of missed demand of customer $i$ costs a penalty $\phi_i$. An emergency facility with index $J$ never fails ($q_J = 0$), costs nothing to open ($f_J = 0$), and "serves" at the penalty cost $d_{iJ} = \phi_i$.

The compact model (RUFL). Each customer is assigned to facilities at levels $r = 0,\dots,R$, with $R \ge 1$; a level-$r$ facility serves her exactly when all her facilities at levels $0,\dots,r-1$ have failed. The variables are $X_j\in\{0,1\}$ (site $j$ open), $Y_{ijr}\in\{0,1\}$ (facility $j$ is customer $i$'s level-$r$ facility) and $P_{ijr}$, the probability that $j$ serves $i$ at level $r$. Constraints (1b)–(1d) say that each level holds one facility until the emergency facility appears, which happens exactly once, and that only open sites are used, each at most once. Constraints (1e)–(1f) fix $P$: $P_{ij0} = 1-q_j$ and $P_{ijr} = (1-q_j)\sum_{k<J}\frac{q_k}{1-q_k}P_{i,k,r-1}Y_{i,k,r-1}$. The objective is

$$\Phi(X,Y,P) = \sum_{j=0}^{J-1} f_jX_j + \sum_{i=0}^{I-1}\sum_{j=0}^{J}\sum_{r=0}^{R}\lambda_i d_{ij}P_{ijr}Y_{ijr}.$$

The scenario model (SSP). A scenario is $\omega\in\Omega=\{0,1\}^J$, with $\delta_{j\omega}=1$ if site $j$ operates in $\omega$ (and $\delta_{J\omega}=1$ always). Its probability is $p_\omega = \prod_{j<J}(1-q_j)^{\delta_{j\omega}}q_j^{1-\delta_{j\omega}}$. With $Y_{ij\omega}\in\{0,1\}$ (customer $i$ served by $j$ in $\omega$), (SSP) minimizes

$$\Psi(X,Y) = \sum_{j=0}^{J-1} f_jX_j + \sum_{i=0}^{I-1}\sum_{j=0}^{J}\sum_{\omega\in\Omega}\lambda_i d_{ij}p_\omega Y_{ij\omega}$$

subject to $\sum_{j=0}^{J}Y_{ij\omega}=1$ and $Y_{ij\omega}\le\delta_{j\omega}X_j$ for regular $j$.

Formalization targets

Goal: Proposition 1 (p. 10)

If $R = J$, both programs have optimal solutions, and for every optimal $(X,Y,P)$ of (RUFL) and every optimal $(X',Y')$ of (SSP),

$$\Phi(X,Y,P) = \Psi(X',Y').$$

The paper states this as "formulation (1a)–(1g) is equivalent to the stochastic programming formulation that covers all failure scenarios"; its proof shows equality of optimal values, and that is the reading formalized.

Milestones (Appendix A.1, pp. 32–34)

1. Scenario-probability identity. For distinct regular sites $j(0),\dots,j(r-1)$ and a facility $j(r)$ not among them, the scenarios in which $j(0),\dots,j(r-1)$ fail and $j(r)$ operates have total probability

$$\sum_{\omega\in\Omega(i,r)}p_\omega = (1-q_{j(r)})\prod_{\ell=0}^{r-1}q_{j(\ell)}.$$

2. (RUFL) → (SSP). Every (RUFL)-feasible $(X,Y,P)$ has an (SSP)-feasible $(X,Y')$ with $\Psi(X,Y') = \Phi(X,Y,P)$, for any $R$.
3. Normalization. Some optimal solution of (SSP), with the same $X$, serves every customer in every scenario by her closest operating open facility, ties to the lowest index.
4. (SSP) → (RUFL). If $R = J$, every normalized (SSP)-feasible $(X,Y)$ has a (RUFL)-feasible $(X,Y',P')$ with $\Phi(X,Y',P') = \Psi(X,Y)$.

Significance

The proposition says that, with enough levels, the compact model is exact: its optimal value is the true minimum expected cost over all failure patterns. The compact model has $O(IJR)$ variables and constraints against $O(IJ2^J)$ for (SSP), and after the standard linearization of the products $P_{ijr}Y_{ijr}$ it is an ordinary mixed-integer program to which the paper's Lagrangian relaxation applies. For $R<J$ the paper notes the compact model is in general not equivalent; Proposition 1 is the anchor that says what the parameter $R$ trades away.

The result is proved in the paper; as far as we know it has not been formalized. The work here is to formalize that proof: the product-measure computation behind milestone 1, the two explicit solution maps, and the exchange argument behind the normalization. The probability identity for "first success in an ordered list of independent trials" and the scenario-to-level bookkeeping are reusable for other reliability models built on Snyder–Daskin levels.

Difficulty

The cost identities are routine once the right bookkeeping is in place; the obstacle is that bookkeeping. In direction (RUFL) → (SSP) one must show that the recursion (1e)–(1f), which multiplies by $q_k/(1-q_k)$ and sums over all regular $k$, collapses to the closed product $(1-q_{j(r)})\prod_{\ell<r}q_{j(\ell)}$ along a customer's assigned list, using that (1b)–(1d) put exactly one facility on each level until the emergency facility and that (1c) forbids repeats. Then the sum over the $2^J$ scenarios must be regrouped by the first operating facility on that list, which needs milestone 1 on a product space.

The converse direction is where $R=J$ enters, and the naive attempt fails: an arbitrary optimal (SSP) solution may serve a customer by different, equally good facilities in different scenarios in ways that no single ordered list reproduces. The normalization step removes this freedom; only after it can one read off a single list per customer, and $R=J$ guarantees that this list (all open sites no farther than the penalty, then $J$) fits on the available levels.

Formalization scope

Indices are 0-based; facilities are Fin (J+1) with Fin.last J the emergency facility, and levels are Fin (R+1). The extended data $d_{iJ}=\phi_i$, $q_J=0$, and the conventions $\delta_{J\omega}=1$, $X_J=1$ are definitions, not hypotheses. Scenarios are Fin J → Bool and $p_\omega$ is the explicit finite product; no measure theory is involved, and (17a) is a finite sum. Binary variables are real numbers constrained to $\{0,1\}$. $R = J$ is expressed by taking the instance type Instance I J J; with the standing $R\ge 1$ this means $J\ge1$. Demand rates are nonnegative; no sign is assumed on $d$, $\phi$ or $f$.

Two printed constraints are corrected, and the formalization commits to the corrections. The first sum in (1b) is printed over the regular sites $j\le J-1$; it runs here over all $J+1$ facilities, as the paper's own reading of (1b) requires. Constraint (17c) is printed as $\sum_i Y_{ij\omega}\le\delta_{j\omega}X_j$, which gives each site a capacity of one customer per scenario and makes Proposition 1 false; it is stated here per customer, $Y_{ij\omega}\le\delta_{j\omega}X_j$.

The goal asserts existence of optimal solutions of both programs before comparing their values, so it cannot hold vacuously. Contributions welcome: proofs of the milestones, a general lemma on the probability that the first success among independent Bernoulli trials occurs at a given position, and an existence-of-optimum lemma for finite binary programs.

Selected references

• X. Cui, Y. Ouyang, Z.-J. M. Shen, Reliable Facility Location Design under the Risk of Disruptions, UCTC-FR-2010-02, University of California Transportation Center, 2010; Operations Research 58(4):998–1011, 2010. https://doi.org/10.1287/opre.1090.0801
• L. V. Snyder, M. S. Daskin, Reliability models for facility location: the expected failure cost case, Transportation Science 39(3):400–416, 2005. https://doi.org/10.1287/trsc.1040.0107
• H. D. Sherali, A. Alameddine, A new reformulation-linearization technique for bilinear programming problems, Journal of Global Optimization 2(4):379–410, 1992. https://doi.org/10.1007/BF00122429

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]).

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).