1 Introduction↩︎

Consider a stochastic convex optimization problem, where for a collection \(\{\loss_z\}_{z \in \mc{Z}}\) of convex functions, we wish to minimize the population objective \[\begin{align} \label{eqn:objective} \poploss_P(x) \defeq \E_P[\loss_Z(x)] = \int \loss_z(x) dP(z). \end{align}\tag{1}\] Upon being asked what it means to solve such a problem, an optimizer might say “finding a point near the minimizer” or “with approximately optimal function value” [1]–[4]. Yet these notions of optimality can become unsatisfactory: consider finding a median of a distribution \(P\), i.e., \(x \in \R\) such that \[\begin{align} \half \le P(Z \le x) ~~ and ~~ \half \le P(Z \ge x) ~~~ for~ Z \sim P. \end{align}\] Finding a median is equivalent to minimizing the expected absolute error \(\loss_z(x) = |x - z|\).

Two examples highlight issues with these optimality criteria here. First, even if the median \(x\opt\) is unique, seeking a point near \(x\opt\) may be a fool’s errand.

Faced with this difficulty, our hypothetical optimizer instead seeks small excess objective:

These examples show that proximity to a minimizer and small excess risk fail to capture the distributional properties defining a median. In this case, the relevant criterion becomes stationarity for the absolute loss: for \(\loss_z(x) = |x - z|\) we have \(\partial \loss_z(x) = \{\sign(x - z)\}\), where \(\sign(0) = [-1,1]\), giving population subdifferential \[\partial \poploss_P(x) = P(x > Z) - P(x < Z) + [-1, 1] \cdot P(Z = x).\] Therefore, for any \(x \in \R\), the subdifferential \(\partial \poploss_P(x)\) contains \(g\) with \(|g| \le 2\epsilon\) if and only if \[\begin{align} P(Z \ge x) \ge \half - \epsilon ~~ and ~~ P(Z \le x) \ge \half - \epsilon, \end{align}\] addressing the issues Examples 1 and 2 highlight.

More generally, we consider finding near stationary points of convex losses. Let \(X \subset \R^d\) be closed convex and \(\{\loss_z\}_{z \in \mc{Z}}\) be a collection of convex functions \(\loss_z : \R^d \to \R \cup \{+\infty\}\), where \(\relint \dom \loss_z \supset X\), so that each is subdifferentiable over \(X\). For a distribution \(P\) on \(Z \in \mc{Z}\), recall 1 the population loss \(\poploss_P(x) = \E_P[\loss_Z(x)]\). Then the standard first-order conditions for optimality [1] state that \(x\opt\) minimizes \(\poploss_P(x)\) over \(x \in X\) if and only if \[0 \in \partial \poploss_P(x\opt) + \normalcone_X(x\opt) ~~~ where ~~~ \normalcone_X(x) = \{v \mid \<v, y - x\> \le 0 ~ for~all~ y \in X\}\] denotes the normal cone to \(X\) at \(x\). We therefore seek estimators achieving \(\epsilon\)-stationarity, that is, \(\what{x}\) for which \[\label{eqn:stationary-point-goal} \dist\left(0, \partial \poploss_P(\what{x}) + \normalcone_X(\what{x})\right) \le \epsilon.\tag{2}\] We develop estimators \(\what{x}_n\) based on a sample \(Z_1, \ldots, Z_n \simiid P\) satisfying \[\dist\left(0, \partial \poploss_P(\what{x}_n) + \normalcone_X(\what{x}_n) \right) \to 0 ~~ in probability~as~ n \to \infty.\]

1.1 Motivation↩︎

The problem of finding \(\epsilon\)-stationary points 2 , and the more abstract variants we consider later of finding near zeros of maximal monotone operators and variational inequalities, arise in many areas of statistics and optimization. For us, a main motivation comes from a direction within statistics and machine learning that seeks “distribution-free” predictive inference [6]–[12]: algorithms for prediction problems that output prediction sets in which true outcomes are guaranteed to lie with a prescribed probability, without making any assumptions on the underlying data generating distribution. The somewhat extended examples here show how these guarantees are equivalent to obtaining the approximate stationary guarantee 2 .

With these examples as motivation, we see that we truly wish to obtain (approximately) stationary points \(\what{x}\): not points that approximately minimize \(\poploss_P\), or are close to a minimizer, but truly satisfy \(\dist(0, \partial \poploss_P(\what{x})+\normalcone_X(\what{x})) \to 0\). We thus return to the general problem of finding an \(\epsilon\)-stationary point of a convex objective given a sample \(Z_1, \ldots, Z_n \simiid P\).

1.2 Problem framing and background↩︎

1.2.1 Monotone operator formulation↩︎

The stationarity target 2 fits within the more general setting of monotone operators and variational inequalities [15], [16]. A set-valued operator \(A:\R^d \toto \R^d\) is monotone on \(X \subset \R^d\) if every \(x, y \in X\), \(u \in A(x)\), and \(v\in A(y)\) satisfies \(\< u-v, x-y\> \geq 0\), with \(\dom A = \{x \mid A(x) \neq \emptyset\}\) and graph \(\gr A = \{(x, u) \mid u \in A(x)\}\). It is maximally monotone if no other monotone operator has graph strictly containing \(\gr A\). When \(F\) is a proper closed convex function, its subdifferential \(\partial F\) is maximal monotone [15] Theorem 12.17. The more general framing for problem 2 , and that in which we state our results, thus follows: let \(X \subset \R^d\) be a closed convex set. Then for a maximal monotone operator \(A\), we seek \(x\opt \in X\) for which there exists \(u \in A(x\opt)\) such that \(\<u, x - x\opt\> \ge 0\) for all \(x \in X\). For the normal cone \(\normalcone_X(x) = \{v \mid \<v, y - x\> \le 0, ~ all~y \in X\}\), we equivalently [15] seek \(x\opt\) satisfying \[\begin{align} 0 \in A(x\opt) + \normalcone_X(x\opt). \end{align}\] Within deterministic optimization, Rockafellar [17] shows that for maximal monotone operators, the proximal point algorithm in quite some generality generates iterates \(x_k\) satisfying \(\dist(0, A(x_k) + \normalcone_X(x_k)) \to 0\).

The statistical setting appears to introduce substantial difficulty. We assume we have a collection of maximal monotone mappings \(A_z(\cdot) : \R^d \toto \R^d\) for \(z \in \mc{Z}\). We will not overly concern ourselves with measurability, but follow standard practice [18] to assure that integrals have appropriate definitions. The standard Aumann integral [19] defines \[A_P(x) \defeq \int A_z(x) dP(z) = \left\{\int a_z(x) dP(z) \mid a_z(x) ~ is~a~measurable~selection~of~ A_z(x)\right\},\] meaning that \(z \mapsto a_z(x)\) is measurable and \(a_z(x) \in A_z(x)\) for \(P\)-almost-all \(z\). We assume \(A_P\) is maximal monotone over \(X\), that is, \(A_P + \normalcone_X\) is maximal monotone. (When \(A_z(x) = \partial \loss_z(x)\) is the subdifferential of a convex function, these conditions hold essentially automatically [20].)

1.2.2 The failure of uniform subdifferential convergence↩︎

Of course, we only have sample access to \(P\) via \(Z_i \simiid P\). Letting \(P_n = \frac{1}{n} \sum_{i = 1}^n \pointmass_{Z_i}\) denote the empirical distribution and \[\begin{align} \dhaus(X, Y) = \max\Big\{\sup_{x \in X} \inf_{y \in Y} \ltwo{x - y}, \sup_{y \in Y} \inf_{x \in X} \ltwo{x - y}\Big\} \end{align}\] denote the Hausdorff distance between sets \(X\) and \(Y\), the most natural notion [21], [22] of convergence of \(A_{P_n}\) to \(A_P\) would be that \[\label{eqn:uniform-convergence-subdifferentials} \sup_{x\in X}\dhaus(\popop(x),\empop(x)) \to 0\tag{3}\] (in probability, expectation, or almost surely). If we had such convergence, the empirical solution \(\what{x}_n \in X\) satisfying \(0 \in \empop(\what{x}_n) + \normalcone_X(\what{x}_n)\) would satisfy \(\dist(0, \popop(\what{x}_n) + \normalcone_X(\what{x}_n)) \to 0\). (Indeed, let \(w_n \in \normalcone_X(\what{x}_n)\) and \(a_n \in \popop(\what{x}_n)\) satisfy \(0 = a_n + w_n\). Then there is \(a_n\opt \in \popop(\what{x}_n)\) satisfying \(\norms{a_n\opt - a_n} \to 0\), and so \(a_n\opt + w_n = a_n + w_n + o(1) \to 0\).)

Unfortunately, the convergence 3 cannot generally hold. In all dimensions \(d \ge 2\), even subdifferentials of convex functions \(\loss_z : \R^d \to \R\) exhibit quite spectacular failures of convergence [22]. Indeed, [22] exhibit a collection of \(1\)-Lipschitz convex functions \(\loss_z : \R^2 \to \R\) and points \(x_n\), which are functions of the sample \(P_n\), satisfying that \[\begin{align} \dhaus(\partial \poploss_{P_n}(x_n), \partial \poploss_P(x_n)) \ge \half ~~~ with~probability~1 \end{align}\] for all sample sizes \(n\). In Section 1.3, we give a simple example, motivated by theirs, for which an empirical minimizer satisfies the same failure.

1.2.3 Existing stationarity guarantees↩︎

These difficulties inspire several alternative approaches. One restricts the class of functions under consideration. For example, [21] show that if the population objective \(\poploss_P\) is differentiable (which, by convexity, implies continuous differentiability), then uniform convergence 3 holds. [23] explores a related approach, showing that classes of functions whose subdifferentials have finite Vapnik-Chervonenkis dimension enjoy uniform convergence, even in weakly convex settings. In this paper, we make no structural assumptions on the function class \(\{\loss_z\}_{z \in \mc{Z}}\), and so cannot achieve uniform convergence 3 , but show how a randomized proximal-point sampling method achieves stationarity 2 .

The other main approach focuses on weaker guarantees. These rely on one of several essentially equivalent conditions, which roughly show that an algorithm generates iterates near approximately stationary points. This viewpoint follows classically from Ekeland’s variational principle [24] and the Brøndsted-Rockafellar theorem [25]; in keeping with our development, we describe the approach via proximal operators.

For \(\lambda > 0\) and a convex function \(F\) or, respectively, maximal monotone operator \(A\), let \[\begin{align} \proxmap_{\lambda F}(x) \defeq \argmin_y \left\{\lambda F(y) + \half \ltwo{y - x}^2 \right\} ~~ and ~~ \proxmap_{\lambda A}(x) \defeq (I + \lambda A)^{-1}(x). \end{align}\] Standard calculations show that if \(x_\lambda = (I + \lambda A)^{-1}(x_0)\), i.e. \(\frac{1}{\lambda} (x_0 - x_\lambda) \in A(x_\lambda)\), then \(\frac{1}{\lambda} \ltwo{x_\lambda - x_0} \le \epsilon\) implies that \(\dist(0, A(x_\lambda)) \le \epsilon\). Conversely, if \(\dist(0, A(x_0)) \le \epsilon\), then using the monotonicity of \(A\) and letting \(v_0 \in A(x_0)\) satisfy \(\ltwo{v_0} \le \epsilon\), the containment \(\frac{1}{\lambda} (x_0 - x_\lambda) = v_\lambda \in A(x_\lambda)\) implies that \[\begin{align} \ltwo{x_0 - x_\lambda}^2 = \lambda \<x_0 - x_\lambda, v_\lambda\> \le \lambda \<x_0 - x_\lambda, v_0\> \le \lambda \ltwo{x_0 - x_\lambda} \epsilon ~~~ so ~~~ \frac{1}{\lambda} \ltwo{x_0 - x_\lambda} \le \epsilon. \end{align}\] Relatedly, the Moreau envelope \(F_\lambda(x) \defeq \inf_y \{F(y) + \frac{1}{2\lambda} \ltwo{x - y}^2\}\) of \(F\) satisfies \[\begin{align} \nabla F_\lambda(x) = \frac{1}{\lambda} (x - \proxmap_{\lambda F}(x)). \end{align}\] (See [26] or [16].) So if \(\ltwo{\nabla F_\lambda(x)} \le \epsilon\) then \(x^+ = \proxmap_{\lambda F}(x)\) satisfies \[\ltwos{x^+ - x} \le \lambda \epsilon, ~~ \dist(0, \partial F(x^+)) \le \epsilon, ~~ and ~~ F(x^+) \le F(x),\] and conversely, if \(\dist(0, \partial F(x)) \le \epsilon\) then \(\ltwo{\nabla F_\lambda(x^+)} \le \epsilon\).

These implications motivate a substantial literature in stochastic and non-stochastic optimization, with a special focus on weakly convex functions, where researchers have provided quantitative guarantees on the rates of convergence of \(\lambda^{-1} \ltwos{x_k - \proxmap_{\lambda \poploss_P}(x_k)}\) of iterates \(x_k\) from, for example, Gauss-Newton and prox-linear methods, stochastic gradient methods, and stochastic model-based minimization [27]–[30]. [31] extend these techniques to demonstrate the empirical convergence \[\sup_{x \in X} \ltwo{\nabla (\poploss_P)_\lambda(x) - \nabla (\poploss_{P_n})_\lambda(x)} = O\left(\sqrt{\frac{d}{n} \log \frac{n \cdot \diam(X)}{\lambda}}\right),\] of the Moreau envelopes of \(\poploss_P\) and \(\poploss_{P_n}\) with high probability, which in turn implies [31] the graphical convergence that \(\dhaus(\gr \partial \poploss_{P_n}, \gr \partial \poploss_P) \to 0\). Related approaches [18], [32] show convergence of iterates of stochastic approximation-type methods to stationary points, that is, they show results of the form \(\dist(x_k, \popop^{-1}(0)) \to 0\) in quite general settings. These guarantees differ fundamentally from the stronger notion 2 of stationarity we seek, which requires procedures to return a point with small stationary residual: they only ensure that a neighborhood of a chosen point contains a nearly stationary point. As the absolute value \(F(x) = |x|\) makes clear, even for \(1\)-Lipschitz functions \(\ltwo{x - \proxmap_F(x)}\) can be arbitrarily small while \(\ltwos{\nabla F(x)} = 1\).

When the random operators \(A_z\) are single-valued and appropriately Lipschitz continuous, additional guarantees are possible. Particular references include [33] and [34], who provide convergence guarantees for stochastic mirror-prox and operator extrapolation methods. Similarly, [35] provide several results for finding stationary points of smooth convex functions, with finite-sample guarantees so long as \(\nabla \loss_z\) is Lipschitz continuous for each \(z\). Our results assume none of this structure.

1.3 An exemplar challenge in finding nearly stationary points↩︎

To highlight the difficulty of finding stationary points, we find it instructive to exhibit convex functions for which empirical risk minimizers are non-stationary at all sample sizes. We take [22] as motivation and construct an example with the property that \(\poploss_{P_n}\) has an empirical minimizer whose population gradient is large with probability 1.

Consider a countable collection \(\{u_j\}_{j \in \N} \subset \sphere^{d-1}\) satisfying the non-alignment condition \[\begin{align} \label{eqn:non-alignment} \Delta_i \defeq \sup\left\{ \<u_i, u_j\> \mid j \neq i\right\} < 1. \end{align}\tag{4}\] (See Appendix 6.2 for a construction.) Then for \(z \in \{0, 1\}^\N\), define the \(1\)-Lipschitz function \[\label{eqn:bad-losses} \loss_z(x) \defeq \sup_{j \in \N} z_j \hinge{\<x, u_j\> - 1},\tag{5}\] which is convex, nonnegative, and \(z\)-measurable. Empirical samples of these losses fail to have uniform convergence of the subdifferential and have highly non-stationary minimizers:

Because the functions 5 are \(1\)-Lipschitz, the non-stationarity Proposition [proposition:bad-losses] exhibits is essentially as bad as possible. Moreover, even though along each line \(\{t u\}_{t \in \R}\), \(u \in \sphere^{d-1}\) the directional derivatives exhibit uniform convergence \(\sup_{t \in \R} |\poploss_{P_n}'(t u; u) - \poploss_P'(t u; u)| \cp 0\), this convergence fails along many directions. This suggests some of the difficulties in showing that “standard” algorithms, such as regularized empirical minimization or proximal point methods, could generate iterates satisfying \(\dist(0, \nabla \poploss_P(\what{x}_n)) \cp 0\): any such analysis (we have tried many) would need to rely on structural properties of convex functions beyond those that are, at least to us, apparent. Accordingly, we take a different approach.

2 The main result↩︎

Our main technical result relies on the proximal point mapping for the sampled objective \(\poploss_{P_n}\). In brief, we show that the set of points \(x \in X\) for which the proximal point update \[\begin{align} x_\lambda \defeq \proxmap_{\lambda \poploss_{P_n}}(x) = \argmin_{y \in X} \left\{\poploss_{P_n}(y) + \frac{1}{2 \lambda} \ltwo{y - x}^2 \right\} \end{align}\] falls into the set of \(\epsilon\)-stationary points has strictly positive (Lebesgue) measure for any \(\lambda > 0\). More precisely, denote the set of \(\epsilon\)-stationary points for the population objective \(\poploss_P\) by \[\begin{align} S_P(\epsilon) \defeq \left\{x \in X \mid \dist(0, \partial \poploss_P(x) + \normalcone_X(x)) \le \epsilon \right\}. \end{align}\] Then as a corollary of Theorem 1 to follow, we have the following convergence result, where \(\mathsf{Leb}\) denotes (\(d\)-dimensional) Lebesgue measure.

Corollary 1 leads to a simple and naive convergent algorithm: assume we have two independent samples \(P_n\) and \(P_n'\), drawn i.i.d. \(P\). Draw \(u^{(i)}\) uniformly from \(X + 4 \epsilon \ball\), \(i = 1, \ldots, m\), and set \(x^{(i)} = \proxmap_{\poploss_{P_n}}(u^{(i)})\). Then report the index \(i\) minimizing the empirical subdifferential on the second sample, \[\begin{align} \what{x}_n = x^{(\what{i})} ~~~ for ~~~ \what{i} \defeq \argmin_i \dist\left(0, \partial \poploss_{P_n'}(x^{(i)}) + \normalcone_X(x^{(i)}) \right). \end{align}\] Because empirical subdifferentials converge pointwise (see Lemma 4 to come), a quick argument via the union bound shows that, so long as \(1 \ll m \ll n\), we have \[\begin{align} \dist\left(0, \partial \poploss_P(\what{x}_n) + \normalcone_X(\what{x}_n)\right) \cp 0. \end{align}\] In the next section, we provide a slightly more sophisticated variant of this procedure, filling out the details of the argument.

2.1 The stochastic monotone operator case↩︎

Our results extend beyond the convex function case to more general (stochastic) maximal monotone operators, which provide a natural setting that coincides with the analytic techniques to prove the results. As in Section 1.2, assume we have a collection of maximal monotone mappings \(A_z(\cdot) : \R^d \toto \R^d\) for \(z \in \mc{Z}\) and \(A_P(x) = \E_P[A_Z(x)]\) for \(Z \sim P\). We assume \(A_P\) is maximal monotone over \(X\), that is, \(A_P + \normalcone_X\) is maximal monotone. We wish to find \(x \in X\) solving \[\label{eqn:vi} 0 \in A_P(x) + \normalcone_X(x),\tag{6}\] which we assume has solutions \(X\opt = \{x \in X \mid 0 \in A_P(x) + \normalcone_X(x)\}\).

We require moment bounds on \(A_z(x)\), and so say that \(A_z\) has \(p \ge 2\) moments near \(x^0\) if there exists \(b > 0\) and envelope \(\lipconst(z) < \infty\) such that for \(\ltwo{x - x^0} \le b\), each (measurable) selection \(a_z(x) \in A_z(x)\) satisfies \(\ltwo{a_z(x)} \le \lipconst(z)\) and \(\E[\lipconst^p(Z)] < \infty\). Summarizing, to state our results, we make the following standing assumptions:

The subdifferential of a convex objective provides the standard example satisfying Assumption 1. In this case, the subdifferential \(A_z(x) = \partial \loss_z(x)\), and the measurability assumptions and subdifferential equality \(\partial \poploss_P(x) = \E_P[\partial \loss_z(x)]\) follow essentially immediately, so that maximal monotonicity holds [15], [20], [36] and part [item:vi-max-monotone] follows. Part [item:vi-solution] trivially holds if \(\poploss_P\) has a minimizer over \(X\), and then so long as \(\loss_z\) is \(\lipconst(z)\)-Lipschitz continuous near \(x^0\) with \(\E[\lipconst^2(Z)] < \infty\), part [item:vi-moments] also holds. Other examples include the situation that \(X \subset \interior \dom A_z\) for each \(z\) and \(P\) has finite support, which guarantees \(A_P\) is maximal monotone [37]. If \(X\) is compact, then \(A_P + \normalcone_X\) is surjective [16] Corollary 21.25 so 6 has a solution.

With this assumption, we give the natural generalization of stationary points, letting \[\begin{align} S_P(\epsilon) = \{x \in X \mid \dist(0, A_P(x) + \normalcone_X(x)) \le \epsilon\} \end{align}\] denote the \(\epsilon\)-stationary points in the problem 6 . Let \(x^0\) denote the minimum norm solution to 6 , as in Assumption 1.[item:vi-moments]. We then have the following theorem.

A simple corollary is cleaner to state, though the more general statement in Theorem 1 allows us to demonstrate that an algorithm that localizes around empirical minimizers converges in probability.

2.2 Proof sketch↩︎

The proof of Theorem 1 reposes on a dimensionality analysis of the graph of the monotone operator \(A_P: \R^d \toto \R^d\) as a subset of \(\R^d \times \R^d\). As we shall see, for an appropriate notion of dimension, the graph has dimension \(d\), and by using the Minty parameterization of the graph, we may divide the graph into pieces, at least one of which is \(d\)-dimensional. Then via a careful analysis using geometric measure theory, we may prove that the measure of any \(d\)-dimensional piece is appropriately preserved.

We provide the formal proof in Section 4, but include an overview of the ideas here, along with the main technique we adopt from convex analysis and monotone operator theory. To begin, we let \(A : \R^d \toto \R^d\) be a maximal monotone operator, and let \[\begin{align} \gr A \defeq \left\{(x, v) \mid x \in \dom A, v \in A(x)\right\} \subset \R^d \times \R^d \end{align}\] be its graph; we will follow [15] to develop the Minty parameterization of this subset, which provides a bi-Lipschitz homeomorphism between \(\R^d\) and \(\gr A\). To that end, for \(\lambda > 0\) define the proximal mapping, or resolvent, and Yosida approximation of \(A\) by \[\begin{align} \res{\lambda}{A} \defeq (I + \lambda A)^{-1} ~~~ and ~~~ \yosida{\lambda}{A} \defeq (\lambda + A^{-1})^{-1} = \frac{1}{\lambda} (I - \res{\lambda}{A}). \end{align}\] Then \(\res{\lambda}{A}\) is \(1\)-Lipschitz, while \(\yosida{\lambda}{A}\) is \(\sqrt{1 + \lambda^{-2}}\)-Lipschitz continuous [16] Remark 23.23. The Minty parameterization [15] \[\begin{align} \label{eqn:minty-param} \minty{\lambda}{A}(x) \defeq \left(\res{\lambda}{A}(x), \yosida{\lambda}{A}(x)\right) ~~ and ~~ \minty{\lambda}{A}^{-1}(y, v) \defeq y + \lambda v \end{align}\tag{7}\] parameterizes the graph \(\gr A\) in that \[\begin{align} \gr A = \left\{\minty{\lambda}{A}(x) \mid x \in \R^d \right\} \end{align}\] and \[(y, v) \in \gr A ~~ if and only if ~~ (y, v) = \minty{\lambda}{A}(x) ~ for ~ x = y + \lambda v.\] As \(\minty{\lambda}{A}\) is \(\sqrt{1 + \lambda^{-2}}\)-Lipschitz and \(\minty{\lambda}{A}^{-1}\) is \(\sqrt{1 + \lambda^2}\)-Lipschitz, this provides a bi-Lipschitz homeomorphism between \(\R^d\) and \(\gr A\).

We divide the proof of Theorem 1 into four parts; the first two hold for any maximal monotone operator \(A\), while the last two specialize to the sampling case that we consider:

1. Use the Minty parameterization to identify full-dimensional compact convex subsets \(K\) of the pre-image \(\res{\lambda}{A}^{-1}(S(\epsilon))\) of the \(\epsilon\)-stationary points, so \(\res{\lambda}{A}(K) \subset S(\epsilon)\).

2. Argue that the subset \(G = \gr A \cap (K \times \R^d)\) of the graph of \(A\) has a particular full-dimensional parameterization, in that there is a further (compact) \(X_0 \subset K\) for which \(X_0\) has dimension \((d - k)\) while \(A(x)\) contains a \(k\)-dimensional but “large enough” simplex for each \(x \in X_0\).

Step [item:full-dim-parameterization] is the most technically new part of the argument, and provides the key to convergence, as pieces of the graph decompose into product-like subsets. It also is the least quantitative, in that the measures of \(X_0\) and \(A(x)\) are only guaranteed to be positive and finite. We therefore provide discussion and examples in Section 4.7 after the proof of the theorem, highlighting the main challenge in providing quantitative guarantees.

3. Taking \(A_P = \int A_z(\cdot) dP(z)\), we argue that most of the sets \(A_P(x)\) from step [item:full-dim-parameterization] survive under sampling, that is, \(A_{P_n}(x)\) still contains a large enough set.

4. Argue that the \(d\)-dimensional Hausdorff measure, which we define formally in the sequel, of the empirical graph \(G_n = \gr A_{P_n} \cap (\res{\lambda}{A_P}(K) \times \R^d)\) is positive because of part [item:fiber-survival], and so the transformation \(\minty{\lambda}{A_{P_n}}^{-1}(G_n) \subset \R^d\) has positive \(d\)-dimensional Lebesgue measure because the Minty parameterization is a bi-Lipschitz homeomorphism.

Given steps [item:minty-to-full-subset]–[item:large-hausdorff], we can then argue that the empirical graph has good behavior: by step [item:large-hausdorff], the set of points \[\begin{align} \minty{\lambda}{A_{P_n}}^{-1}(G_n) = \left\{x \mid \minty{\lambda}{A_{P_n}}(x) \in G_n \right\} \subset \res{\lambda}{A_{P_n} + \normalcone_X}^{-1}(S(\epsilon)). \end{align}\] Because step [item:large-hausdorff] demonstrates that \(\mathsf{Leb}(\minty{\lambda}{A_{P_n}}^{-1}(G_n)) \ge c\) (for some constant \(c > 0\)) with high probability, this then implies the result.

3 Algorithms↩︎

Our algorithms, which rely on Theorem 1, proceed in two phases: the first uses a regularized empirical risk minimizer to find a point near the minimum norm solution \(x^0\) to the variational inequality 6 (recall Assumption 1 part [item:vi-moments]). Given such an initial point, we can then sample in a small neighborhood around \(\what{x}_n\), and Theorem 1 shows that one of these points likely has small subdifferential, and we can check this.

To set the stage, observe that strong monotonicity implies there is a unique solution \(x^\alpha\in X\) to the problem \[0 \in A_P(x) + \normalcone_X(x) + \alpha x.\] A quick argument (see Appendix 6.3) shows that \(x^\alpha\to x^0\) as \(\alpha\downarrow 0\):

Applying a standard convergence argument, this lemma in turn implies that the empirical regularized solution \(\what{x}_n \in X\) uniquely solving \[\begin{align} \label{eqn:l2-empirical-vi} 0 \in A_{P_n}(x) + \alpha x + \normalcone_X(x), ~~~ i.e. ~~~ \what{x}^\alpha= \proxmap_{\alpha^{-1} A_{P_n} + \normalcone_X}(0) \end{align}\tag{8}\] also converges to \(x^0\). We have the following proposition, which holds in arbitrary Hilbert spaces, whose proof we defer to Appendix 6.4.

Proposition [proposition:sampled-vi-convergence] shows that if we have access to a proximal point operator, then obtaining the “classical” guarantees of convergence, i.e., that \(\what{x}_n\) is close to nearly stationary point, is essentially immediate.

3.1 The randomized sampling algorithm↩︎

To actually guarantee stationarity, we arrive at a simple algorithm:

3.2 The convex minimization randomized sampling procedure↩︎

In the convex case where \(A_z(x)=\partial \loss_z(x)\) our Algorithm [alg:resample] becomes:

4 Proof of Theorem 1↩︎

This section provides the full proof of Theorem 1. We follow the outline in steps [item:minty-to-full-subset]–[item:large-hausdorff] in Sections 4.3 through 4.6. Before continuing, however, we provide some preliminary background on convergence of empirical subdifferentials (and monotone operators) to their population counterparts, as well as a brief review of the dimension theory we will use.

4.1 Convergence of empirical subdifferentials and monotone maps↩︎

We begin with guarantees on the (pointwise) convergence of subdifferentials of convex functions and, more generally, set-valued random functions. Recall the Hausdorff distance \[\dhaus(A, B) = \max\Big\{\sup_{a \in A} \dist(a, B), \sup_{b \in B} \dist(b, A)\Big\}\] between sets, where \(\dist(x, A) = \inf_{a \in A} \ltwo{x - a}\). Now, let \(S_z \subset \R^d\) be compact-convex-valued sets, measurable in \(z\), and let \(S_P = \E_P[S_Z] = \int S_z dP(z)\) be the (Aumann) expectation of \(S_Z\) for \(Z \sim P\) and \(S_{P_n} = \frac{1}{n} \sum_{i = 1}^n S_{Z_i}\) its empirical counterpart.

See Appendix 7.1 for a proof. The lemma also has the converse that there exist \(1\)-Lipschitz convex loss functions for which the empirical subdifferential has expected Hausdorff distance from the population subdifferential at least \(\min\{1, \sqrt{d/n}\}\), though we omit this result.

Lemma 4 extends to the monotone operator cases we consider: if \(A_z : \R^d \toto \R^d\) satisfy Assumption 1, then we immediately obtain the following result.

In some cases, we will truncate the random monotone operators we consider. Accordingly, we have the following result, whose proof we defer to Appendix 7.2.

For any \(t > 0\), in the setting of the lemma, \(\dhaus(C \cap 2B, D \cap 2B) \ge t\) therefore implies that either \(\dhaus(C, D) \ge \frac{c}{2}\) or \(\dhaus(C, D) \ge t/6\). Thus, if \(S_n\) is a random sequence of convex sets and \(B = c \cdot \ball\) satisfies \(S \cap B \neq \emptyset\), we obtain for any \(p < \infty\) and \(k \ge 0\) that \[\begin{align} \E\left[\dhaus^p(S_n \cap 2B, S \cap 2B) \wedge k \right] & = \int_0^{k^{1/p}} \P(\dhaus(S_n \cap 2B, S \cap 2B) \ge t^{1/p}) dt \\ & \le \int_0^{k^{1/p}} \left(\P(\dhaus(S_n, S) \ge t^{1/p} / 6) + \P(\dhaus(S_n, S) \ge c/2)\right) dt \\ & \le \E\left[6^p \dhaus^p(S_n, S) \wedge k\right] + k^{1/p} \cdot \P(\dhaus(S_n, S) \ge c/2). \end{align}\] The following lemma is then immediate.

4.2 A brief review of dimension theory↩︎

We provide definitions relating to dimension and Hausdorff measure that play a role in our arguments. For background, see the books on dimension theory by [38] and [39], and geometric measure theory by [40] and [41].

Defining the dimensional constant \(\omega_s \defeq \frac{\pi^{s/2}}{\Gamma(\frac{s}{2} + 1)}\) for \(s \ge 0\), the \(s\)-Hausdorff measure of a set \(X\) at diameter \(\delta\) is \[\begin{align} \hausmeasure^s_\delta(X) \defeq \inf\left\{\sum_{i = 1}^\infty \omega_s \left(\frac{\diam(U_i)}{2} \right)^s \mid \diam(U_i) < \delta, \bigcup_i U_i \supset X \right\} \end{align}\] the infimum taken over covers of \(X\) with diameters less than \(\delta\). Because the infima are over smaller collections, \(\hausmeasure^s_\delta\) is non-decreasing as \(\delta \downarrow 0\), so one defines \[\begin{align} \hausmeasure^s(X) = \sup_{\delta > 0} \hausmeasure_\delta^s(X) = \lim_{\delta \downarrow 0} \hausmeasure_\delta^s(X). \end{align}\] The \(0\)-dimensional Hausdorff measure \(\hausmeasure^0(X) = \card(X)\). For Borel \(X \subset \R^d\), we always have \(\hausmeasure^d(X) = \mathsf{Leb}(X)\), where \(\mathsf{Leb}\) denotes Lebesgue measure. The Hausdorff dimension of a set is \[\begin{align} \mathrm{dim}_{\hausmeasure}(X) \defeq \inf\left\{s \ge 0 \mid \hausmeasure^s(X) = 0 \right\}. \end{align}\] For convex sets \(X\), this coincides with the usual affine dimension of the set [42], and if \(X \subset \R^d\) has interior, then \(\mathrm{dim}_{\hausmeasure}(X) = d\) [41].

We also use the little inductive dimension \(\mathop{\mathrm{ind}}\), where recalling [38], [39], \(\mathop{\mathrm{ind}}X\) has the inductive definition that

1. \(\mathop{\mathrm{ind}}\emptyset = -1\)

2. \(\mathop{\mathrm{ind}}X \le d\) if for all points \(x \in X\) and all neighborhoods \(V \ni x\), there exists a neighborhood \(U\) of \(x\) with \(x \in \cl U \subset V\) and \(\mathop{\mathrm{ind}}\bd U \le d - 1\), where \(\bd U = \cl U \setminus U\).

3. \(\mathop{\mathrm{ind}}X = d\) if \(\mathop{\mathrm{ind}}X \le d\) and \(\mathop{\mathrm{ind}}X > d - 1\).

4. \(\mathop{\mathrm{ind}}X = +\infty\) if there is no \(d\) for which \(\mathop{\mathrm{ind}}X \le d\).

Standard results [38], [39] imply that \(\mathop{\mathrm{ind}}\R^d = d\), and we always have \(\mathop{\mathrm{ind}}X \le \mathrm{dim}_{\hausmeasure}X\) for all sets \(X\) [38]. The inductive dimension coincides with the topological dimension (the Lebesgue covering dimension) when \(X\) is a separable metric space [39]. As we work exclusively with separable metric spaces (subsets of \(\R^d\) or \(\R^{2d}\) inheriting the Euclidean metric structure), these notions of dimension are equivalent.

4.3 Step [item:minty-to-full-subset]: the pre-image of stationary points↩︎

We first consider a generic maximal monotone operator \(A\) with stationary points, where the Minty parameterization demonstrates that subsets of the graph \(\gr A\) have Hausdorff and inductive dimension \(d\). Let \(A : \R^d \toto \R^d\) have a solution \(x_0\) with \(0 \in A(x_0)\), and for \(\epsilon > 0\) let \[\begin{align} S(\epsilon) \defeq \left\{x \mid \dist(0, A(x)) \le \epsilon \right\}. \end{align}\] To state the result, recall the Minty parameterization 7 of the graph \(\gr A\).

4.4 Step [item:full-dim-parameterization]: full dimensional graph subsets↩︎

The second step of the proof demonstrates that if \(A\) is maximal monotone, we may decompose the graph \(\gr A\) into pieces that maintain full-dimensionality. Let \(0 < r< \infty\) be a (constant) radius of interest, \(X \subset \R^d\) be a compact subset of \(\dom A\), and define the graph piece \[\begin{align} G(X, r) \defeq \gr A \cap \left(X \times r\ball\right). \end{align}\] Because \(A\) is maximal monotone, \(\gr A\) is closed, so \(G(X, r)\) is a compact set. Define the projection \(\pi : G \to \R^d\), \(\pi(x, v) \defeq x\), onto the first \(d\) coordinates, and the slices (fibers) \[\begin{align} G_x \defeq A(x) \cap r\ball ~~ and ~~ \pi^{-1}(x) = \{x\} \times G_x, \end{align}\] where \(G_x\) is compact convex because \(A\) is maximal monotone and \(r< \infty\). The following lemma shows that we can decompose pieces of \(G \subset \R^d \times \R^d\) into parts of the form \((x, G_x)_{x \in X_k}\), where each slice \(G_x\) is \(k\)-dimensional, and for at least one \(k\), the base set \(X_k\) is \((d - k)\)-dimensional.

4.5 Step [item:fiber-survival]: the empirical slices survive under sampling↩︎

We specialize to the randomized sampling setting. Let \(x^0\) be the minimum norm solution to \(0 \in A_P(x) + \normalcone_X(x)\), so \(x^0\) also solves \(0 \in A_P(x) + \normalcone_{X \cap C}(x)\) for any closed convex set \(C\) containing \(x^0\). Without loss of generality, we may thus take \(X = X \cap B'\) for a ball \(B' \supset B = \{x^0 + b \ball\}\), where \(b\) is the value in the statement of Theorem 1 (i.e., Assumption 1.[item:vi-moments]); we slightly abuse notation to write \(X\) for the remainder of the proof. We also tacitly and without loss of generality assume \(\lambda\epsilon \le b\) (otherwise replace \(\epsilon\) with \(\epsilon' = b / \lambda\); the argument applies, as this can only shrink \(S_P(\epsilon)\)).

Lemma 8 shows that \(K = x^0 + \lambda \epsilon \ball\) satisfies \(\res{\lambda}{A}(K) \subset S_P(\epsilon) \cap K\), \(\res{\lambda}{A}(K) \subset K\) is compact, and that for \(\epsilon \le r/ 2 < \infty\), the piece \(\gr (A_P + \normalcone_X) \cap (\res{\lambda}{A}(K) \times (r/2) \ball)\) is compact with Hausdorff and inductive dimension \(d\), as is \[\begin{align} G_P \defeq \gr (A_P + \normalcone_X) \cap \left(\res{\lambda}{A}(K) \times r\ball\right). \end{align}\] Consider the empirical analogue \[\begin{align} \label{eqn:empirical-graph} A_n \defeq A_{P_n} + \normalcone_X ~~ and ~~ G_n \defeq \gr A_n \cap (\res{\lambda}{A}(K) \times r\ball) \end{align}\tag{11}\] of \(G_P\), where we let \(G_n(x) = A_n(x) \cap r\ball\) (we use \(\epsilon \le r/2\) to be able to apply Lemma 7). Lemma 9 guarantees the existence of \(0 \le k \le d\), \(\epsilon > 0\), and a compact set \(X_0\) such that \[\begin{align} X_0\subset \left\{\res{\lambda}{A}(K) \mid \mathrm{dim}_{\hausmeasure}G_P(x) = k ~ and ~ G_P(x) ~ contains~a~k-simplex~ \simplex_k ~ with ~ \epsilon \le \hausmeasure^k(\simplex_k) \right\} \end{align}\] with \(0 < \hausmeasure^{d - k}(X_0) < \infty\). We then define its empirical analogue \[\begin{align} \label{eqn:construct-empirical-goodset} X_{0,n}\defeq \left\{x \in X_0\mid G_n(x) ~ contains a ~ k-simplex~\simplex_k ~with~ \hausmeasure^k(\Delta_k) \ge \frac{\epsilon}{2}\right\}. \end{align}\tag{12}\] Then with high probability, the set \(X_{0,n}\) has large \(d-k\)-dimensional Hausdorff measure:

4.6 Step [item:large-hausdorff]: the measure of graph pieces↩︎

We complete the proof of Theorem 1, continuing in the setting of the previous section. We first demonstrate that the empirical graph piece 11 , \(G_n = \gr(A_{P_n} + \normalcone_X) \cap (\res{\lambda}{A}(K) \times r\ball)\), has reasonable \(d\)-dimensional Hausdorff measure for \(K = x^0 + \lambda \epsilon \ball\). Recall the projection mapping \(\pi(x, v) = x\), which is \(1\)-Lipschitz and whose inverse relative to \(G_n\), \(\pi|_{G_n}^{-1}(x) = \{x\} \times G_n(x)\), is bounded. Then as in the proof of Lemma 9, [40], implies \[\begin{align} \int_{X_{0,n}}^* \hausmeasure^k\left(G_n \cap \pi|_{G_n}^{-1}(x)\right) d \hausmeasure^{d - k}(x) \le C_{k,d-k} \cdot \hausmeasure^d(G_n). \end{align}\] By the construction 12 of \(X_{0,n}\), \(\hausmeasure^k(G_n \cap \pi_{G_n}^{-1}(x)) \ge \epsilon/2\) for each \(x \in X_{0,n}\), and so \[\begin{align} \frac{\epsilon}{2} \hausmeasure^{d-k}(X_{0,n}) \le C_{k,d-k} \hausmeasure^d(G_n). \end{align}\]

We may now combine the pieces: Lemma 11 shows that as \(n \to \infty\), for any \(\gamma > 0\) we have \(\hausmeasure^{d - k}(X_{0,n}) \ge (1 - \gamma) \hausmeasure^{d - k}(X_0)\) with high probability, and the application of Lemma 9 in step [item:fiber-survival] guarantees \(0 <\hausmeasure^{d - k}(X_0) < \infty\). Summarizing, we have shown that there exists a constant \(c > 0\) such that \[\begin{align} \label{eqn:hausdorff-of-empirical-graph} \P\left(\hausmeasure^d(G_n) \ge c \right) \to 1, \end{align}\tag{13}\] where \(G_n\) is the empirical graph piece 11 , which for \(K = \{x^0 + \lambda \epsilon \ball\}\) satisfies \[\begin{align} G_n \subset \gr A_n \cap \left(S_P(\epsilon) \cap K \times r\ball\right). \end{align}\]

As the final step, we use the Minty parameterization eq. ¿eq:eqn:minty-param?

we observe that \[\begin{align} \minty{\lambda}{A_n}^{-1}(G_n) = \{x + \lambda v \mid (x, v) \in G_n\} \subset \left\{x \mid \res{\lambda}{A_n}(x) \in S_P(\epsilon) \cap K \right\} \cap \{\res{\lambda}{A}(K) + \lambda r\ball\} \end{align}\] by the construction 11 as \((x, v) \in G_n\) implies \(\ltwo{v} \le r\). Lemma 8 part [item:resolvent-near-stationary] shows that \(\res{\lambda}{A}(K) \subset K = x^0 + \lambda \epsilon \ball\), and so we obtain \[\begin{align} \minty{\lambda}{A_n}^{-1}(G_n) \subset \res{\lambda}{A_n}^{-1}(S_P(\epsilon) \cap K) \cap \{x^0 + 2 \lambda r\ball\}. \end{align}\] Thus \[\begin{align} \hausmeasure^d\left(\minty{\lambda}{A_n}^{-1}(G_n)\right) = \mathsf{Leb}\left(\minty{\lambda}{A_n}^{-1}(G_n)\right) \le \mathsf{Leb}\left(\res{\lambda}{A_n}^{-1}(S_P(\epsilon) \cap K) \cap \{x^0 + 2 \lambda r\ball\} \right) \end{align}\] as \(\hausmeasure^d = \mathsf{Leb}\) on \(\R^d\). Because the Minty parameterization is a bi-Lipschitz homeomorphism, [41] then shows that \[\begin{align} \hausmeasure^d(G_n) \le \left(1 + \lambda^{-2}\right)^{d/2} \hausmeasure^d\left(\minty{\lambda}{A_n}^{-1}(G_n)\right). \end{align}\] The choice \(A_n = A_{P_n} + \normalcone_{X \cap B'}\) in the proof of step [item:fiber-survival] in Section 4.5 means that inequality 13 implies the theorem, because \(r\ge 2 \epsilon\) is otherwise arbitrary.

4.7 Discussion of step [item:full-dim-parameterization] and graph decomposition↩︎

The decomposition of \(\gr \partial f\) into pieces \(K_i = \{(x, g) \mid g \in \partial f(x), \dim \partial f(x) = i\}\), which sampling approximately preserves, constitutes the central tool in our results, and it guarantees \(\proxmap_{\poploss_{P_n}}^{-1}(S_P(\epsilon))\) has positive Lebesgue measure for \(\epsilon > 0\). The lack of quantitative guarantees makes it challenging to use as an ingredient to obtain rates of convergence. Part of this is that the “bases” \(\pi(K_i) = \{x \mid \dim \partial f(x) = i\}\) take many forms, as they can be countable, non-convex, and disconnected, making the application of convex-analytic tools non-obvious. We thus provide examples illustrating the decomposition in Lemma 9, after which we discuss the challenges of making the result more quantitative.

With these examples, we discuss the final point [item:base-graph-set] in Lemma 9, which provides a non-quantitative guarantee that \(\hausmeasure^{d-k}(X_k) > 0\) and that the images \(\partial f(x)\) for \(x \in X_k\) have positive \(k\)-dimensional Hausdorff measure \(\hausmeasure^k(\partial f(x))\). The challenging part of this result is to demonstrate that the base set \(X_k = \pi(K)\) for the projection \(\pi(x, v) = x\) has positive Hausdorff measure \(\hausmeasure^{d-k}(\pi(K)) > 0\). To do so, we use inductive dimension (leveraging the dimension of mappings), to (essentially) derive two-sided bounds on the Fubini-type integral 10 relating \(\hausmeasure^d(K_i)\) and \(\hausmeasure^{d - i}(\pi(K_i))\) for the partitioning sets \(K_i\) in part [item:partition-by-hausdorff] of the lemma.

The co-area formula [40] analogizes Fubini’s theorem for Hausdorff measures, and states that for sufficiently regular sets \(K \subset \R^N\) and \(u : K \to X \subset \R^m\) with \(m \le n\), \[\label{eqn:co-area} \int_{K} \mathsf{J}_m^K u(y) d\hausmeasure^n(y) = \int_{X} \hausmeasure^{n-m}\left(u^{-1}(x)\right) d\hausmeasure^m(x),\tag{14}\] where \(\mathsf{J}_m^K u\) denotes the \(m\)-dimensional Jacobian of \(u\) along the tangents to \(K\). In the integral 10 , we have \(u = \pi\), the (orthogonal) projection onto the first \(d\) components, \(N = 2d\), \(n = d\), and \(m = d - k\). The projection satisfies \(0 \le \mathsf{J}_{d-k}^K \pi \le 1\), and so (with regularity) \[\begin{align} \int_K \mathsf{J}_{d-k}^K \pi(y) d\hausmeasure^d(y) = \int_{\pi(K)} \hausmeasure^{k}(G \cap \pi^{-1}(x)) d\hausmeasure^{d-k}(x). \end{align}\] One might hope to guarantee a lower bound \(\mathsf{J}_{d-k}^K \pi \ge c > 0\) to develop a two-sided inequality \[\begin{align} c \cdot \hausmeasure^d(K) \le \int_{\pi(K)} \hausmeasure^k(G \cap \pi^{-1}(x)) d\hausmeasure^{d-k}(x) \le \hausmeasure^d(K), \end{align}\] which would allow quantitative bounds on \(\hausmeasure^{d-k}(\pi(K))\) in (e.g.) Lemma 11.

Even for \(1\)-Lipschitz \(\mc{C}^\infty\) convex functions, such a quantitative lower bound cannot generally hold. For example, for \(\epsilon > 0\), consider \(f_\epsilon(x) = \epsilon \log(e^{x/\epsilon} + e^{-x/\epsilon})\), which is \(1\)-Lipschitz with \(f_\epsilon'(x) = (e^{2x/\epsilon} - 1) / (e^{2x/\epsilon} + 1)\) and second derivative \(f_\epsilon''(0) = \epsilon^{-1}\). Because the \(\mc{C}^\infty\) map \(x \mapsto (x, f_\epsilon'(x))\) parametrizes the \(1\)-dimensional set \(K_0 = \{(x, f_\epsilon'(x))\}_{-1 \le x \le 1}\), this has \(1\)-dimensional Jacobian \(\mathsf{J}_1^{K_0} \pi(0) = 1 / \sqrt{1 + f''_\epsilon(0)^2} = (1 + o(1)) \cdot \epsilon\) as \(\epsilon \to 0\). This somewhat trivial example shows that, at least following our proof approach, new insights may be necessary to provide quantitative guarantees on the characteristics of the set \(\proxmap_{\poploss_{P_n}}^{-1}(S(\epsilon))\).

5 Discussion↩︎

Under local moment and maximal monotonicity assumptions, the randomized resampling estimator \(\what{x}_n\) that Algorithm [alg:resample] defines satisfies \[\dist\left(0, \popop(\what{x}_n) + \normalcone_X(\what{x}_n) \right) \to 0 ~~ in probability.\] For the subdifferential \(\partial \poploss_P\), this gives the convergence \(\dist(0, \partial \poploss_P(\what{x}_n) + \normalcone_X(\what{x}_n)) \cp 0\). Our guarantee does not require the operator to be single-valued or Lipschitz continuous, and does not rely on a smoothed surrogate as a proxy for stationarity. This distinction is central to our contribution: we control the stationary residual \(\dist\left(0, \popop(\what{x}_n) + \normalcone_X(\what{x}_n)\right)\) directly. Many statistical and optimization questions remain open.

A natural strengthening of our results would show that the random resampling step in Algorithm [alg:resample] is unnecessary, and that the empirical proximal mapping maps a neighborhood of \(x^0\) onto the \(\epsilon\)-stationary set. We record this as our first open question:

If the answer to this is affirmative, then the regularized empirical minimizer \[\what{x}^{\lambda_n} = \argmin_{x \in X} \left\{\poploss_{P_n}(x) + \frac{\lambda_n}{2} \ltwo{x}^2 \right\},\] which converges to \(x^0\) (Proposition [proposition:sampled-vi-convergence]), would eventually lie in this neighborhood with high probability. A single empirical proximal point step, initialized at \(\what{x}^{\lambda_n}\), would then produce a population \(\epsilon\)-stationary point. This would replace the resampling-and-validation scheme by a simpler algorithm while preserving the same stationarity guarantees.

Two special cases suggest the answer to the question should be true, and they allow a much simpler argument than the dimensionality results we develop: when the population objective \(\poploss_P\) is differentiable at the minimum norm minimizer \(x^0\) and when \(\partial \poploss_P(x^0)\) contains a neighborhood of \(0\) (i.e., \(\poploss_P\) has a sharp minimizer). In each of these, regularized minimizers exhibit asymptotic stationarity. For simplicity, we consider cases when Assumption 1.[item:vi-moments] holds with \(p\) moments near \(x^0\), meaning that \(\loss_z\) is \(\lipconst(z)\)-Lipschitz at \(x^0\) with \(\E[\lipconst^p(Z)] < \infty\).

One additional question regards uniformity of convergence: our results are pointwise in the distribution \(P\). This contrasts with classical optimization. To make the contrast apparent, suppose that the functions \(\{\loss_z\}_{z \in \mc{Z}}\) are \(\lipconst\)-Lipschitz and the domain \(X\) is bounded. Then stochastic gradient methods achieve (worst-case optimal) convergence rate \(\E[\poploss_P(x_k)] - \poploss_P\opt \lesssim \lipconst \cdot \diam(X) / \sqrt{k}\), independent of \(P\) and dimension [43], [44], and letting \(\mc{P}\) be the collection of distributions on \(\mc{Z}\), the distribution-free guarantee \(\inf_{\what{x}_n} \sup_{P \in \mc{P}} \E_{P^n}[\poploss_P(\what{x}_n) - \poploss_P\opt] \to 0\) holds if and only if the losses \(\{\loss_z\}_{z \in \mc{Z}}\) are uniformly Lipschitz on each compact \(X_0 \subset \relint X\) [45]. In a companion paper we are preparing, we show that in distinction, procedures achieving stationarity 2 must exhibit some dependence on the underlying dimension.

We nonetheless have a natural open question:

That is, we ask whether it is possible to derive a convergence rate to stationarity for \(\lipconst\)-Lipschitz convex objectives. One avenue would be to adapt [46]’s stability-based arguments [46], which show that stochastic optimization is possible if and only if estimators exist that are appropriately stable to replacing an individual example \(Z_i\). Thus far, the natural instability of subdifferential mappings has made it difficult for us to generalize replace-one stability to Hausdorff or other set distances, leaving Question [question:uniform] quite open.

6 Miscellaneous deferred proofs↩︎

6.1 Proof of Observation [observation:conditional-coverage]↩︎

Writing the subdifferential of \(\loss_{x,z}(\theta) = \alpha \hinge{\theta^\top \phi(x) - z} + (1 - \alpha) \hinge{z - \theta^\top \phi(x)}\) and using the shorthand \(\tau_\theta(x) = \theta^\top \phi(x)\), we obtain \[\begin{align}\\ & = \E_P\big[\alpha \phi(X) \indic{Z < \tau_\theta(X)} - (1 - \alpha) \phi(X) \indic{Z > \tau_\theta(X)} + [-(1 - \alpha), \alpha] \cdot \phi(X) \indic{Z = \tau_\theta(X)} \big] \\ & = -\E\left[\phi(X) \indic{Z > \tau_\theta(X)}\right] + \alpha \E[\phi(X)] + \E\left[[-1, 0] \cdot \phi(X) \indic{Z = \theta^\top \phi(X)}\right]. \end{align}\] So if there is a vector \(u \in \partial \poploss_P(\theta)\), there is a functional \(0 \le s(x, z) \le 1\) for which \[\begin{align} u = \alpha \E[\phi(X)] - \E\left[\phi(X) \indic{Z > \tau_\theta(X)}\right] - \E\left[s(X, Z) \cdot \phi(X) \indic{Z = \tau_\theta(X)}\right]. \end{align}\] Rewriting this by conditioning on \(X\), we have \[\begin{align} u = \E\left[\phi(X) \cdot (\alpha - P(Z > \tau_\theta(X) \mid X))\right] - \E\left[s(X, Z) \cdot \phi(X) \indic{Z = \tau_\theta(X)}\right]. \end{align}\] In particular, for all vectors \(v \in \R^d\) and nonnegative weighting functions of the form \(w(x) = v^\top \phi(x)\), we have \[\begin{align} u^\top v \le \E\left[w(X) \cdot \left(\alpha - P(Z > \tau_\theta(X) \mid X)\right) \right] \end{align}\] or, renormalizing, \(P_w(Z > \tau_\theta(X)) \le \alpha + \epsilon \ltwo{v} / \E[w(X)]\).

6.2 Constructing a separated subset of the sphere↩︎

For \(\epsilon > 0\) and \(x\in \sphere^{d-1}\), let \(C_\epsilon = \{y \in \sphere^{d-1} \mid \ltwo{y - x} \le \epsilon\}\) be the spherical cap centered at \(x\). Then for \(\epsilon\) small the surface area \(\mathrm{SA}_d(C_\epsilon) \le c_d \epsilon^{d-1}\) while \(\mathrm{SA}_d(\sphere^{d-1}) = \frac{2 \pi^{d/2}}{\Gamma(d/2)} = C_d\), where \(c_d\) and \(C_d\) are dimension dependent constants. Let \(\varepsilon_i > 0\) be a non-increasing sequence satisfying \(\frac{c_d}{C_d} \sum_{i = 1}^\infty \varepsilon_i^{d-1} \le 1\), and recursively define a set of points as follows. Take \(U_1 = \sphere^{d-1}\), \(u_1 \in U_1\), and construct \[\begin{align} U_{i + 1} = U_i \setminus \{u_i + \varepsilon_i \ball_2^d\}, \end{align}\] taking \(u_i \in U_i\). At each step \(i\), any point \(u_j\) for \(j < i\) satisfies \(\ltwo{u_i - u_j} \ge \varepsilon_j\), while each point \(j > i\) satisfies \(\ltwo{u_i - u_j} \ge \varepsilon_i\); we can never have \(U_i = \emptyset\) by a volume argument because \(\frac{c_d}{C_d} \sum_i \varepsilon_i^{d-1} \le 1\). So \(\ltwo{u_i - u_j} \ge \varepsilon_i\) for each \(i \neq j\), and thus we have inequality 4 : \[\begin{align} \sup_{j \neq i} \<u_i, u_j\> = \sup_{j \neq i} \left\{1 - \half \ltwo{u_i - u_j}^2 \right\} \le 1 - \half \varepsilon_i^2 < 1. \end{align}\]

6.3 Proof of Lemma 1↩︎

First, we claim that \(\ltwos{x^0} \ge \ltwos{x^\lambda}\) for all \(\lambda > 0\). By definition of \(x^0\) and \(x^\lambda\) as solutions to their respective V.I.s, we have \[\begin{align} 0 \in A(x^0) + \normalcone_X(x^0) ~~ and ~~ -\lambda x^\lambda \in A(x^\lambda) + \normalcone_X(x^\lambda). \end{align}\] Then by monotonicity, for the elements \(a^0 \in A(x^0)\) and \(w^0 \in \normalcone_X(x^0)\) and \(a^\lambda \in A(x^\lambda)\) and \(w^\lambda \in \normalcone_X(x^\lambda)\) realizing \(0 = a^0 + w^0\) and \(a^\lambda + w^\lambda + \lambda x^\lambda = 0\), \[\begin{align} \<0 +\lambda x^\lambda, x^0 - x^\lambda\> = \<a^0 + w^0 - a^\lambda - w^\lambda, x^0 - x^\lambda\> \ge 0. \end{align}\] Rearranging and applying Cauchy-Schwarz, \[\begin{align} \lambda \ltwos{x^\lambda}^2 & \le \lambda \<x^0, x^\lambda\> \le \lambda \ltwos{x^0} \ltwos{x^\lambda}, \end{align}\] and dividing by \(\ltwos{x^\lambda}\) gives the claim (except if \(\ltwos{x^\lambda} = 0\), in which case it is trivial).

We now argue that \(x^\lambda \to x^0\). Consider any potential limit of a (subsequence) of \(x^\lambda\) as \(\lambda \downarrow 0\), i.e., \(x^\lambda \to x\). We claim \(x\) solves 6 . The outer semicontinuity of \(A+\normalcone_X\) yields \(\limsup_{\lambda \to 0} A(x^\lambda)+\normalcone_X(x^\lambda) \subset A(x)+\normalcone_X(x)\); let \(u^\lambda \in A(x^\lambda)+\normalcone_X(x^\lambda)\) be the vector satisfying \(u^\lambda+\lambda x^\lambda=0\). Then as \(\lambda x^\lambda \to 0\), we obtain that \(\lim_{\lambda \to 0} u^\lambda = 0 \in A(x)+\normalcone_X(x)\). Hence, since \(\ltwos{x^0} \ge \ltwos{x^\lambda}\), \(x\) must be the minimum norm element of \(X\opt\).

6.4 Proof of Proposition [proposition:sampled-vi-convergence]↩︎

We use a typical stochastic convergence argument.

The final step uses Assumption 1. We see that for small enough \(\alpha\), \(x^\alpha\) belongs to the neighborhood of \(x^0\) Assumption 1 specifies, and so \[\begin{align} \ltwobig{\what{x}^{\alpha} - x^0} \le \ltwobig{x^\alpha- x^0} + \frac{1}{\alpha} \ltwobig{a_{P_n}(x^\alpha) - a_P(x^\alpha)}. \end{align}\] Taking an expectation yields that \(\E[\ltwos{a_{P_n}(x^\alpha) - a_P(x^\alpha)}^2] \lesssim \frac{1}{n}\), and so long as \(n \alpha_n^2 \to \infty\) while \(\alpha_n \to 0\) we have the first result. For the almost-sure convergence, observe by the vector Khinchine-Kahane inequalities and a symmetrization argument [47] that for any mean zero vectors \(W_i \in \R^d\), we have for \(p > 2\) that \[\begin{align} \E\left[\ltwobigg{\sum_{i = 1}^n W_i}^p\right] \le 2^{p} \E\left[\ltwobigg{\sum_{i = 1}^n \varepsilon_i W_i}^p\right] & \lesssim \E\left[\E\left[\ltwobigg{\sum_{i = 1}^n \varepsilon_i W_i}^2 \mid W_1^n \right]^{p/2}\right] \\ & = \E\left[\bigg(\sum_{i = 1}^n \ltwo{W_i}^2\bigg)^{p/2}\right] \le n^{p/2 - 1} \sum_{i = 1}^n \E[\ltwo{W_i}^p], \end{align}\] where \(\varepsilon \simiid \uniform\{-1, 1\}\) are an independent Rademacher sequence. Accordingly, given \(p\)th moment, we obtain \[\begin{align} \E\left[\ltwobig{a_{P_n}(x^\alpha) - a_P(x^\alpha)}^p\right] \lesssim n^{-p/2} \E\left[\ltwobig{a_Z(x^\alpha)}^p\right], \end{align}\] so that \[\begin{align} \sum_{n = 1}^\infty \P\left(\frac{1}{\alpha_n} \ltwobig{a_{P_n}(x^{\alpha_n}) - a_P(x^{\alpha_n})} > \epsilon\right) \lesssim \frac{1}{\epsilon^p} \sum_{n = 1}^\infty \frac{1}{\alpha_n^p \cdot n^{p/2}} \cdot \E\left[\ltwobig{a_Z(x^\alpha)}^p\right]. \end{align}\] So long as \(\sum_n \frac{1}{\alpha_n^p n^{p/2}} < \infty\), the Borel-Cantelli lemma implies almost sure convergence.

7 Deferred proofs for Theorem 1↩︎

7.1 Proof of Lemma 4↩︎

Using the standard identification [1] of the Hausdorff distance with support functions of sets, that is, that \(\dhaus(A, B) = \sup_{\ltwo{v} \le 1} |\sigma_A(v) - \sigma_B(v)|\) for compact convex \(A, B\) with \(\sigma_A(v) = \sup_{a \in A} \<a, v\>\), we see that \[\begin{align} \dhaus\left(S_{P_n}, S_P\right) & = \sup_{\ltwo{v} \le 1} \left|\frac{1}{n} \sum_{i = 1}^n \sigma_{S_{Z_i}}(v) - \sigma_{S_P}(v) \right|. \end{align}\] Introducing the random signs \(\varepsilon_i \simiid \uniform\{-1, 1\}\), a symmetrization argument [47] gives \[\begin{align} \E\left[\dhaus^p\left(S_{P_n}, S_P\right)\right] & \le C_p \cdot \E\left[\sup_{\ltwo{v} \le 1} \left|\frac{1}{n} \sum_{i = 1}^n \varepsilon_i \sigma_{S_{Z_i}}(v) \right|^p\right], \end{align}\] where \(C_p\) depends only on \(p\).

Define the function \(H_n(\varepsilon) \defeq \sup_{v \in \ball} |\sum_{i = 1}^n \varepsilon_i \sigma_{S_{Z_i}}(v)|\), which is a \(\sqrt{\sum_{i = 1}^n \lipconst^2(Z_i)}\)-Lipschitz function of \(\varepsilon \in \{-1, 1\}^n\). Then the convex concentration inequality [5] states that \[\begin{align} \P\left( \left|H_n(\varepsilon) - \E[H_n(\varepsilon) \mid P_n]\right| \ge t \mid P_n \right) \le 2 \exp\left(-\frac{t^2}{2 n \E_{P_n}[\lipconst^2(Z)]} \right) \end{align}\] for all \(t \ge 0\), where \(\E_{P_n}[\lipconst^2(Z)] = \frac{1}{n} \sum_{i = 1}^n \lipconst^2(Z_i)\). Using that \(\E[X] = \int_0^\infty P(X \ge t) dt\) for \(X \ge 0\), we thus obtain \[\begin{align} \E\left[n^p \dhaus^p\left(S_{P_n}, S_P\right)\right] & \le C_p \cdot \E\left[|H_n(\varepsilon)|^p\right] \\ & \le 2^p C_p \cdot \left(\E\left[ \E\left[ \left|H_n(\varepsilon) - \E[H_n(\varepsilon) \mid P_n]\right|^p \mid P_n\right]\right] + \E\left[\left|\E[H_n(\varepsilon) \mid P_n]\right|^p \right]\right). \end{align}\] For the first inner conditional expectation, we observe that \[\begin{align} \E\left[\left|H_n(\varepsilon) - \E[H_n(\varepsilon) \mid P_n]\right|^p \mid P_n\right] & = \int_0^\infty \P\left(|H_n(\varepsilon) - \E[H_n(\varepsilon) \mid P_n] | \ge t^{1/p}\right) dt \\ & \le 2 \int_0^\infty \exp\left(-\frac{t^{2/p}}{2n \E_{P_n}[\lipconst^2(Z)]} \right) dt \\ & = p \Gamma(p) n^{p/2} \E_{P_n}[\lipconst^2(Z)]^{p/2} \end{align}\] by a standard gamma integral computation. To upper bound the expectation \(\E[|H_n(\varepsilon)| \mid P_n]\) we use a standard entropy integral [5], [48] to obtain \(\E[|H_n(\varepsilon)| \mid P_n] \lesssim \sqrt{d n \E_{P_n}[\lipconst^2(Z)]}\). Substituting above yields \[\begin{align} \E\left[n^p \dhaus^p(S_{P_n}, S_P)\right] & \le C_p \cdot (dn)^{p/2} \cdot \E\left[\E_{P_n}[\lipconst^2(Z)]^{p/2}\right] \le C_p \cdot (dn)^{p/2} \cdot \E[\lipconst^p(Z)] \end{align}\] by Jensen’s inequality.

7.2 Proof of Lemma 6↩︎

Let \(\delta = \dhaus(C, D) \le \frac{c}{2}\). Then each \(x \in D\) satisfies \(\dist(x, C) \le \delta\), and similarly, each \(y \in C\) satisfies \(\dist(y, D) \le \delta\). We proceed in three steps: first, arguing that \(C \cap 2B\) is non-empty, then arguing that \(\dist(x, C \cap 2B) \lesssim \delta\) for all \(x \in D \cap 2B\), and finally that \(\dist(y, D \cap 2B) \lesssim \delta\) for all \(y \in C \cap 2B\).

1. We first argue that \(C \cap 2B\) is non-empty. Take \(x_0 \in D \cap B\) and the \(y_0 \in C\) satisfying \(\ltwo{x_0 - y_0} \le \delta\), so \(\ltwo{y_0} \le \ltwo{x_0} + \delta < 2c\), and \(y_0 \in 2 B\). Thus \(C \cap 2B\) is non-empty.

2. Take \(x \in D \cap 2B\), and let \(y \in C\) satisfy \(\ltwo{x - y} \le \delta\), so \(\ltwo{y} \le 2c + \delta\). Let \(0 < t < 1\) to be chosen and let \[\begin{align} z = (1 - t) y + t y_0 \in C, \end{align}\] where the inclusion follows by convexity of \(C\). Computing norms, we have \[\begin{align} \ltwo{z} & \le (1 - t) \ltwo{y} + t \ltwo{y_0} \le (1 - t) (\ltwo{x} + \delta) + t (\ltwo{x_0} + \delta) \\ & \le (1 - t) (2c + \delta) + t(c + \delta) = 2 c - t c + \delta, \end{align}\] which satisfies \(\ltwo{z} \le 2c\) whenever \(\delta \le t c\), that is, \(z \in 2B\). Taking \(t = \frac{\delta}{c}\) then yields \[\begin{align} \ltwo{z - x} \le \ltwo{z - y} + \ltwo{x - y} \le t \ltwo{y - y_0} + \delta \le \frac{\delta}{c}(c + \delta + 2c + \delta) + \delta \le 6 \delta. \end{align}\]

3. Take \(y \in C \cap 2B\). Let \(x \in D\) be the closest point to \(y\), so \(\ltwo{x - y} \le \delta\) and \(\ltwo{x} \le 2c + \delta\). Take any point \(x' \in D \cap B\), and for \(0 \le t \le 1\) to be chosen let \(z = (1 - t) x + t x' \in D\). Evidently, if \((1 - t) (2c + \delta) + t c = 2c - tc + (1 - t) \delta \le 2 c\), i.e., \(\delta \le t (c + \delta)\), then \(z \in D \cap 2B\); we may take \(t = \delta / c\) to achieve this because \(c \ge 2\delta\). In this case \[\begin{align} \ltwo{y - z} \le \ltwo{y - x} + t \ltwo{x - x'} \le \delta + \frac{\delta}{c} (2c + \delta + c) \le 5 \delta. \end{align}\]

As \(\delta = \dhaus(C, D) \le c/2\), we see that \[\begin{align} \dhaus(C, D) \le \frac{c}{2} ~~ implies ~~ \dhaus(C \cap 2B, D \cap 2B) \le 6 \cdot \dhaus(C, D). \end{align}\]

References↩︎