Probability

Why lower bounds for kernel regression matter

Regularized least squares (RLS), also called kernel ridge regression, is one of the standard estimators of statistical learning: given a sample of input–output pairs, it fits a function from a reproducing kernel Hilbert space by minimizing the empirical squared error plus a multiple of the squared norm. Caponnetto and De Vito (FoCM 2007) proved that, over a class of distributions described by two parameters — the decay of the eigenvalues of the kernel's covariance operator and the smoothness of the regression function relative to it — RLS with a well-chosen regularization parameter converges at the rate $\ell^{-bc/(bc+1)}$ in the sample size $\ell$ (their Theorem 1). An upper rate on its own does not say whether another method could do better. This mission formalizes the matching minimax lower rate (their Theorem 2): when the output space is finite dimensional, no learning algorithm whatsoever converges faster than $\ell^{-bc/(bc+1)}$ uniformly over the class. Together the two theorems say that RLS is rate-optimal, and they are the reference point for later work on spectral regularization, early stopping and distributed kernel methods.

The paper's lower bound adapts the minimax analysis of DeVore, Kerkyacharian, Picard and Temlyakov (FoCM 2006) to the vector-valued RKHS setting.

Setting

Inputs $x$ lie in a Polish space $X$ and outputs $y$ in a real Hilbert space $Y$ of finite dimension $d$. The hypothesis space $\mathcal H$ is a separable Hilbert space of functions $f : X \to Y$ in which evaluation is continuous; $K_x : Y \to \mathcal H$ is the adjoint of evaluation at $x$, so $f(x) = K_x^* f$. Hypothesis 1 asks that $(x,t) \mapsto \langle K_t v, K_x w\rangle_{\mathcal H}$ be measurable and that $\operatorname{Tr}(K_x^*K_x) \le \kappa$ for every $x$.

A distribution $\rho$ on $Z = X \times Y$ has risk $\mathcal E[f] = \int \|f(x) - y\|^2\, d\rho$. Hypothesis 2 asks that $\int\|y\|^2 d\rho < \infty$, that the risk has a minimizer $f_{\mathcal H} \in \mathcal H$ (taken of minimal norm), and that the noise $y - f_{\mathcal H}(x)$ satisfies a Bernstein moment condition with constants $M, \Sigma$. With $\rho_X$ the marginal of $\rho$, the operator $T = \int K_x K_x^*\, d\rho_X(x)$ on $\mathcal H$ has the quadratic form $\langle Tf, g\rangle = \int \langle f(x), g(x)\rangle\, d\rho_X$ and a spectral decomposition $T = \sum_n t_n \langle\cdot, e_n\rangle e_n$.

The prior $\mathcal P(b,c)$ (Definition 1, with fixed positive $M, \Sigma, R, \alpha, \beta$, $1 < b < \infty$ and $1 \le c \le 2$) is the set of probability measures $\rho$ satisfying Hypothesis 2 with $M, \Sigma$ such that

• $T$ has infinitely many positive eigenvalues $t_1 \ge t_2 \ge \dots$ with $\alpha \le n^b t_n \le \beta$ (a capacity condition);
• $f_{\mathcal H} = T^{(c-1)/2} g$ for some $g$ with $\|g\|^2 \le R$ (a source condition).

Formalization targets

Goal: Theorem 2 (p. 11)

$$\lim_{\tau \to 0}\ \liminf_{\ell \to \infty}\ \inf_{f_\ell}\ \sup_{\rho \in \mathcal P(b,c)}\ \mathbb P_{\mathbf z \sim \rho^\ell}\Big[\mathcal E[f^\ell_{\mathbf z}] - \mathcal E[f_{\mathcal H}] > \tau\, \ell^{-\frac{bc}{bc+1}}\Big] = 1,$$

the infimum over all learning algorithms $f_\ell : Z^\ell \to \mathcal H$. The constant in front of the rate is left free, so the goal asserts only the exponent.

Milestones

1. Proposition 4 (p. 21): for $f = T^{(c-1)/2}g$, $\|g\|^2 \le R$, the explicit distribution $\rho_f$ with marginal $\nu$ (the marginal of some $\rho_0 \in \mathcal P(b,c)$) and $2d$-point conditional law is a probability measure with regression function $f$, and lies in $\mathcal P(b,c)$ when $\min(M,\Sigma) \ge 2(4d+1)\sqrt{\kappa^c R}$.
2. Proposition 4, (54): $\mathcal K(\rho_f, \rho_{f'}) \le \frac{16}{15 d L^2}\|\sqrt T(f - f')\|^2$ with $L = 4\sqrt{\kappa^c R}$.
3. Proposition 6 (p. 24): for $m > 16$ there are $N \ge e^{m/24}$ sign vectors in $\{-1,+1\}^m$ with pairwise $\sum_n(\sigma_i^n - \sigma_j^n)^2 \ge m$.
4. Proposition 5 (pp. 22–23): for small $\epsilon$ there are $N_\epsilon \ge e^{\gamma\epsilon^{-1/(bc)}}$ functions in the source class with $\epsilon \le \|\sqrt T(f_i - f_j)\|^2 \le 4\epsilon$.
5. Theorem 5 (p. 24): for every algorithm some $\rho_* \in \mathcal P(b,c)$ has $\mathbb P[\text{excess risk} > \epsilon/4] \ge \min\{N^*/(N^*+1),\ e^{-3/e}\sqrt{N^*}\, e^{-4\ell\epsilon/(15 d\kappa^c R)}\}$, $N^* = e^{\gamma\epsilon^{-1/(bc)}}$.

Significance

The result. Theorem 2 shows that the exponent $bc/(bc+1)$ attained by RLS cannot be improved by any estimator over $\mathcal P(b,c)$ when $\dim Y < \infty$; for $c = 1$ RLS is optimal up to a logarithmic factor. It separates what is a property of the problem class from what is a property of the algorithm: improvements to kernel methods must change the class (stronger assumptions) rather than the rate. The construction — a packing of the source class measured in the $\sqrt T$-norm, combined with a KL bound for an explicit noise model — is the template reused in many later lower bounds for kernel and inverse-problem estimators.

Formalizing it. The result has been proved since 2007; no machine-checked version exists. A formal proof needs the information-theoretic lower-bound machinery (a Fano-type inequality for many hypotheses), a Varshamov–Gilbert-type packing of the Hamming cube, KL divergence between explicit mixtures, and the spectral description of the covariance operator of a vector-valued RKHS. Each of these is reusable well beyond this paper. The mission also makes precise the paper's implicit conventions (see below), which a pen-and-paper reader fills in silently.

Difficulty

The upper half of the argument is not the hard part; the hard part is that the lower bound is uniform over all measurable algorithms, which no direct computation reaches. The step that does not follow from the paper alone is Theorem 5: its proof invokes Lemma 3.3 and Eq. 3.12 of DeVore et al., a Fano-type inequality bounding the probability of correct identification among $N$ hypotheses with pairwise KL divergence at most a given level. That lemma is not stated in the paper and is not in Mathlib; it must be formalized. A second obstacle is the packing (Proposition 6), whose proof is a probabilistic union bound with Hoeffding's inequality. A naive attempt to prove Theorem 2 by exhibiting a single bad distribution fails: for any fixed $\rho$ some algorithm (the constant one returning $f_{\rho}$) has zero excess risk, so the bad distribution must depend on the algorithm, and the order of quantifiers is essential.

Formalization scope

Mathlib's RKHS ℝ H X Y provides the function space with continuous evaluation, RKHS.kerFun is $K_x$, and InformationTheory.klDiv is the Kullback–Leibler information. The operator $T$ is never built as an operator-valued Bochner integral: it is recorded through its quadratic form $\int\langle f(x), g(x)\rangle d\rho_X$ and an eigen-system indexed by $\mathbb N$ from $0$ (the paper's $t_n$ is t (n-1)). The trace in Hypothesis 1 is a series over a Hilbert basis of $Y$; the conditional law in Hypothesis 2 is Measure.condKernel, and the moment integral there is a lintegral. $\mathcal E[f_{\mathcal H}]$ in the goal is $\inf_{f \in \mathcal H}\mathcal E[f]$.

Conventions and added hypotheses, each implicit on the page:

• $\mathcal P(b,c)$ is assumed nonempty; the proof fixes $\rho_0 \in \mathcal P(b,c)$, and over an empty prior the supremum is over the empty set and the statement is false.
• Algorithms are measurable maps $Z^\ell \to \mathcal H$; this is the reading under which the probability in the goal is defined, and it restricts the infimum relative to "all mappings".
• The constants $M, \Sigma, R, \alpha, \beta, \kappa$ are positive; the basis $(v_j)$ of $Y$ in Proposition 4 is orthonormal.
• Proposition 4's "$\|g\|^2 \le R$" for $f'$ is read as $\|g'\|^2 \le R$; Proposition 6's "$i \ne, j$" as $i \ne j$; (56) is required for $i \ne j$. The proof's variance display on p. 22 is wrong for $d \ge 2$, but the conclusion of Proposition 4 holds.

The goal is the $\varepsilon$–$\tau$–$L$ unfolding of the limit and mentions neither $\rho_f$, the KL bound nor the packing, so it cannot be discharged by any of the milestones' constructions in isolation; the distribution $\rho$ is chosen after the algorithm, never before it.

Contributions welcome: a general Fano/DeVore-type lemma for finitely many hypotheses, the Varshamov–Gilbert bound, KL ≤ χ² for finite mixtures, and lemmas relating covForm to the excess risk.

Selected references

• A. Caponnetto, E. De Vito, Optimal rates for the regularized least-squares algorithm, Found. Comput. Math. 7 (2007) 331–368. https://doi.org/10.1007/s10208-006-0196-8
• R. DeVore, G. Kerkyacharian, D. Picard, V. Temlyakov, Approximation methods for supervised learning, Found. Comput. Math. 6 (2006) 3–58. https://doi.org/10.1007/s10208-004-0158-6
• L. Györfi, M. Kohler, A. Krzyżak, H. Walk, A Distribution-Free Theory of Nonparametric Regression, Springer, 2002. https://doi.org/10.1007/b97848
• A. B. Tsybakov, Introduction to Nonparametric Estimation, Springer, 2009. https://doi.org/10.1007/b13794

Motivation

Regularized least squares (RLS, also called kernel ridge regression or Tikhonov regularization) is the simplest learning algorithm built on a reproducing kernel Hilbert space. Given $\ell$ examples $(x_i,y_i)$ drawn independently from an unknown distribution $\rho$, it returns the function in a hypothesis space $\mathcal H$ that minimizes the empirical squared error plus $\lambda$ times the squared norm. It is used in regression, in multi-task learning with vector-valued outputs, and as the reference case for spectral regularization methods.

The basic statistical question is how fast the excess risk of the RLS estimator goes to zero as $\ell$ grows, and how to choose $\lambda=\lambda_\ell$ to get that speed. Caponnetto and De Vito (FoCM 2007) answered it for a family of priors $\mathcal P(b,c)$ described by two numbers: the decay rate $b$ of the eigenvalues of the covariance operator of the input distribution and the regularity $c$ of the target function. Their Theorem 1 gives the upper rate $\ell^{-bc/(bc+1)}$; Theorems 2 and 3 show that no algorithm does better when the output space is finite dimensional. This mission formalizes the upper rate.

Timeline. Cucker and Smale (2002) and De Vito, Caponnetto and Rosasco (2005) gave rates for RLS that did not depend on the eigenvalue decay of $\rho_X$. Zhang (2005) introduced the effective dimension $\mathcal N(\lambda)$ as the complexity measure. Caponnetto and De Vito (authors' copy dated 2006, published 2007) combined it with a source condition to obtain rates that are optimal over $\mathcal P(b,c)$, for vector-valued outputs. Steinwart, Hush and Scovel (2009) and Fischer and Steinwart (2020) later extended the analysis to other norms and to $c<1$.

Setting

The input space $X$ is a Polish space and the output space $Y$ is a real separable Hilbert space. The hypothesis space $\mathcal H$ is a real separable Hilbert space of functions $f:X\to Y$ in which evaluation at each point is continuous. For $x\in X$, $K_x:Y\to\mathcal H$ is the adjoint of evaluation at $x$, so $f(x)=K_x^*f$. Hypothesis 1 adds a measurability condition and a uniform trace bound $\operatorname{Tr}(K_x^*K_x)\le\kappa$.

A distribution $\rho$ on $Z=X\times Y$ has marginal $\rho_X$ and conditional laws $\rho(\cdot\mid x)$. The expected risk of $f\in\mathcal H$ is $\mathcal E[f]=\int\|f(x)-y\|_Y^2\,d\rho$. Hypothesis 2 asks that $\mathbb E\|y\|^2<\infty$, that $\mathcal E$ attains its infimum over $\mathcal H$ at some $f_{\mathcal H}$ (the minimizer of minimal norm is used), and that the noise $y-f_{\mathcal H}(x)$ satisfies a Bernstein moment condition with constants $M,\Sigma$.

The covariance operator $T=\int_XK_xK_x^*\,d\rho_X$ is positive and trace class, with $\langle Tf,f\rangle_{\mathcal H}=\int_X\|f(x)\|_Y^2\,d\rho_X$ and eigen-decomposition $T=\sum_nt_n\langle\cdot,e_n\rangle e_n$, $t_1\ge t_2\ge\dots>0$. The prior $\mathcal P(b,c)$, for $1<b<\infty$ and $1\le c\le2$, consists of the $\rho$ satisfying Hypothesis 2, with $f_{\mathcal H}=T^{(c-1)/2}g$ for some $\|g\|^2\le R$ (source condition), and with $\alpha\le n^bt_n\le\beta$ for all $n$ (eigenvalue decay).

The RLS estimator $f_{\mathbf z}^\lambda$ minimizes $\frac1\ell\sum_i\|f(x_i)-y_i\|_Y^2+\lambda\|f\|_{\mathcal H}^2$ over $\mathcal H$.

Formalization targets

Goal: Theorem 1, $1<b<+\infty$

With $\lambda_\ell=\ell^{-b/(bc+1)}$ and $a_\ell=\ell^{-bc/(bc+1)}$ for $c>1$, and $\lambda_\ell=a_\ell=(\log\ell/\ell)^{b/(b+1)}$ for $c=1$,

$$\lim_{\tau\to\infty}\limsup_{\ell\to\infty}\sup_{\rho\in\mathcal P(b,c)}\mathbb P_{\mathbf z\sim\rho^\ell}\Big[\mathcal E[f_{\mathbf z}^{\lambda_\ell}]-\mathcal E[f_{\mathcal H}]>\tau a_\ell\Big]=0.$$

The statement fixes the rate, not the constants: the threshold $\tau$ absorbs every constant of the prior.

Milestones

1. Proposition 1 iii)–v): the excess risk is $\|\sqrt T(f-f_{\mathcal H})\|^2$; the regularized expected and empirical risks have unique minimizers $f^\lambda=(T+\lambda)^{-1}Tf_{\mathcal H}$ and $f_{\mathbf z}^\lambda=(T_{\mathbf x}+\lambda)^{-1}g_{\mathbf z}$.
2. Proposition 2: a Bernstein inequality for means of i.i.d. Hilbert-space-valued variables.
3. Theorem 4: with probability $\ge1-\eta$,

$$\mathcal E[f_{\mathbf z}^\lambda]-\mathcal E[f_{\mathcal H}]\le3C_\eta\Big(\mathcal A(\lambda)+\frac{\kappa^2\mathcal B(\lambda)}{\ell^2\lambda}+\frac{\kappa\mathcal A(\lambda)}{\ell\lambda}+\frac{\kappa M^2}{\ell^2\lambda}+\frac{\Sigma^2\mathcal N(\lambda)}{\ell}\Big),$$

provided $\ell\ge2C_\eta\kappa\mathcal N(\lambda)/\lambda$ and $\lambda\le\|T\|$, where $C_\eta=32\log^2(6/\eta)$, $\mathcal A(\lambda)=\mathcal E[f^\lambda]-\mathcal E[f_{\mathcal H}]$, $\mathcal B(\lambda)=\|f^\lambda-f_{\mathcal H}\|^2$ and $\mathcal N(\lambda)=\operatorname{Tr}[(T+\lambda)^{-1}T]$. 4. Proposition 3: on $\mathcal P(b,c)$, $\mathcal A(\lambda)\le\lambda^c\|T^{(1-c)/2}f_{\mathcal H}\|^2$, $\mathcal B(\lambda)\le\lambda^{c-1}\|T^{(1-c)/2}f_{\mathcal H}\|^2$ and $\mathcal N(\lambda)\le\frac b{b-1}\beta^{1/b}\lambda^{-1/b}$.

Significance

The result. Theorem 1 says that RLS with $\lambda_\ell$ chosen from $(b,c)$ achieves the rate $\ell^{-bc/(bc+1)}$ uniformly over the prior, and the companion lower bounds show this is the minimax rate for $1<c\le2$ and finite-dimensional $Y$. The rate interpolates between the parametric rate $1/\ell$ (fast eigenvalue decay, smooth target) and slower nonparametric rates, and it identifies the effective dimension, rather than the dimension of $\mathcal H$, as the quantity that governs complexity. Theorem 4 is a non-asymptotic bound of independent use; it is the template for later analyses of spectral regularization, gradient descent with early stopping, and random-feature approximations of kernel methods.

The formalization. The result is proved on paper; no machine-checked proof of a kernel ridge regression rate is known to exist. A formal development requires a Bernstein inequality in Hilbert spaces, spectral calculus for a trace-class operator defined from a measure, and the operator-perturbation argument of Theorem 4. Two printed constants are corrected here: Theorem 4's proof yields $3C_\eta$, not $C_\eta$, and Proposition 3's bound on $\mathcal N(\lambda)$ has $\beta^{1/b}$ in place of $\beta$. The paper also states Theorem 1 for $b=+\infty$; that branch fails for $1\le c<2$ and is not posed.

Difficulty

The obvious argument bounds $\|f_{\mathbf z}^\lambda-f^\lambda\|_{\mathcal H}$ by uniform concentration of $T_{\mathbf x}$ around $T$ and multiplies by $\|\sqrt T\|$. That gives a variance term that ignores the eigenvalue decay of $T$, and hence a rate that does not improve with $b$. The optimal rate needs the variance measured through $\mathcal N(\lambda)$, which requires controlling $\sqrt T(T_{\mathbf x}+\lambda)^{-1}$ in operator norm with high probability. This in turn needs concentration of $(T+\lambda)^{-1/2}(T-T_{\mathbf x})$ in Hilbert–Schmidt norm, and the condition $\ell\gtrsim\mathcal N(\lambda)/\lambda$ under which the empirical operator is close enough to $T$. The noise is unbounded, so only the moment condition (9) is available, and the concentration step must use moment bounds rather than boundedness.

Formalization scope

$\mathcal H$ is Mathlib's RKHS ℝ H X Y, and $K_x$ is RKHS.kerFun H x. The trace in Hypothesis 1 is computed in a fixed Hilbert basis of $Y$. $\rho_X$ is the first marginal and $\rho(\cdot\mid x)$ is condKernel. The operator $T$ is not built as an operator-valued integral. It enters through its quadratic form $\int\langle f(x),g(x)\rangle\,d\rho_X$ and through an eigen-system $(e_n,t_n)$, indexed from $0$, so (17) reads $\alpha\le(n+1)^bt_n\le\beta$. The effective dimension is $\sum_nt_n/(t_n+\lambda)$. Samples are Fin ℓ → X × Y under the product measure, and probabilities of events are outer measures. "With probability at least $1-\eta$" is stated as "the bad event has measure at most $\eta$". The noise condition and the moments of Proposition 2 are integrals of nonnegative functions with values in $[0,\infty]$.

The hypotheses the paper uses but does not display are added: positivity of $M,\Sigma,R,\alpha,\beta,\kappa$, and integrability of the random variable in Proposition 2. The RLS estimator of the goal is any family of minimizers of (18) for $\ell\ge2$; at $\ell=1$ and $c=1$ the parameter $\lambda_1=0$ is degenerate. The goal quantifies uniformly: $L$ is chosen before $\rho$, and the estimator is fixed before $\rho$. The goal must not be replaced by a statement about $\mathcal A,\mathcal B,\mathcal N$ or by Theorem 4's event: it is the uniform rate for the estimator itself.

Reusable infrastructure: Hilbert-space Bernstein inequalities (Proposition 2 alone is a valuable target), spectral calculus for compact positive operators given by a quadratic form, and the representer/normal equation for vector-valued RLS. Contributions of proofs for any milestone, and of supporting lemmas about trace-class operators and effective dimension, are welcome.

Selected references

• A. Caponnetto, E. De Vito, Optimal rates for the regularized least-squares algorithm, Found. Comput. Math. 7 (2007) 331–368. https://doi.org/10.1007/s10208-006-0196-8
• F. Cucker, S. Smale, On the mathematical foundations of learning, Bull. Amer. Math. Soc. 39 (2002) 1–49. https://doi.org/10.1090/S0273-0979-01-00923-5
• E. De Vito, A. Caponnetto, L. Rosasco, Model selection for regularized least-squares algorithm in learning theory, Found. Comput. Math. 5 (2005) 59–85. https://doi.org/10.1007/s10208-004-0134-1
• T. Zhang, Learning bounds for kernel regression using effective data dimensionality, Neural Comput. 17 (2005) 2077–2098. https://doi.org/10.1162/0899766054323008
• I. Pinelis, Optimum bounds for the distributions of martingales in Banach spaces, Ann. Probab. 22 (1994) 1679–1706. https://doi.org/10.1214/aop/1176988477
• I. Steinwart, D. Hush, C. Scovel, Optimal rates for regularized least squares regression, COLT 2009. https://www.cs.mcgill.ca/~colt2009/papers/038.pdf
• S. Fischer, I. Steinwart, Sobolev norm learning rates for regularized least-squares algorithms, J. Mach. Learn. Res. 21 (2020) 1–38. https://jmlr.org/papers/v21/19-734.html

Motivation

Many design problems in control, finance and engineering are convex programs whose constraints depend on an uncertain parameter $\delta$: a solution must satisfy $x\in\mathcal X_\delta$ for every $\delta$ in a possibly infinite set $\Delta$. Enforcing all constraints (robust optimization) is often intractable or overly conservative. The scenario approach draws $N$ independent samples of $\delta$, solves the convex program with those $N$ constraints only, and asks how likely it is that the resulting solution violates a fresh constraint. The question matters wherever a randomized design is certified by a confidence statement, from robust control to chance-constrained portfolio selection.

Timeline.

• Calafiore and Campi (Math. Program. 2005; IEEE TAC 2006) introduced the method and bounded the probability that the violation exceeds $\varepsilon$ by a quantity of order $\binom Nd(1-\varepsilon)^{N-d}$. The bound is valid but loose.
• Campi and Garatti (SIAM J. Optim. 2008, this mission's source) proved the bound $\sum_{i=0}^{d-1}\binom Ni\varepsilon^i(1-\varepsilon)^{N-i}$ for every convex problem satisfying existence and uniqueness of solutions. They showed it is attained with equality by every fully-supported problem, so it cannot be improved without further assumptions.
• Later work extended the result to non-unique solutions, constraint removal, and non-convex decisions (Campi and Garatti, Introduction to the Scenario Approach, SIAM 2018).

Setting

Let $(\Delta,\mathcal D,\mathbb P)$ be a probability space, $c\in\mathbb R^d$ with $d\ge1$, and let $\mathcal X\subseteq\mathbb R^d$ and $\mathcal X_\delta\subseteq\mathbb R^d$ ($\delta\in\Delta$) be convex closed sets. The violation probability of a point $x$ is

$$V(x)=\mathbb P\{\delta\in\Delta:\ x\notin\mathcal X_\delta\}.$$

For a multi-extraction $(\delta^{(1)},\dots,\delta^{(m)})\in\Delta^m$, the program $P_m$ minimises $c^\top x$ over $x\in\mathcal X\cap\bigcap_{i=1}^m\mathcal X_{\delta^{(i)}}$. It is assumed that every $P_m$ has a unique solution $x^*_m$. A constraint $\delta^{(r)}$ is a support constraint of $P_m$ if its removal changes the solution. A convex $P_m$ has at most $d$ support constraints (Proposition 2.2). The problem is fully-supported if, for every $m\ge d$, the program $P_m$ built from $m$ independent samples has exactly $d$ support constraints with $\mathbb P^m$-probability one.

Two further objects carry the argument. For $\mathcal I\subseteq\{1,\dots,m\}$ of cardinality $d$, $S_{\mathcal I}$ is the set of multi-extractions whose support constraints have exactly the indexes in $\mathcal I$. The violation law is

$$F(\alpha)=\mathbb P^d\{V(x^*_d)\le\alpha\},$$

the distribution of the violation of the solution built from $d$ samples.

Formalization targets

Goal: Theorem 2.4, equation (2.3)

For a fully-supported problem, every $N\ge d$ and every $\varepsilon\in[0,1]$,

$$\mathbb P^N\{V(x^*_N)>\varepsilon\}=\sum_{i=0}^{d-1}\binom Ni\varepsilon^i(1-\varepsilon)^{N-i}.$$

Milestones (PART 1 of §3)

• Proposition 2.2: at most $d$ support constraints.
• $S_{\bar{\mathcal I}}\subseteq\widetilde S_{\bar{\mathcal I}}$ for $\bar{\mathcal I}=\{1,\dots,d\}$, where $\widetilde S_{\bar{\mathcal I}}$ is the set where $\delta^{(d+1)},\dots,\delta^{(m)}$ are not violated by the solution generated by $\delta^{(1)},\dots,\delta^{(d)}$; and $\widetilde S_{\bar{\mathcal I}}\subseteq S_{\bar{\mathcal I}}$ up to a probability-zero set.
• (3.3): $\mathbb P^m\{S_{\mathcal I}\}=\int_0^1(1-\alpha)^{m-d}F(\mathrm d\alpha)$ for every $\mathcal I$ of cardinality $d$.
• (3.4): $\binom md\int_0^1(1-\alpha)^{m-d}F(\mathrm d\alpha)=1$ for all $m\ge d$.
• Moment uniqueness: $F(\alpha)=\alpha^d$ is the only distribution on $[0,1]$ satisfying (3.4).
• (3.2): $F(\alpha)=\alpha^d$.
• Partition chain: $\mathbb P^N\{V(x^*_N)>\varepsilon\}=\binom Nd\int_{(\varepsilon,1]}(1-\alpha)^{N-d}F(\mathrm d\alpha)$.
• Integration by parts: $\binom Nd\int_\varepsilon^1(1-\alpha)^{N-d}\,d\,\alpha^{d-1}\,\mathrm d\alpha=\sum_{i=0}^{d-1}\binom Ni\varepsilon^i(1-\varepsilon)^{N-i}$.

Significance

The result. Equation (2.3) shows that the scenario bound (2.2) is tight: no bound that depends only on $N$, $d$ and $\varepsilon$ can be smaller, because a fully-supported problem attains it. The distribution of $V(x^*_N)$ is then a Beta law, $\mathbb P^N\{V(x^*_N)\le\varepsilon\}$ being the probability that a $\mathrm{Binomial}(N,\varepsilon)$ variable is at least $d$, the same for every fully-supported problem. This is what fixes the sample sizes used in practice: $N$ is chosen so that the binomial tail is below a confidence level $\beta$. Fact (3.2), that $V(x^*_d)$ has distribution function $\alpha^d$ whatever the problem, is a distribution-free statement of independent interest.

Formalizing it. The result is proved in the source. As far as is known it has no machine-checked proof. The goal statement is already posed on the platform, and this mission supplies the paper's proof structure as milestones. Two milestones are reusable outside the scenario approach: the uniqueness of a distribution on $[0,1]$ given the moments $\int(1-\alpha)^k\,\mathrm dF=1/\binom{d+k}d$, and the incomplete-beta identity for binomial tails.

Difficulty

The obvious route would compute the law of $V(x^*_N)$ directly, but it depends on the geometry of the constraints. The paper never computes it. It obtains the law of $V(x^*_d)$ only implicitly, through the infinite family of identities (3.4), and recovers it by a uniqueness theorem for moment problems. Two points need care. First, full support holds only almost surely: duplicated samples, for instance, produce programs with fewer than $d$ support constraints, so every set identity holds only up to null sets. Second, the claim that removing a non-support constraint keeps the first $d$ constraints as the only support constraints uses Proposition 2.2. Two identical non-support constraints show that a constraint can become a support constraint after another is removed, unless the count is bounded by $d$.

Formalization scope

Goal. The goal is the already-posed platform statement ScenarioApproach.Generalization.violation_tail_eq_binomial_sum_of_fullySupported (theorem id cffaa932-832c-42ca-9e81-1848ffab7e34), referenced as it stands and not restated. Proposition 2.2 is the platform statement card_support_constraints_le_dim (f70e8aa3-…). This mission adds the PART 1 steps as milestones under ScenarioExact.PartOne.

Representation. Decisions are vectors in EuclideanSpace ℝ (Fin d). A multi-extraction is ω : Fin m → Δ, with 0-based indexes, so $\bar{\mathcal I}$ is $\{i : i<d\}$ and "$\delta^{(d+1)},\dots,\delta^{(m)}$" are the indexes $j\ge d$. $\mathbb P^m$ is Measure.pi (fun _ : Fin m => P). $V$, the feasible set, solutions, support constraints and full support are the published definitions violation, feasibleSet, IsSolution, IsSupportConstraint and FullySupported. A support constraint is one whose removal admits a feasible point of strictly smaller cost, which under uniqueness is the paper's "its removal changes the solution". Full support is almost sure, not pointwise.

Hypotheses made explicit. Assumption 1 is entered as existence and uniqueness of the solution for every number of constraints and every sample, together with a family of solution maps θs k, each assumed to solve $P_k$ and to be measurable. Under uniqueness, θs N is the goal's solution map. The paper's "measurability ... is assumed for granted" (p. 4) is replaced by joint measurability of $\{(x,\delta):x\in\mathcal X_\delta\}$ and measurability of the solution maps, the same two hypotheses as the goal. No set $S_{\mathcal I}$ is assumed measurable. The nonempty-interior clause of Assumption 1 is unused in PART 1 and is not assumed, so the milestones compose with the goal.

Conventions. $F$ is the push-forward measure violationLaw on $\mathbb R$, with $F(\alpha)$ = violationLaw … (Set.Iic α). Integrals against $F$ are lower Lebesgue integrals of nonnegative integrands, as extended nonnegative reals: over $[0,1]$ for $\int_0^1$, and over $(\varepsilon,1]$ for $\int_\varepsilon^1$ in the partition chain, since that integral comes from the event $V>\varepsilon$. The integration-by-parts identity is a real interval integral. Ranges are $1\le d$, $d\le m$, $d\le N$ and $0\le\varepsilon\le1$.

Ruled out. A pointwise "exactly $d$ support constraints for every sample" would be unsatisfiable for many problems (repeated samples) and would trivialise the probabilistic content, so it is not used. Assuming measurability of the event $\{V(x^*_N)>\varepsilon\}$ or of $S_{\mathcal I}$, or the identity $\mathbb P^m\{S_{\mathcal I}\}=\mathbb P^m\{\widetilde S_{\mathcal I}\}$, as a hypothesis would assume part of the conclusion, so none of these is a hypothesis.

Infrastructure. A complete development needs: the support-constraint count (Proposition 2.2, a Helly-type argument), invariance of product measures under coordinate permutations, the change-of-variables formula for push-forward measures, the Hausdorff moment uniqueness theorem on $[0,1]$, and the binomial–incomplete-beta identity. The last two are general results, and contributions of them are welcome independently.

Selected references

• M. C. Campi, S. Garatti, The exact feasibility of randomized solutions of uncertain convex programs, SIAM J. Optim. 19(3) (2008) 1211–1230. https://doi.org/10.1137/07069821X
• G. Calafiore, M. C. Campi, Uncertain convex programs: randomized solutions and confidence levels, Math. Program. 102 (2005) 25–46. https://doi.org/10.1007/s10107-003-0499-y
• G. Calafiore, M. C. Campi, The scenario approach to robust control design, IEEE Trans. Automat. Control 51(5) (2006) 742–753. https://doi.org/10.1109/TAC.2006.875041
• M. C. Campi, S. Garatti, Introduction to the Scenario Approach, SIAM, 2018. https://doi.org/10.1137/1.9781611975444
• A. N. Shiryaev, Probability, 2nd ed., Springer, 1996, Chapter II, §12. https://doi.org/10.1007/978-1-4757-2539-1

Motivation

MAXCUT asks for a partition of the vertices of a weighted graph into two sets that maximizes the total weight of the edges between them. It is one of Karp's original NP-hard problems, so no polynomial-time exact algorithm is expected, and the natural question is how close a polynomial-time algorithm can come to the optimum. Sampling a uniformly random partition already achieves, in expectation, half of the optimal value. For two decades this factor $1/2$ was essentially the best known.

Goemans and Williamson (J. ACM 42(6), 1995) replaced the combinatorial problem by a semidefinite relaxation, solvable in polynomial time by interior point methods, and rounded its solution with a random Gaussian hyperplane. They proved that the resulting cut has expected weight at least $0.878$ times the maximum. The technique founded the use of semidefinite programming in approximation algorithms. Khot, Kindler, Mossel and O'Donnell (SIAM J. Comput. 37(1), 2007) showed that, assuming the Unique Games Conjecture, no polynomial-time algorithm achieves a better constant. Nesterov (Optim. Methods Softw. 9, 1998) extended the rounding analysis to maximizing any positive semidefinite quadratic form over the hypercube, with the constant $2/\pi$.

This mission formalizes the presentation of these results in §6.6 of S. Bubeck, Convex Optimization: Algorithms and Complexity (arXiv:1405.4980v2), pp. 343–347.

Setting

Let $n\ge 0$ and let $A\in\mathbb R^{n\times n}$ be a symmetric matrix with non-negative entries; $A_{i,j}$ is the weight between points $i$ and $j$. The graph Laplacian is $L=D-A$, where $D$ is the diagonal matrix with entries $\sum_{j=1}^n A_{i,j}$. For $x\in\{-1,1\}^n$ the vector $x$ encodes a partition, and MAXCUT is (6.7)

$$\max_{x\in\{-1,1\}^n} x^\top L x .$$

Write $\langle M,X\rangle=\operatorname{Tr}(M^\top X)$ for the Frobenius inner product and $\mathbb S^n_+$ for the symmetric positive semidefinite matrices. Since $x^\top Lx=\langle L,xx^\top\rangle$ and $xx^\top\in\mathbb S^n_+$ has unit diagonal, MAXCUT is bounded above by the SDP relaxation

$$\max\bigl\{\langle L,X\rangle : X\in\mathbb S^n_+,\ X_{i,i}=1,\ i\in[n]\bigr\}.$$

A solution $\Sigma$ of the relaxation is any feasible matrix attaining this maximum. The rounding draws $\xi\sim\mathcal N(0,\Sigma)$, a centered Gaussian vector with covariance $\Sigma$, and outputs $\zeta=\operatorname{sign}(\xi)\in\{-1,1\}^n$ coordinatewise.

Formalization targets

Goal: Theorem 6.11 (Goemans–Williamson)

For $A$ symmetric with non-negative entries, $L=D-A$, $\Sigma$ any solution of the relaxation, $\xi\sim\mathcal N(0,\Sigma)$ and $\zeta=\operatorname{sign}(\xi)$:

$$\mathbb E\,\zeta^\top L\zeta\ \ge\ 0.878\max_{x\in\{-1,1\}^n}x^\top Lx.$$

Milestones

1. Bounded entries. If $\Sigma\in\mathbb S^n_+$ and $\Sigma_{i,i}=1$, then $|\Sigma_{i,j}|\le 1$ (remark in the proof of Lemma 6.12).
2. Lemma 6.12 (Sheppard's formula). If $\xi\sim\mathcal N(0,\Sigma)$ with $\Sigma_{i,i}=1$ and $\zeta=\operatorname{sign}(\xi)$, then $\mathbb E\,\zeta_i\zeta_j=\frac{2}{\pi}\arcsin(\Sigma_{i,j})$.
3. Inequality (6.8). $1-\frac{2}{\pi}\arcsin(t)\ge 0.878(1-t)$ for all $t\in[-1,1]$.
4. Relaxation inequality. $\max_{x}x^\top Lx=\max_x\langle L,xx^\top\rangle\le\langle L,\Sigma\rangle$ for every solution $\Sigma$.

The separately stated Laplacian identity on p. 346 is also included as a theorem item: if $X_{i,i}=1$ for all $i$, then $\langle L,X\rangle=\sum_{i,j}A_{i,j}(1-X_{i,j})$; for $x\in\{-1,1\}^n$, $x^\top Lx=\sum_{i,j}A_{i,j}(1-x_ix_j)$.

Companion: Theorem 6.13 (Nesterov)

For $B\in\mathbb S^n_+$, $\Sigma$ a solution of $\max\{\langle B,X\rangle : X\in\mathbb S^n_+,\ X_{i,i}=1\}$, $\xi\sim\mathcal N(0,\Sigma)$ and $\zeta=\operatorname{sign}(\xi)$:

$$\mathbb E\,\zeta^\top B\zeta\ \ge\ \frac{2}{\pi}\max_{x\in\{-1,1\}^n}x^\top Bx.$$

Significance

The result. Theorem 6.11 is a polynomial-time randomized $0.878$-approximation for MAXCUT: the relaxation is a semidefinite program, and sampling a Gaussian vector and taking signs is cheap. Repeated sampling turns the bound in expectation into a cut of value close to $0.878$ times the optimum with high probability. The same scheme of relaxation followed by randomized rounding underlies approximation algorithms for MAX-2SAT, correlation clustering and quadratic programs over the hypercube, and Nesterov's Theorem 6.13 is the version for an arbitrary positive semidefinite objective.

Formalizing it. Both theorems were proved long ago. To our knowledge neither has a machine-checked proof in Mathlib. The platform has related statements from other books, in different forms: Grothendieck's identity for a standard Gaussian and two unit vectors, and the relaxation guarantee with a Grothendieck constant. This mission states the textbook's results for a Gaussian with a possibly singular covariance matrix, which is the form the rounding uses. A complete development needs Sheppard's formula for a degenerate bivariate Gaussian, an elementary but careful real-variable inequality, and a link between Mathlib's multivariate Gaussian and Gram factorizations of $\Sigma$. All three are reusable.

Difficulty

The algebra (the Laplacian identity and milestone 4) is routine. The probabilistic core is Lemma 6.12. The textbook argument reduces it to the probability that a uniformly random direction separates two unit vectors, which is "a quick picture" on paper. In Lean this requires showing that the pair $(\xi_i,\xi_j)$ has the law of $(\langle V_i,\varepsilon\rangle,\langle V_j,\varepsilon\rangle)$ for a standard Gaussian $\varepsilon$, and then computing an angular measure in the plane, including the degenerate cases $\Sigma_{i,j}=\pm1$, where the pair is supported on a line. A density-based argument fails there, because $\mathcal N(0,\Sigma)$ has no density when $\Sigma$ is singular, and singular solutions of the relaxation occur (for instance $\Sigma=xx^\top$). Inequality (6.8) is a statement about a transcendental function on a closed interval with a tight constant ($0.878$ against the true minimum $\approx0.87856$), so crude estimates do not suffice near the minimizer $t\approx-0.689$.

Formalization scope

• Matrices are Matrix (Fin n) (Fin n) ℝ, vectors Fin n → ℝ. $\mathbb S^n_+$ is Matrix.PosSemidef, which includes symmetry, and $\langle M,X\rangle$ is trace (Mᵀ * X).
• $\mathcal N(0,\Sigma)$ is Mathlib's ProbabilityTheory.multivariateGaussian 0 Σ on EuclideanSpace ℝ (Fin n), defined for every positive semidefinite $\Sigma$, singular ones included. Expectations are Bochner integrals against it, and each theorem also asserts integrability of its (bounded) integrand.
• The sign is $\{-1,1\}$-valued: $\operatorname{sign}(r)=1$ for $r\ge0$ and $-1$ for $r<0$. Mathlib's Real.sign would give $\operatorname{sign}(0)=0$, which takes $\zeta$ out of $\{-1,1\}^n$; the two agree almost surely because $\Sigma_{i,i}=1$.
• The maximum over the hypercube is a finite maximum (Finset.sup') over the $2^n$ Boolean vectors read as $\pm1$ vectors, so it is never a junk value. "The solution" of the relaxation means any maximizer, and maximizers exist since the feasible set is compact and contains the identity.
• Standing hypotheses: in Theorem 6.11, $A$ symmetric with non-negative entries (the book's MAXCUT setting); in Lemma 6.12, $\Sigma$ positive semidefinite (implicit in "$\xi\sim\mathcal N(0,\Sigma)$"); in Theorem 6.13, $B$ positive semidefinite. The identities of milestones 4 and 5 hold for every real matrix $A$ and are stated without hypotheses on $A$.
• Ruled out: tying $\xi$'s law to anything other than $\Sigma$, or dropping optimality of $\Sigma$, would make the goal false or vacuous; here the law is exactly $\mathcal N(0,\Sigma)$ and $\Sigma$ is a maximizer.
• Welcome contributions: Sheppard's formula in Mathlib's multivariate Gaussian language, a proof of (6.8), and the Schur product theorem ($A,B\succeq0\Rightarrow A\circ B\succeq0$) used in Theorem 6.13.

Selected references

• S. Bubeck, Convex Optimization: Algorithms and Complexity, Foundations and Trends in Machine Learning 8(3–4):231–358, 2015. arXiv:1405.4980v2
• M. X. Goemans, D. P. Williamson, Improved approximation algorithms for maximum cut and satisfiability problems using semidefinite programming, J. ACM 42(6):1115–1145, 1995. doi:10.1145/227683.227684
• Yu. Nesterov, Semidefinite relaxation and nonconvex quadratic optimization, Optim. Methods Softw. 9(1–3):141–160, 1998. doi:10.1080/10556789808805690
• S. Khot, G. Kindler, E. Mossel, R. O'Donnell, Optimal inapproximability results for MAX-CUT and other 2-variable CSPs?, SIAM J. Comput. 37(1):319–357, 2007. doi:10.1137/S0097539705447372
• W. F. Sheppard, On the application of the theory of error to cases of normal distribution and normal correlation, Phil. Trans. R. Soc. A 192:101–167, 1899. doi:10.1098/rsta.1899.0003

Motivation

When a problem has millions of variables, even one full gradient can be too expensive to compute, while a single partial derivative $\partial f/\partial x_i$ is often cheap: in regularized regression, support vector machines and many structured problems, updating one coordinate costs a small fraction of a full gradient step. Coordinate descent methods exploit this by moving along one coordinate at a time. They are among the oldest optimization schemes and were for a long time analysed only for cyclic orders and only asymptotically.

Nesterov (2012) showed that choosing the coordinate at random, with probabilities depending on the coordinate-wise smoothness constants, gives global, non-asymptotic rates that can beat full gradient descent in total work. This mission formalizes that analysis as presented in §6.4 of S. Bubeck, Convex Optimization: Algorithms and Complexity (arXiv:1405.4980v2), pp. 338–342, and in particular its linear rate for strongly convex functions (Theorem 6.8).

Setting

Let $f:\mathbb R^n\to\mathbb R$ be differentiable, write $\nabla_i f(x)=\frac{\partial f}{\partial x_i}(x)$ and let $e_i$ be the $i$-th standard basis vector. The function is directionally smooth with constants $\beta_1,\dots,\beta_n>0$ if

$$|\nabla_i f(x+ue_i)-\nabla_i f(x)|\le\beta_i|u|\qquad\text{for all } i\in[n],\ x\in\mathbb R^n,\ u\in\mathbb R,$$

equivalently, each one-variable restriction $u\mapsto f(x+ue_i)$ is $\beta_i$-smooth.

For a real exponent $c$, the weighted norms are

$$\|x\|_{[c]}=\sqrt{\textstyle\sum_{i}\beta_i^{c}x_i^2},\qquad \|x\|^*_{[c]}=\sqrt{\textstyle\sum_{i}\beta_i^{-c}x_i^2}.$$

For $\alpha>0$, $f$ is $\alpha$-strongly convex w.r.t. a norm $\|\cdot\|$ if $f(x)-f(y)\le\nabla f(x)^\top(x-y)-\frac{\alpha}{2}\|x-y\|^2$ for all $x,y$. The point $x^*$ is a minimizer of $f$.

For $\gamma\ge0$, RCD($\gamma$) starts at $x_1\in\mathbb R^n$ and iterates

$$x_{s+1}=x_s-\frac{1}{\beta_{i_s}}\nabla_{i_s}f(x_s)\,e_{i_s},$$

where $i_1,i_2,\dots$ are drawn independently from $p_\gamma(i)=\beta_i^\gamma/\sum_{j}\beta_j^\gamma$. The case $\gamma=0$ is uniform sampling; $\gamma=1$ samples proportionally to $\beta_i$.

Formalization targets

Goal: Theorem 6.8 (p. 341)

Let $\gamma\ge0$, let $f$ be $\alpha$-strongly convex w.r.t. $\|\cdot\|_{[1-\gamma]}$ and directionally smooth with constants $\beta_i$, and let $\kappa_\gamma=\sum_i\beta_i^\gamma/\alpha$. Then for every $t\ge0$

$$\mathbb E f(x_{t+1})-f(x^*)\le\Big(1-\frac{1}{\kappa_\gamma}\Big)^t\big(f(x_1)-f(x^*)\big).$$

Milestones

1. Lemma 6.9 (p. 341): for $f$ $\alpha$-strongly convex w.r.t. any norm, $f(x)-f(x^*)\le\frac{1}{2\alpha}\|\nabla f(x)\|_*^2$.
2. One coordinate step (p. 340): $f\big(x-\frac{1}{\beta_i}\nabla_i f(x)e_i\big)-f(x)\le-\frac{1}{2\beta_i}(\nabla_i f(x))^2$.
3. Expected decrease (p. 340): $\mathbb E_{i_s}f(x_{s+1})-f(x_s)\le-\frac{1}{2\sum_i\beta_i^\gamma}\big(\|\nabla f(x_s)\|^*_{[1-\gamma]}\big)^2$.
4. Lemma 6.9 in the weighted norm (p. 342): $\big(\|\nabla f(x)\|^*_{[1-\gamma]}\big)^2\ge2\alpha(f(x)-f(x^*))$.
5. Contraction (pp. 341–342): one step multiplies the expected gap by at most $1-1/\kappa_\gamma$.

Companion: Theorem 6.7 (pp. 339–340)

For $f$ convex and directionally smooth, and $t\ge2$,

$$\mathbb E f(x_t)-f(x^*)\le\frac{2R_{1-\gamma}^2(x_1)\sum_i\beta_i^\gamma}{t-1},\qquad R_{1-\gamma}(x_1)=\sup_{f(x)\le f(x_1)}\|x-x^*\|_{[1-\gamma]}.$$

Significance

Theorem 6.8 says random coordinate descent converges linearly, with a rate governed by $\sum_i\beta_i^\gamma/\alpha$ instead of the global smoothness constant. For $\gamma=1$, directional smoothness implies $f$ is $\beta$-smooth with $\beta\le\sum_i\beta_i$, so for functions whose global smoothness constant is of the order of $\sum_i\beta_i$, RCD(1) attains the accuracy of gradient descent after the same number of iterations (book, p. 340, comparing Theorem 6.7 with Theorem 3.3), while each iteration touches a single coordinate. The same per-step inequalities underlie later accelerated and parallel coordinate methods.

These results are proved in the literature (Nesterov 2012; Bubeck 2015). Their contribution here is a machine-checked version. As far as a search of the Prove2Me catalogue shows, no coordinate descent rate of this kind has been formalized there; a Euclidean-norm special case of Lemma 6.9 exists on the platform as a separate result, but not the arbitrary-norm lemma or the weighted-norm instance used here.

Difficulty

The main obstacle is bookkeeping of the randomness: the per-step inequality holds for each fixed iterate, while the theorem is about the expectation over the whole sequence of draws $i_1,\dots,i_t$, so the pointwise contraction has to be passed through the tower of conditional expectations. In the strongly convex case this is linear and exact; for Theorem 6.7 the recursion on $\delta_s=\mathbb Ef(x_s)-f(x^*)$ is quadratic, and since $\delta_s$ is an expectation while the gradient norm at $x_s$ is random, the pointwise inequality does not transfer to $\delta_s$ verbatim. A second point is geometric: strong convexity, the dual norm and the sampling distribution must use matching weights ($\beta_i^{1-\gamma}$ against $\beta_i^{\gamma}$), and a mismatch silently changes the constant.

Formalization scope

• $\mathbb R^n$ is EuclideanSpace ℝ (Fin n). The gradient is an explicit map $g$ with HasGradientAt f (g x) x, so $\nabla_i f(x)=g(x)_i$. Lemma 6.9 is stated for a finite-dimensional real normed space with the Fréchet derivative and the operator norm as the dual norm.
• Powers $\beta_i^c$ are real powers. The theorems assume $n\ge1$, $\alpha>0$ and $\beta_i>0$, which the book uses implicitly; $\gamma\ge0$ is the book's.
• RCD(γ) is a deterministic function of the drawn coordinates, and the expectation over $t$ independent draws from $p_\gamma$ is the finite sum $\sum_{(i_1,\dots,i_t)\in[n]^t}\prod_s p_\gamma(i_s)\,F(i_1,\dots,i_t)$. No measure theory or integrability conventions are involved.
• The minimizer $x^*$ is assumed to exist, as the book does throughout; its uniqueness, which the book assumes "only for sake of notation", is not used.
• In Theorem 6.7 the supremum $R_{1-\gamma}(x_1)$ is passed as any real upper bound $R$ on the sublevel set, which is equivalent when the supremum is finite and avoids Lean's value $0$ for an unbounded supremum.
• Directional smoothness is required at every $x$ and $u$, and $p_\gamma$ is fixed by the $\beta_i$; neither is weakened to hold only along the iterates, which would change the theorem.

A complete development needs the one-dimensional descent lemma (3.5), weighted Cauchy–Schwarz for the dual pair $\|\cdot\|_{[c]},\|\cdot\|^*_{[c]}$, and a decomposition of the finite expectation over $[n]^{t+1}$ into the last draw and the first $t$. The weighted-norm and finite-expectation lemmas are reusable for other randomized coordinate and sampling methods. Proofs of the milestones and of either theorem are welcome.

Selected references

• S. Bubeck, Convex Optimization: Algorithms and Complexity, Foundations and Trends in Machine Learning 8(3–4):231–358, 2015. arXiv:1405.4980v2, §6.4, pp. 338–342.
• Yu. Nesterov, Efficiency of coordinate descent methods on huge-scale optimization problems, SIAM Journal on Optimization 22(2):341–362, 2012. doi:10.1137/100802001
• P. Richtárik and M. Takáč, Parallel coordinate descent methods for big data optimization, Mathematical Programming 156:433–484, 2016. arXiv:1212.0873

Motivation

Retailers of seasonal goods (fashion, electronics, airline seats) sell a fixed stock over a finite season to customers who arrive over time, and those customers know that prices may fall. A customer who expects a markdown waits, and a seller who ignores this loses revenue. The classical revenue-management literature (Gallego and van Ryzin 1994; Talluri and van Ryzin 2004) models myopic customers who buy on arrival or leave; the literature on forward-looking (strategic) buyers, for instance Aviv and Pazgal (2008), studies particular price paths.

Board and Skrzypacz ask the mechanism-design question: among all selling schemes, which maximizes the seller's expected discounted revenue when buyers arrive over time, have private values and time their purchases strategically? Their answer, published in the Journal of Political Economy in 2016, is that the optimal mechanism has a simple structure: in every period the seller sells to the highest remaining buyer if and only if his value exceeds a cutoff that depends only on the period and the number of units left. When demand is weakly decreasing over time, the cutoffs are characterized by one-period indifference conditions, which in the continuous-time limit can be implemented by posted prices. The source used here is the authors' accepted manuscript of February 6, 2015; all page numbers refer to that manuscript.

Setting

A seller has units of a good and sells them over periods $t\in\{1,\dots,T\}$; unsold units are worth zero after period $T$. Payoffs are discounted by $\delta\in(0,1)$. At the start of period $t$ a random number $N_t$ of buyers arrives, independently across periods, with a law that may depend on $t$. Each buyer wants one unit; his value is drawn independently from a distribution with continuous density $f$, distribution function $F$ and support $[\underline v,\bar v]$. The marginal revenue of a buyer with value $v$ is

$$m(v)=v-\frac{1-F(v)}{f(v)},$$

assumed strictly increasing and continuously differentiable, with $m(\underline v)<0$.

By the standard mechanism-design reduction (§2.1, eq. (2.5)), the seller's problem is to choose when to serve each buyer so as to maximize the expected discounted sum of the served buyers' marginal revenues. The state in period $t$, after the period-$t$ entrants have arrived, is the number $k$ of units left and the values $y^1\ge y^2\ge\cdots$ of the buyers present. The value $\Pi^k_t$ and the pre-entry value $\tilde\Pi^k_{t}$ satisfy the Bellman equation (4.3):

$$\Pi^k_t(\mathbf y)=\max_{0\le j\le k}\Big[\sum_{i=1}^j m(y^i)+\delta\,\tilde\Pi^{k-j}_{t+1}(\mathbf y^{-j})\Big],\qquad \tilde\Pi^k_{t+1}(\mathbf y)=E_{t+1}\big[\Pi^k_{t+1}(\mathbf y\cup\mathbf v_{t+1})\big],$$

where $\mathbf y^{-j}$ is the set of buyers left after the $j$ highest are served and $\mathbf v_{t+1}$ the next period's entrants. Selling one unit to $y^1$ today rather than none gives the difference function

$$\Delta\Pi^k_t(y^1,\mathbf y^{-1})=m(y^1)+\delta\tilde\Pi^{k-1}_{t+1}(\mathbf y^{-1})-\delta\tilde\Pi^k_{t+1}(y^1,\mathbf y^{-1}),$$

and the cutoff $x^k_t$ is the smallest $y\in[\underline v,\bar v]$ with $\Delta\Pi^k_t(y,\varnothing)\ge 0$. Comparing selling to $y^1$ today with waiting and selling at least one unit tomorrow (to the best of $y^1$ and the entrants) gives $D\Pi^k_t(y^1)$ (p. 17). Demand is weakly decreasing in the usual stochastic order if $P(N_{t+1}>x)\le P(N_t>x)$ for all $x$ and $t$.

Formalization targets

Goal: Theorem 2 (p. 17)

If $N_t$ is weakly decreasing in the usual stochastic order then, for every $k\ge 1$,

$$x^k_{t+1}\le x^k_t\quad(1\le t\le T-1),\qquad D\Pi^k_t(x^k_t)=0,$$

and $x^k_t$ is the unique root of $D\Pi^k_t$ in $[\underline v,\bar v]$ for $t\le T-1$: the seller is indifferent between selling to the cutoff type today and waiting one period to sell that unit tomorrow (the one-period-look-ahead property).

Central milestone: Theorem 1 (p. 15)

For every $t$ and $k\ge 1$, the optimal rule sells to the highest buyer iff $y^1\ge x^k_t$, whatever the values of the lower buyers; $x^{k+1}_t\le x^k_t$; and $x^k_t$ is the unique root of $\Delta\Pi^k_t$.

Milestones

In attack order:

1. Lemma 1: allocations are monotone in values.
2. Lemma 2: with cutoffs decreasing in the unit index, units can be treated one at a time.
3. Equation (A.1): increasing differences of $\Pi$.
4. Lemma 3: $\Delta\Pi$ is independent of lower buyers, continuous and strictly increasing in $y^1$, and increasing in $k$.
5. Footnote 12: the boundary values of $\Delta\Pi$.
6. Theorem 1.
7. Strict monotonicity of $D\Pi$ in $y^1$ (p. 18).
8. Lemma 4: $D\Pi^k_{t+1}\ge D\Pi^k_t$.

After the goal, (4.7) gives the period-$(T-1)$ cutoff equation $m(x^k_{T-1})=\delta E_T[\max\{m(x^k_{T-1}),m(v^k_T)\}]$.

Significance

Theorem 1 says that the optimal allocation does not depend on how many buyers are present or what their values are, only on time and inventory. This is what makes the optimal mechanism implementable without eliciting values from buyers as they arrive. Theorem 2 turns the global dynamic program into local indifference conditions. In the continuous-time limit (§5 of the paper) these become differential equations, and the optimum is implemented by posted prices with an auction at the end of the season. Under weakly decreasing demand, therefore, the classical revenue-management practice of posting prices loses nothing against the best possible mechanism.

The paper's results are proved, in prose, with envelope-theorem and coupling arguments. They have not been machine-checked. A formal development would produce a verified backward-induction model of multi-unit dynamic allocation with random arrivals, with the structural results (monotonicity, deterministic cutoffs, monotone comparative statics in inventory and time) that recur across dynamic pricing and optimal stopping. Two printed gaps are recorded below: the positivity of $m(\bar v)$, and the restriction of footnote 12 to $t\le T-1$.

Difficulty

The obvious argument fails at "deterministic". A priori the cutoff for the highest buyer depends on the values of the lower buyers, because selling a unit today changes which of them will be served later and when. Lemma 3(a) holds only under the induction hypothesis that all future cutoffs are already deterministic and decreasing in inventory, so Lemma 3, Theorem 1 and (A.1) form a single backward induction over periods and units, and none of them can be proved in isolation. The value function is an expectation, over a random number of i.i.d. entrants, of a maximum over sorted values, so continuity and strict monotonicity in $y^1$ (Lemma 3(b), and the same for $D\Pi$) are not available from general facts. The tempting argument for Theorem 2, that cutoffs fall over time simply because fewer buyers arrive later, is incomplete: Lemma 4 has to compare two periods with different arrival laws and different future cutoffs at once.

Formalization scope

Periods are natural numbers $1,\dots,T$ with $T\ge1$, and units are natural numbers. Values, marginal revenues and profits are real numbers. The buyers present form a finite multiset of reals, and an absent buyer is absent, never a value $0$. The value law is the measure with density $f$; the density is positive and continuous on $[\underline v,\bar v]$ and zero outside, and $m$ is defined from $f$ and $F$. Expectations over a cohort are lower Lebesgue integrals against $\sum_n P(N_t=n)\,\mu^{\otimes n}$ of nonnegative bounded quantities, so no non-measurable or non-integrable integrand can silently become $0$. The value function is defined by the Bellman equation (4.3), with $\Pi_{T+1}\equiv0$. The sequence problem (4.1) over purchase times is not formalized; the paper says either may be used (footnote 16).

Standing assumptions and handled gaps:

• $\delta\in(0,1)$;
• $N_t$ independent across periods (only the marginal laws enter);
• $m$ strictly increasing and $C^1$ on $[\underline v,\bar v]$ with $m(\underline v)<0$;
• added: $m(\bar v)>0$, which footnote 12 uses without stating; without it no unit is ever sold and no cutoff exists;
• added: $f>0$ on the closed support, needed for $m$ to be defined there;
• footnote 12's equality $\Delta\Pi^k_t(\bar v)=(1-\delta)m(\bar v)$ is stated for $t\le T-1$ only, since $\Delta\Pi^k_T=m$;
• $D\Pi^k_t$ is used only for $t\le T-1$, and Lemma 4 needs $t+1\le T-1$;
• "decreasing" and "increasing" are weak except in Lemma 3(b) and for $D\Pi$;
• at $y^1=x^k_t$ both selling and waiting are optimal.

The cutoff is defined from $\Delta\Pi$, never as the threshold of an optimal policy, and "optimal" always means maximal in (4.3) over every number of units sold. A formalization that postulates a threshold policy, replaces the random cohort by its mean, or sets absent buyers to value $0$ would trivialize or change the statements and is ruled out. The mechanism-design reduction (IC/IR to (2.5)) and the continuous-time results of §5 are out of scope.

The usual stochastic order is the published platform definition StochasticOrders.Usual.UsualOrder. A complete development needs finite-horizon dynamic programming over multisets, expectations of functions of sorted i.i.d. samples, envelope arguments, and monotone coupling for the usual stochastic order on $\mathbb N$; these parts are reusable beyond this mission. Proofs of any milestone are welcome, as is a proof that (4.3) agrees with the sequence problem (4.1).

Selected references

• S. Board and A. Skrzypacz, Revenue Management with Forward-Looking Buyers, Journal of Political Economy 124(4), 2016. https://doi.org/10.1086/686713
• R. B. Myerson, Optimal Auction Design, Mathematics of Operations Research 6(1), 1981. https://doi.org/10.1287/moor.6.1.58
• G. Gallego and G. van Ryzin, Optimal Dynamic Pricing of Inventories with Stochastic Demand over Finite Horizons, Management Science 40(8), 1994. https://doi.org/10.1287/mnsc.40.8.999
• Y. Aviv and A. Pazgal, Optimal Pricing of Seasonal Products in the Presence of Forward-Looking Consumers, Manufacturing & Service Operations Management 10(3), 2008. https://doi.org/10.1287/msom.1070.0183
• K. T. Talluri and G. J. van Ryzin, The Theory and Practice of Revenue Management, Springer, 2004. https://doi.org/10.1007/b139000
• M. Shaked and J. G. Shanthikumar, Stochastic Orders, Springer, 2007. https://doi.org/10.1007/978-0-387-34675-5