Stochastic L-convex Function Minimization
Year: 2021

Abstract:
We study an extension of the stochastic submodular minimization problem, namely, the stochastic L-convex minimization problem. We develop the first polynomialtime algorithms that return a near-optimal solution with high probability. We design a novel truncation operation to further reduce the computational complexity of the proposed algorithms. When applied to a stochastic submodular function, the computational complexity of the proposed algorithms is lower than that of the existing stochastic submodular minimization algorithms. In addition, we provide a strongly polynomial approximate algorithm. The algorithm execution also does not require any prior knowledge about the objective function except the L-convexity. A lower bound on the computational complexity that is required to achieve a high probability error bound is also derived. Numerical experiments are implemented to demonstrate the efficiency of our theoretical findings.

heading: 1 Introduction
text: As an important branch of optimization, discrete stochastic optimization problems (Futschik & Pflug, 1995; Gutjahr & Pflug, 1996; Futschik & Pflug, 1997; Kleywegt et al., 2002; Semelhago et al., 2020) appear in a wide range of applications such as image processing (Jegelka & Bilmes, 2011), pricing optimization (Ito & Fujimaki, 2016) and simulation-based optimization (Altman et al., 2003; Freund et al., 2017). Since a discrete optimization problem is generallyNP-hard to solve, specific structures are necessary to guarantee the existence of polynomial-time algorithms. In the special case when the feasible set is {0, 1}d for some positive integer d, submodular functions are viewed as the discrete analogy of convex functions in the field of combinatorial optimization. The submodular function minimization (SFM) problem has wide applications in computer vision, economics, game theory and is well-studied in literature (Lee et al., 2015; Axelrod et al., 2020; Zhang et al., 2020). However, for certain problems, because of the stochastic structure, the exact value of the objective function is not available but only noisy samples of the objective function value are available. These scenarios, for example, include the optimization of expected performances for complex stochastic systems. In contrast to the deterministic problem, the stochastic SFM problem is less understood (Ito, 2019; Chakrabarty et al., 2017; Blais et al., 2019; Axelrod et al., 2020).
In this work, we consider the stochastic minimization problem of a class of generalized submodular functions, named the L\-convex functions (Murota, 2003). More concretely, a function
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
f : {1, . . . , N}d 7→ R is called an L\-convex function if it satisfies the discrete midpoint convexity:
f(x) + f(y) ≥ f(d(x+ y)/2e) + f(b(x+ y)/2c), ∀x, y ∈ {1, . . . , N}d, where N, d are positive integers and d·e, b·c are the ceiling and flooring functions applied componentwisely to vectors. We note that L\-convex functions are reduced to submodular functions if N = 2. As an extension to the stochastic SFM problem, the stochastic L\-convex minimization problem also appears in a variety of decision making problems. For example, the user satisfaction function of a bike-sharing system (Freund et al., 2017) and the expected total waiting time of a service network (Altman et al., 2003) are proved to be multimodular functions, which is equivalent to an L\-convex function under a linear transformation. Since the evaluation of both functions includes the simulation of random user behaviors, we only have access to noisy objective function value oracles and thus both problems can be formulated as stochastic L\-convex minimization problems. We refer the readers to the appendix and Murota (2003) for more concrete examples of L\-convex minimization problems.
The contributions of this work are two-fold. For the algorithm design, we provide the first polynomialtime algorithms to achieve high-probability bounds for the stochastic L\-convex function minimization problem. We note that we only consider the number of objective function value evaluations in the computational complexity, since the function value evaluation oracle is usually the major contributor to the running time of algorithms. The first algorithm is based on the observation that the Lovász extension of an L\-convex function transforms the original problem to an equivalent continuous convex optimization problem. Therefore, we can apply the stochastic subgradient method to solve the continuous convex problem and then round the solution to an integral solution. The idea of combining the Lovász extension and the stochastic subgradient method has already been considered in Hazan & Kale (2011); Ito (2019) for stochastic online convex optimization and stochastic submodular minimization, respectively. However, the direct application of the stochastic subgradient method will lead to a large computational complexity. We add an extra truncation step to the stochastic subgradient and the computational complexity is reduced by O(d). We also implement numerical experiments to show the significance of the truncation step. In the special case when N = 2, our problem reduces to the stochastic SFM problem and the computational complexity of the truncated subgradient method is O(d2) smaller than the stochastic cutting-plane method in Blais et al. (2019). In addition, we propose the dimension reduction method and get a strongly polynomial approximate algorithm (see the definition in Section 4.2) with computational complexity O(d3N2(d + log(N)) −2 log(1/δ)). Finally, the uniform sampling method in Zhang et al. (2021) is extended to the multi-dimensional case and the computational complexity is independent of N . Table 1 summarizes the upper bounds on the computational complexities to find a solution whose objective value is at most larger than the optimum with probability at least 1− δ. On the other hand, we derive lower bounds on the computational complexity for achieving a highprobability error bound. We construct a set of L\-convex functions whose landscapes are close to each other while their optima are different with each other. Intuitively, the algorithm requires sufficiently many evaluations to differentiate those “similar”L\-convex functions, where the similarity between functions is characterized by the KL-divergence. More rigorously, we utilize the information inequality in Kaufmann et al. (2016) and get the lower bound O(d −2 log(1/δ)) on the complexity. Recently, Ito (2019) proved a lower bound for achieving an -expected error. There is a difference in terms of d and δ regarding the lower bound results between ours and Ito’s. We note that these two lower bounds are established under different settings, i.e., the high-probability setting and the expectation setting. Hence, a directly comparison of the two lower bounds is unavailable; see the discussion in Section 5.

heading: 1.1 Notation
text: For N ∈ N, we define [N ] := {1, 2, . . . , N}. For a given set S and an integer d ∈ N, the product set Sd is defined as {(x1, x2 . . . , xd) : xi ∈ S, i ∈ [d]}. For two vectors x, y ∈ Rd, we use (x ∧ y)i := min{xi, yi} and (x ∨ y)i := max{xi, yi} to denote the component-wise minimum and maximum. Similarly, the ceiling function d·e and the flooring function b·c round each component to an integer when applied to vectors. The failing probability of algorithms is denoted as δ. The notation f = O(g) (resp. f = Θ(g)) means that there exist absolute constants c1, c2 > 0 such that f ≤ c1g (resp. c1g ≤ f ≤ c2g). Similarly, the notation f = Õ(g) (resp. f = Θ̃(g)) means that there exist absolute constants c1, c2 > 0 and constant c3 > 0 independent of δ such that f ≤ c1g + c3 (resp. c1g ≤ f ≤ c2g + c3).

heading: 2 Problem Formulation
text: We have in mind a complex stochastic system and consider the following discrete stochastic optimization problem
min x
f(x) := E[F (x, ξx)] s.t. x ∈ X , (1)
in which X is the discrete feasible set, ξx summarizes all the randomness supported on some proper probability space (Y,BY) and F : X×Y 7→ R is a deterministic function. Specifically, ξx summarizes all the random variables and stochastic processes involved in the system for the decision variable x. The function F maps the decision and associated randomness into a random system performance. For example, in service systems such as the bike-sharing service, this function F captures the full operational logic regarding how services are arranged. We write ξx,1, ξx,2, . . . , ξx,n as independent and identically distributed (i.i.d.) copies of ξx. We use F̂n(x) := 1n ∑n j=1 F (x, ξx,j) to denote the empirical mean of the n independent evaluations for the choice of decision variables x. The triad of the set X , the space of randomness (Y,BY) and the function F (·, ·) is called the model of problem (1). In this work, we assume the L\-convexity of the objective function. For the ease of presentation, we assume the feasible set X = [N ]d, where d denotes the dimension and N represents the number of choices on each dimension for the decision variable x.
Assumption 1. The objective function f(·) is an L\-convex function on X = [N ]d, i.e.,
f(x) + f(y) ≥ f(d(x+ y)/2e) + f(b(x+ y)/2c), ∀x, y ∈ X ,
Remark 1. In general, the definition of L\-convexity also requires the feasible set X to be an L\convex set, namely, the indicator function of X should be L\-convex on Zd. In our case, the feasible set X = [N ]d is an L\-convex set. Hence, we only need the discrete mid-point convexity property to define L\-convex functions on X .
Moreover, we assume that the distribution of the stochastic output F (x, ξx) is sub-Gaussian.
Assumption 2. The distribution of F (x, ·) is sub-Gaussian with a known parameter σ2 for all x ∈ X .
The knowledge of σ2 represents the maximum level of stochasticity involved in the system performance. For scenarios where this σ2 is not known, it has to be adaptively learned and adds additional complexity to the problem, but is not the focus of this work. We further note that our algorithm design and analysis can be generalized to any probability distribution with provable concentration inequalities. The set of all models satisfying both Assumption 1 and Assumption 2 is denoted as MC(X ), or simplyMC. The main goal is to design algorithms that guarantee the selection of a near-optimal solution with high probability. Formally, this criterion is called the ( , δ)-Probably Approximately Correct (PAC) guarantee; see Even-Dar et al. (2002); Kaufmann et al. (2016); Ma & Henderson (2017). We say that a solution x̂, as the output of an algorithm, satisfies the ( , δ)-PAC guarantee if
P [ f(x̂)−min
x∈X f(x) ≥
] ≤ δ,
where the probability is calculated over the distribution of x̂.
An algorithm is called an [( , δ)-PAC,MC]-algorithm if it always returns an ( , δ)-PAC solution for any model inMC. The computational complexity of an algorithm is defined as the number of evaluations of F (·, ·) before termination. The difficulty of satisfying the PAC guarantee for any model inMC is characterized by the worst-case performance of an optimal algorithm. Definition 1. The computational complexity of the ( , δ)-PAC guarantee and classMC is defined as
T ( , δ,MC) := inf A is (( ,δ)-PAC,MC) sup M∈MC E [τA] ,
where the random variable τA is the computational complexity for each implementation of A.
The reminder of this paper is focused on providing upper bounds on T ( , δ,MC) by designing algorithms and providing lower bounds by an information-theoretical analysis. In this paper, we consider the behavior of algorithms in the high-probability regime, namely, in the setting when the failing probability δ is small enough.
Assumption 3. The failing probability δ is sufficiently small.
In general, the sub-Gaussian randomness leads to
T ( , δ,M) = C( ,M) · log (1/δ) 0 < δ 1,
where C( ,M) is some function only depending on the precision and the model classM. For example, it is known that the sample complexity is Θ(Nd −2) log(1/δ) for the general best-arm selection problem (Even-Dar et al., 2002) when no convex structure is in presence.

heading: 3 Truncated Stochastic Subgradient Algorithms
text: In this section, we develop gradient-based algorithms for problem (1). We first show that the properties of L\-convex functions make it possible to construct a continuous convex function on the continuous cube [1, N ]d that has the same set of optima as the original objective function. Then, a rounding process is designed to round the solution of the continuous problem to an integral solution. Finally, the truncated stochastic subgradient algorithm is proposed and its computational complexity is analyzed. Moreover, in this section, we make the assumption that the `∞-Lipschitz constant of the objective function can be estimated a priori. For example, this estimation is possible for service systems optimization, where it represents an upper bound on how much impact an additional server over a certain time range contributes to the system performance.
Assumption 4. The `∞-Lipschitz constant of f(x) is known to be L a priori. Namely, we know beforehand that
|f(x)− f(y)| ≤ L, ∀x, y ∈ X s.t. ‖x− y‖∞ ≤ 1.

heading: 3.1 Reduction to a continuous convex optimization problem
text: We first show that the original discrete minimization problem (1) can be equivalently transformed to a continuous convex optimization problem. Based on the fact that the restriction of an L\-convex function to a hypercube is a submodular function (Murota, 2003), we show that we can piece together the Lovász extension of the objective function in each hypercube to form a continuous convex function on conv(X ) = [1, N ]d. More importantly, it can be shown that the continuous convex function has exactly the same global minima as the original L\-convex function and an unbiased subgradient estimator is available at each point at the cost of 2d evaluations to F (·, ·). Theorem 1. Suppose that function f(·) is an L\-convex function on [N ]d. Then, there exists a continuous convex function f̃(·) on [1, N ]d, which we also name as the Lovász extension, such that
• It holds f̃(x) = f(x) for any x ∈ X .
• The minimizers of f̃ satisfy arg miny∈[1,N ]d f̃(y) = arg miny∈[N ]d f(y).
• For any point x ∈ [1, N ]d, there exists an unbiased estimator ĝx to a subgradient of f̃(·) at point x, which can be evaluated via 2d evaluations to F (·, ·).
Algorithm 1 Rounding process to a feasible solution Input: Model X ,BY, F (x, ξx), optimality guarantee parameters , δ, ( /2, δ/2)-PAC solution x̄ to
problem (2). Output: An ( , δ)-PAC solution x∗ to problem (1).
1: Compute the neighbouring points (defined in the appendix) of x̄, denoted as S0, . . . , Sd. 2: Simulate at Si until the 1− δ/(4d) confidence half-width of F̂n(Si) is smaller than /4 for all i. 3: Return the optimum x∗ ← arg minS∈{S0,...,Sd} F̂n(S).
Algorithm 2 Projected and truncated SSGD method for the PAC guarantee Input: Model X ,BY, F (x, ξx), optimality guarantee parameters , δ, number of iterations T , step
size η, truncation threshold M . Output: An ( , δ)-PAC solution x∗ to problem (1).
1: Choose the initial point x0 ← ((N + 1)/2, . . . , (N + 1)/2). 2: for t = 0, . . . , T − 1 do 3: Generate a stochastic subgradient ĝt at xt. 4: Truncate the stochastic subgradient g̃t ← TM (ĝt). 5: Update xt+1 ← P (xt − ηg̃t). 6: end for 7: Compute the averaging point x̄← (∑T−1 t=0 x t ) /T . 8: Round x̄ to an integral point by Algorithm 1.
The proof along with a detailed discussion on L\-convex functions is provided in the appendix.
Next, we describe the rounding process that takes an approximate solution to the relaxed problem
min x∈[1,N ]d f̃(x) (2)
to an approximate integral solution to the original problem (1). We list the pseudo-code in Algorithm 1. The following theorem proves the correctness of Algorithm 1 and estimates the computational complexity. Note that all upper bounds on computational complexity in this paper are proved to hold almost surely and thus also in expectation. Theorem 2. Suppose that Assumptions 1-4 hold. The solution returned by Algorithm 1 satisfies the ( , δ)-PAC guarantee. The computational complexity of Algorithm 1 is at most
O
[ d
2 log
( d
δ
) + d ] = Õ [ d
2 log
( 1
δ
)] .
Although we transform the original problem to a continuous convex problem, the algorithm design in this work relies heavily on other two important properties of L\-convex functions. The first property is that the global minima are integral points, which makes it possible to design strongly polynomial algorithms; see Section 4.2. The second property is that the Lovász extension is a piecewise linear function. When the sub-optimality gap is known, we prove that the Lovász extension satisfies the Weak Sharp Minima property and the subgradient descent method can be accelerated to achieve a better dependence on the scale N ; see the appendix.

heading: 3.2 Algorithm and complexity analysis
text: Utilizing the tools derived in the last subsection, we can first apply the stochastic subgradient (SSGD) method to find an approximate solution to problem (2) and then round the solution to get an approximate solution to problem (1). The main difficulty of giving sharp upper bounds lies in the fact that the Lovász extension is neither smooth nor strongly-convex. This property of the Lovász extension prohibits the application of Nesterov’s acceleration method (Nesterov, 1983) and variance reduction methods (Johnson & Zhang, 2013; Defazio et al., 2014).
Now, we propose the projected and truncated SSGD method for the ( , δ)-PAC guarantee. The orthogonal projection onto the feasible set, which is defined as
P(x) := (x ∧N1) ∨ 1, ∀x ∈ Rd,
is applied after each iteration to ensure the feasibility of iteration point. In addition to the projection, the component-wise truncation of stochastic subgradient is critical in reducing expected computational complexities. The truncation operator with threshold M > 0 is defined as
TM (g) := (g ∧M1) ∨ (−M1), ∀g ∈ Rd. The pseudo-code of projected and truncated SSGD method is listed in Algorithm 2. With a suitable choice of the step size, the truncation threshold and the number of iterations, Algorithm 2 returns an ( , δ)-PAC solution and the computational complexity scales as O(d2N2). Theorem 3. Suppose that Assumptions 1-4 hold. If we choose
T = Θ̃
[ dN2
2 log(
1 δ )
] , M = Θ̃ [√ log( dNT ) ] , η = N
M √ T ,
then Algorithm 2 returns an ( , δ)-PAC solution. Furthermore, we have
T ( , δ,MC) = O
[ d2N2
2 log(
1 δ ) +
d3N2
2 log
( d2N
3
) + d3N2L2
2
] = Õ [ d2N2
2 log(
1 δ )
] .
The intuition behind the truncation operation is that the truncation makes the norm of the subgradient smaller, which in turn allows the choice of a larger step size. With the larger step size, the convergence rate can be improved. Using the truncation operation to accelerate the algorithms is novel in the context of the stochastic L\-convex function minimization problem. The truncation operation has been considered in other contexts, such as the sparse phase retrieval problem (Wang et al., 2017). We can prove that the computational complexity without the truncation step (i.e., M =∞) is
Õ
[ d3N2
2 log(
1 δ )
] .
Hence, the truncation of stochastic subgradient is necessary for reducing the asymptotic expected computational complexity. This relation has also been observed in our numerical experiments. While the error of the normal SSGD method only contains the optimization residual and the variance terms, the error of the truncated SSGD method has an extra bias term. We note that the bias term can be made arbitrarily small with high probability by choosing large enough M and utilizing the tail bound for sub-Gaussian random variables, and therefore the total error can be controlled similarly as the normal SSGD method. We mention that Algorithm 2 also applies when X is a general L\-convex set. Next, we consider a special case with additional condition that the norm of stochastic subgradients is bounded almost surely. We then show that this additional condition, whenever available, can contribute to further complexity reduction of our proposed algorithm. Assumption 5. There exist a constant G and an unbiased subgradient estimator ĝ such that
P [‖ĝ‖1 ≤ G] = 1. Moreover, the cost of generating each ĝ is at most β evaluations of F (·, ·).
We note that G and β may depend on d and N . This assumption is common when analyzing the high-probability error bound of stochastic subgradient methods (Hazan & Kale, 2011; Xu et al., 2016). We first give two examples where Assumption 5 holds. Example 1. Suppose that the randomness of each feasible x ∈ X shares the same measure space, i.e., there exists a measure space (Z,BZ) such that ξx can be any element in the measure space for all x ∈ X . Moreover, for any fixed ξ ∈ BZ, the function F (·, ξ) is also L\-convex and has `∞-Lipschitz constant L̃. Then, the properties of L\-convex functions imply that there exists an unbiased subgradient estimator ĝ such that the cost of evaluating ĝ is d + 1 and ‖ĝ‖1 ≤ 3L̃/2. Therefore, the Assumption 5 holds with G = 3L̃/2 and β = d + 1. We note that this example is similar to the case when k = d in Ito (2019) while the noise of our problem can be unbounded.
When the distribution at each feasible point is bounded, we show that Assumption 5 also holds. Example 2. Suppose that the distribution at each point x ∈ X is bounded, namely, there exists a constant C ≥ 0 such that Pξx [|F (x, ξx)| ≤ C] = 1, ∀x ∈ X . We note that bounded distributions are sub-Gaussian with parameter 2C2. In this case, there exists an unbiased subgradient estimator such that Assumption 5 holds with G = 2dC and β = 2.
More details about these two examples are provided in the appendix. Since the stochastic subgradient is bounded, the truncation step is unnecessary in Algorithm 2. The computational complexity of Algorithm 2 is estimated in the following theorem. The proof is similar to Lemma 10 in Hazan & Kale (2011) and, since the feasible set is the hypercube [1, N ]d, we replace the `2-norm with the `∞-norm to bound distances between points. Theorem 4. Suppose that Assumptions 1-5 hold. If we skip the truncation step in Algorithm 2 (i.e., set M =∞) and choose
T = Θ̃
[ N2(L+G)2
2 log(
1 δ )
] , η = √ dN2
TG2 ,
then Algorithm 2 returns an ( , δ)-PAC solution. Furthermore, we have T ( , δ,MC) = O [ βN2(L+G)2 + d
2 log(
1 δ ) +
d2G2
2 ] = Õ [ βN2(L+G)2 + d
2 log(
1 δ )
] .
In the case of Example 1, we have β = d+ 1, G = 3L̃/2 and then the computational complexity of Algorithm 2 is at most Õ[dN2(L+ L̃)2 −2 log(1/δ)]. If both Lipschitz constants are independent of d and N , the computational complexity becomes Õ[dN2 −2 log(1/δ)], which is O(d) better than the general case without Assumption 5. Moreover, in the case of Example 2, we have G = d and β = 2. Hence, the computational complexity is at most Õ[d2N2 −2 log(1/δ)], even without the truncation.
Finally, we consider the case when the global optimum is unique and the sub-optimality gap is known to be at least c > 0. In this case, the Lovász extension satisfies the weak sharp minimum condition (Burke & Ferris, 1993) and we can design an acceleration scheme similar to that in Xu et al. (2016). The acceleration scheme reduces the dependence on N from N2 to log(N) and we provide detailed discussions in the appendix.

heading: 4 Stochastic Localization Algorithms
text: Although the algorithms developed in Section 3 are guaranteed to find an approximate solution in polynomial time, gradient-based algorithms require the prior knowledge about the Lipschitz constant of the objective function to choose the optimal step size. In addition, the computational complexity of subgradient methods has a quadratic dependence on the Lipschitz constant. However, for certain applications, the upper bound on the Lipschitz constant may be difficult to estimate a priori or is very large so that the computational complexity is prohibitive for large-scale problems. To deal with those difficult scenarios, we utilize the discrete natural of problem (1) to design algorithms that do not require the prior knowledge about the objective function except the L\-convexity and whose computational complexity has an upper bound that is independent of the Lipschitz constant. More specifically, we propose two approaches to deal with this problem. The first approach is to generalize the uniform sampling algorithm in Zhang et al. (2021) to the multi-dimensional case. In the second approach, we propose a strongly polynomial approximate algorithm by gradually reducing the dimension.

heading: 4.1 Multi-dimensional uniform sampling algorithm
text: We first generalize the uniform sampling algorithm in Zhang et al. (2021) to the multi-dimensional case. The basic idea is that if the computational complexity for finding an ( , δ)-PAC solution has the form C −2 log(1/δ), we can view the algorithm as a (possibly biased) estimator with sub-Gaussian tails to the optimal value, where C > 0 is a constant. For example, the uniform sampling algorithm provides a sub-Gaussian estimator for the optimal value of problem (1) in the one-dimensional case. Now suppose that we have derived a sub-Gaussian estimator for the optimal value for the (d− 1)-dimensional case. Consider the function
f1(x1) := min x2,...,xd∈[N ] f(x1, . . . , xd), ∀x1 ∈ [N ], (3)
which is an one-dimensional convex function. Since the (d− 1)-dimensional algorithm has a subGaussian tail, we can view each implementation of the algorithm as a noisy function value oracle to f1(·) and the uniform sampling algorithm can be modified to solve the convex problem (3). By a similar method as Zhang et al. (2021), we can prove that the computational complexity of the resulting d-dimensional algorithm is a constant multiple of that of the (d− 1)-dimensional algorithm. By the induction method, we get an approximate algorithm for the multi-dimensional problem (1) whose computational complexity is independent of the Lipschitz constant L and the scale N in the high-probability regime. A rigorous proof to the following result is provided in the appendix. Theorem 5. Suppose that Assumptions 1-3 hold. There exists an [( , δ)-PAC,MC]-algorithm with computational complexity
T ( , δ,MC) = Õ [ Cd
2 log(
1 δ )
] ,
where C is an absolute constant.
The exponential dependence on the dimension d makes it impossible to use the generalized algorithm even for moderate dimension problems. However, this result implies that the computational complexity almost does not grow with N in the regime N d in the asymptotic regime with the failing probability δ is very small. Therefore this algorithm can be preferable when d is small and N is large.

heading: 4.2 Dimension reduction algorithm
text: Next, we take an alternative approach and derive a strongly polynomial approximate algorithm for problem (1). In this work, we do not treat the parameters d and N as input parameters. The message we hope to convey is that the sample complexities are upper bounded by a constant independent of the objective function. This property is important for some applications of the stochastic L\-function minimization problem where the objective function has unfavorable properties (e.g. large Lipschitz constant, small sub-optimality gap). Such unfavorable properties may not be known a priori before the problem is solved. If the parameters d andN is treated as an input parameter, the proposed algorithms are pseudo-polynomial. We observe that a subgradient of a convex function provides a separation oracle to the global minima, which can be used to gradually shrink the set that contains global minima. This idea leads to the deterministic cutting-plane methods (Vaidya, 1996; Bertsimas & Vempala, 2004; Lee et al., 2015). In this work, we develop a framework for extending cutting-plane methods to the stochastic case. We prove that if we average Õ[dN2 −2 log(1/δ)] stochastic subgradients, we can get a separation oracle for the set of -approximate solutions with probability at least 1− δ. Utilizing this stochastic separation oracle, we may extend the deterministic cutting-plane methods to the stochastic case to find PAC solutions.
However, the stochastic cutting-plane methods still require the knowledge of the Lipschitz constant of the objective function. To design a strongly polynomial algorithm, we observe that if a convex body P ⊂ Rd has a volume smaller than (d!)−1 = O[exp(−d log(d))], then all integral points in P must lie on a hyperplane. Otherwise, if there exist d+ 1 points x0, . . . , xd ∈ P that are not on the same hyperplane, then the convex body P contains the polytope conv{x0, . . . , xd}, which has the volume
1 d! |det(x1 − x0, . . . , xd − x0)| ≥ 1 d! > vol(P ),
where vol(P ) is the volume of the polytope P . This leads to a contradiction. Hence, we may use the stochastic version of Vaidya’s method (Vaidya, 1996) or the random walk-based (Bertsimas & Vempala, 2004) cutting-plane method to reduce the volume of the search polytope P to O[exp(−d log(d))] and then reduce the problem dimension by projecting the polytope onto the hyperplane. After d − 1 dimension reductions, we have an one-dimensional convex problem and algorithms in Zhang et al. (2021) can be applied. Theorem 6. Suppose that Assumptions 1-3 hold. Then, the dimension reduction algorithm is an [( , δ)-PAC,MC]-algorithm with computational complexity
T ( , δ,MC) = Õ [ d3N2(d+ log(N))
2 log(
1 δ )
] .
The dimension reduction method is the first algorithm for the stochastic L\-convex function minimization that is strongly polynomial and does not require prior knowledge about the objective function
except the L\-convexity. Recently, Jiang (2020) proved that the LLL algorithm (Lenstra et al., 1982) is able to reduce the number of arithmetic operations to be polynomial and make the algorithm more practical, but the LLL algorithm will not reduce the number of evaluations to the stochastic objective values. Their results can be naturally extended to the stochastic case to achieve the PAC guarantee. We note that if we allow exponentially many arithmetic operations, the LLL algorithm is not necessary and the computational complexity can be reduced to Õ[d4N2 −2 log(1/δ)]. Similar to the gradient-based method, the acceleration scheme can reduce the number of required objective function value evaluations if we have a positive lower bound to the sub-optimality gap.

heading: 5 Lower Bounds
text: In this section, we provide a lower bound result on [( , δ)-PAC,MC]-algorithms. This result shows the limit of algorithms that achieve the PAC guarantee for any stochastic L\-convex functions. To prove lower bounds, we make non-trivial extensions to the construction in Graur et al. (2020) to get several convex models that are “similar” to each other but they have distinct optimal solutions, where the similarity between two models is characterized by the KL-divergence between their distributions. Hence, algorithms need a large number of evaluations to differentiate these models. More rigorously, the information-theoretical inequality in Kaufmann et al. (2016) is utilized to get the following lower bound:
Theorem 7. Suppose Assumptions 1-3 hold. We have T ( , δ,MC) ≥ Θ [ d
2 log(
1 δ )
] .
Comparing with the lower bound for the expected error in Ito (2019), our high-probability lower bound is O(d) smaller but includes an extra term log(1/δ). Since the failing probability δ can be sufficiently small, the two lower bounds are not directly comparable. The results in Sections 3-4 imply that
T ( , δ,MC) = Õ [
min{d2N2,Md} 2 log( 1
δ )
] .
Therefore, we conjecture that the optimal computational complexity is asymptotically independent of N and leave the design of optimal algorithms for future works.

heading: 6 Numerical Experiments
text: In this section, we present a brief account of numerical experiments to illustrate several of our theoretical findings. First, we show that the extra new truncation step in the gradient-based method can achieve a better performance, compared to the gradient-based method without the truncation step. We implement Algorithm 2 with and without the truncation operation on a staffing allocation stochastic optimization problem with staffing cost constraint, which can be formulated into the regularized and the constrained forms. Details about the problem setup are provided in the appendix. We summarize the computational complexities and the objective values in Table 2. We can see that the two algorithms achieve similar objective values but the “Truncated" version observes an advantage in the computational cost compared to the "Not Truncated" version. We note that our algorithms do not require the feasible region to be a hypercube and in fact the problems considered in this section have more general feasible set with constraints. Especially, we can see that the truncation plays an important role in reducing the computational complexity, especially when the dimension is high.
Next, we compare the performances of the gradient-based methods and stochastic localization algorithms on multi-dimensional problems. To facilitate the comparison, we construct a separable stochastic convex function, which is an L\-convex function by the results in Murota (2003). We have observed that using (N /4, δ/4)-SO oracles is enough for producing high-probability guarantees using the dimension reduction method within the range 1 ≤ N ≤ 150. More details about the experiment setup are provided in the appendix. We summarize the results in Table 3, where the coverage rate refers to the percentage of implementations that produce an -optimal solution. Since the coverage rates of the algorithms are 100%, the PAC guarantee is satisfied by all of the algorithms except the (d,N) = (20, 50) case for the dimension reduction method. This failure is because of
our loosened choice of the separation oracles. When the dimension d is not large, the dimension reduction method has the best performance and has the advantage of not requiring any knowledge about the Lipschitz constant; when the dimension d is large, the higher order dependence on d makes the dimension reduction method inferior to the stochastic subgradient method. In addition, the computational complexity of the dimension reduction method has little dependence on N . This is because we choose the separation oracle to be (N /4, δ/4)-SO and the computational complexity of each evaluation is independent of N .

heading: 7 Conclusions
text: In this work, we provide the first polynomial-time algorithms for the stochastic L\-convex minimization problem, which serves as an important extension of the stochastic submodular minimization problem. Comparing to existing stochastic submodular minimization algorithms, our algorithms either achieve a lower computational complexity or do not require extra information about the objective function. We also derive a lower bound on the computational complexity to achieve the PAC guarantee. Numerical results are implemented to verify our theoretical findings.
Acknowledgments and Disclosure of Funding
Haixiang Zhang and Javad Lavaei were supported by grants from AFOSR, ARO, ONR, NSF and C3.ai Digital Transformation Institute.

Review summaries:
This paper focuses on the problem of minimizing $f(x) = \mathbb{E}_{\xi}[F(x,\xi_x)]$ where $f:[N]^d\to\mathbb{R}$ is an $L^\natural$ function. This class of functions can be interpreted as a discrete analog of convex functions. The goal is to design polynomial time algorithms that guarantee a solution whose value is $\epsilon$ away from the optimal value with probability at least $1-\delta$ (in the paper referred as $(\epsilon,\delta)$-PAC guarantee).  The authors provide upper bounds and a lower bound for the computational complexity of the problem, which is defined as the worst-case expected time over all possible $(\epsilon,\delta)$-PAC algorithms. First, the authors present the standard Lovasz extension for $L^\natural$ functions. This extension is convex so common techniques in convex optimization can be used. First, they present a truncated stochastic subgradient algorithm that has an appropriate rounding procedure as a final step. They show this method has a complexity of $\tilde{O}(\frac{d^2N^2}{\epsilon^2}\log(1/\delta))$ when $f$ is $L$-Lipschitz. Later, they show that under the assumption that the subgradients of the function are uniformly bounded with probability 1, the same algorithm without truncation has an improved complexity. The authors also study the performance of the algorithms when the suboptimality gap is known. Then, they give two methods for the case when the Lipschitz constant is not known: multidimensional uniform sampling and dimension reduction. Finally, they give an information-theoretical lower bound of $\Theta(\frac{d}{\epsilon^2}\log(1/\delta))$. 

In this paper, the authors studied the problem of stochastic $L^{\natural}$-convex function minimization. This problem is a extension of submodular minimization. Authors come up with an algorithm that uses projected gradient descent with truncated gradient on lovasz extension. They also provide lower bound for the complexity of this problem.

This paper studies the minimization of a stochastic L$^\\natural$-convex function. - The authors analyze the convergence of three algorithms, truncated subgradient, multi-dim uniform sampling, and stochastic cutting-plane. The time complexity of these algorithms is typically $O(d^2N^2/\\epsilon^2)$ or something higher, where $d$ is the dimension, $N$ is the range of the feasible region, and $\\epsilon$ is the error-accuracy. - The authors also proved a lower bound that any algorithm takes $\\Omega(d/\\epsilon^2)$-time to find an $\\epsilon$-optimal solution. - In numerical experiments, the truncated subgradient method rapidly found an $\\epsilon$-optimal solution.