A Proof of Lemma 1

Stochastic Submodular Maximization: The Case of Coverage Functions

Abstract

Stochastic optimization of continuous objectives is at the heart of modern machine learning. However, many important problems are of discrete nature and often involve submodular objectives. We seek to unleash the power of stochastic continuous optimization, namely stochastic gradient descent and its variants, to such discrete problems. We first introduce the problem of stochastic submodular optimization, where one needs to optimize a submodular objective which is given as an expectation. Our model captures situations where the discrete objective arises as an empirical risk (e.g., in the case of exemplar-based clustering), or is given as an explicit stochastic model (e.g., in the case of influence maximization in social networks). By exploiting that common extensions act linearly on the class of submodular functions, we employ projected stochastic gradient ascent and its variants in the continuous domain, and perform rounding to obtain discrete solutions. We focus on the rich and widely used family of weighted coverage functions. We show that our approach yields solutions that are guaranteed to match the optimal approximation guarantees, while reducing the computational cost by several orders of magnitude, as we demonstrate empirically.

\newcolumntype

H>c<@

1 Introduction

Submodular functions are discrete analogs of convex functions. They arise naturally in many areas, such as the study of graphs, matroids, covering problems, and facility location problems. These functions are extensively studied in operations research and combinatorial optimization \citepkrause2012submodular. Recently, submodular functions have proven to be key concepts in other areas such as machine learning, algorithmic game theory, and social sciences. As such, they have been applied to a host of important problems such as modeling valuation functions in combinatorial auctions, feature and variable selection [Krause and Guestrin(2005a)], data summarization [Lin and Bilmes(2011)], and influence maximization [Kempe et al.(2003)Kempe, Kleinberg, and Tardos].

Classical results in submodular optimization consider the oracle model whereby the access to the optimization objective is provided through a black box — an oracle. However, in many applications, the objective has to be estimated from data and is subject to stochastic fluctuations. In other cases the value of the objective may only be obtained through simulation. As such, the exact computation might not be feasible due to statistical or computational constraints. As a concrete example, consider the problem of influence maximization in social networks [Kempe et al.(2003)Kempe, Kleinberg, and Tardos]. The objective function is defined as the expectation of a stochastic process, quantifying the size of the (random) subset of nodes influenced from a selected seed set. This expectation cannot be computed efficiently, and is typically approximated via random sampling, which introduces an error in the estimate of the value of a seed set. Another practical example is the exemplar-based clustering problem, which is an instance of the facility location problem. Here, the objective is the sum of similarities of all the points inside a (large) collection of data points to a selected set of centers. Given a distribution over point locations, the true objective is defined as the expected value w.r.t. this distribution, and can only be approximated as a sample average. Moreover, evaluating the function on a sample involves computation of many pairwise similarities, which is computationally prohibitive in the context of massive data sets.

In this work, we provide a formalization of such stochastic submodular maximization tasks. More precisely, we consider set functions , defined as for , where is an arbitrary distribution and for each realization , the set function is monotone and submodular (hence is monotone submodular). The goal is to maximize subject to some constraints (e.g. the -cardinality constraint) having access only to i.i.d. samples .

Methods for submodular maximization fall into two major categories: (i) The classic approach is to directly optimize the objective using discrete optimization methods (e.g. the Greedy algorithm and its accelerated variants), which are state-of-the-art algorithms (both in practice and theory), at least in the case of simple constraints, and are most widely considered in the literature; (ii) The alternative is to lift the problem into a continuous domain and exploit continuous optimization techniques available therein \citepvondrak-pipage. While the continuous approaches may lead to provably good results, even for more complex constraints, their high computational complexity inhibits their practicality.

In this paper we demonstrate how modern stochastic optimization techniques (such as SGD, AdaGrad [Duchi et al.(2011)Duchi, Hazan, and Singer] and Adam [Kingma and Ba(2015)]), can be used to solve an important class of discrete optimization problems which can be modeled using weighted coverage functions. In particular, we show how to efficiently maximize them under matroid constraints by (i) lifting the problem into the continuous domain using the multilinear extension [Vondrák(2008)], (ii) efficiently computing a concave relaxation of the multilinear extension [Seeman and Singer(2013)], (iii) efficiently computing an unbiased estimate of the gradient for the concave relaxation thus enabling (projected) stochastic gradient ascent-style algorithms to maximize the concave relaxation, and (iv) rounding the resulting fractional solution without loss of approximation quality [Calinescu et al.(2011)Calinescu, Chekuri, Pál, and Vondrák]. In addition to providing convergence and approximation guarantees, we demonstrate that our algorithms enjoy strong empirical performance, often achieving an order of magnitude speedup with less than error with respect to Greedy. As a result, the presented approach unleashes the powerful toolkit of stochastic gradient based approaches to discrete optimization problems.

Our contributions.

In this paper we (i) introduce a framework for stochastic submodular optimization, (ii) provide a general methodology for constrained maximization of stochastic submodular objectives, (iii) prove that the proposed approach guarantees a approximation in expectation for the class of weighted coverage functions, which is the best approximation guarantee achievable in polynomial time unless , (iv) highlight the practical benefit and efficiency of using continuous-based stochastic optimization techniques for submodular maximization, (v) demonstrate the practical utility of the proposed framework in an extensive experimental evaluation. We show for the first time that continuous optimization is a highly practical, scalable avenue for maximizing submodular set functions.

2 Background and problem formulation

Let be a ground set of elements. A set function is submodular if for every , it holds . Function is said to be monotone if for all . We focus on maximizing subject to some constraints on . The prototypical example is maximization under the cardinality constraint, i.e., for a given integer , find , , which maximizes . Finding an exact solution for monotone submodular functions is NP-hard [Feige(1998)], but a -approximation can be efficiently determined [Nemhauser et al.(1978)Nemhauser, Wolsey, and Fisher]. Going beyond the -approximation is NP-hard for many classes of submodular functions [Nemhauser et al.(1978)Nemhauser, Wolsey, and Fisher, Krause and Guestrin(2005b)]. More generally, one may consider matroid constraints, whereby is a matroid with the family of independent sets , and maximize such that . The Greedy algorithm achieves a -approximation [Fisher et al.(1978)Fisher, Nemhauser, and Wolsey], but Continuous Greedy introduced by \citetvondrak2008optimal, calinescu2007maximizing can achieve a -optimal solution in expectation. Their approach is based on the multilinear extension of , , defined as

(1)

for all . In other words, is the expected value of of over sets wherein each element is included with probability independently. Then, instead of optimizing over , we can optimize over the matroid base polytope corresponding to : , where is the matroid’s rank function. The Continuous Greedy algorithm then finds a solution which provides a approximation. Finally, the continuous solution is then efficiently rounded to a feasible discrete solution without loss in objective value, using Pipage Rounding [Ageev and Sviridenko(2004), Calinescu et al.(2007)Calinescu, Chekuri, Pál, and Vondrák]. The idea of converting a discrete optimization problem into a continuous one was first exploited by \citetlovasz1983submodular in the context of submodular minimization and this approach was recently applied to a variety of problems [Vondrák(2007), Iyer and Bilmes(2015), Bach(2010)].

Problem formulation.

The aforementioned results are based on the oracle model, whereby the exact value of for any is given by an oracle. In absence of such an oracle, we face the additional challenges of evaluating , both statistical and computational. In particular, consider set functions that are defined as expectations, i.e. for we have

(2)

where is an arbitrary distribution and for each realization , the set function is submodular. The goal is to efficiently maximize subject to constraints such as the -cardinality constraint, or more generally, a matroid constraint.

As a motivating example, consider the problem of propagation of contagions through a network. The objective is to identify the most influential seed set of a given size. A propagation instance (concrete realization of a contagion) is specified by a graph . The influence of a set of nodes in instance is the fraction of nodes reachable from using the edges . To handle uncertainties in the concrete realization, it is natural to introduce a probabilistic model such as the Independent Cascade [Kempe et al.(2003)Kempe, Kleinberg, and Tardos] model which defines a distribution over instances that share a set of nodes. The influence of a seed set is then the expectation , which is a monotone submodular function. Hence, estimating the expected influence is computationally demanding, as it requires summing over exponentially many functions . Assuming as in (2), one can easily obtain an unbiased estimate of for a fixed set by random sampling according to . The critical question is, given that the underlying function is an expectation, can we optimize it more efficiently?

Our approach is based on continuous extensions that are linear operators on the class of set functions, namely, linear continuous extensions. As a specific example, considering the multilinear extension, we can write where denotes the extension of . As a consequence, the value of , when , is an unbiased estimator for and unbiased estimates of the (sub)gradients may be obtained analogously. We explore this avenue to develop efficient algorithms for maximizing an important subclass of submodular functions that can be expressed as weighted coverage functions. Our approach harnesses a concave relaxation detailed in Section 3.

Further related work. The emergence of new applications, combined with a massive increase in the amount of data has created a demand for fast algorithms for submodular optimization. A variety of approximation algorithms have been presented, ranging from submodular maximization subject to a cardinality constraint [Mirzasoleiman et al.(2015)Mirzasoleiman, Badanidiyuru, Karbasi, Vondrak, and Krause, Wei et al.(2014)Wei, Iyer, and Bilmes, Badanidiyuru and Vondrák(2014)], submodular maximization subject to a matroid constraint [Calinescu et al.(2007)Calinescu, Chekuri, Pál, and Vondrák], non-monotone submodular maximization [Feige et al.(2011)Feige, Mirrokni, and Vondrak], approximately submodular functions [Horel and Singer(2016)], and algorithms for submodular maximization subject to a wide variety of constraints [Kulik et al.(2009)Kulik, Shachnai, and Tamir, Feldman et al.(2011)Feldman, Naor, and Schwartz, Vondrák(2013), Iyer and Bilmes(2013), Ene and Nguyen(2016)]. A closely related setting to ours is online submodular maximization [Streeter and Golovin(2008)], where functions come one at a time and the goal is to provide time-dependent solutions (sets) such that a cumulative regret is minimized. In contrast, our goal is to find a single (time-independent) set that maximizes the objective (2). Another relevant setting is noisy submodular maximization, where the evaluations returned by the oracle are noisy [Hassidim and Singer(2016), Singla et al.(2016)Singla, Tschiatschek, and Krause]. Specifically, [Singla et al.(2016)Singla, Tschiatschek, and Krause] assumes a noisy but unbiased oracle (with an independent sub-Gaussian noise) which allows one to sufficiently estimate the marginal gains of items by averaging. In the context of cardinality constraints, some of these ideas can be carried to our setting by introducing additional assumptions on how the values vary w.r.t. to their expectation . However, we provide a different approach that does not rely on uniform convergence and compare sample and running time complexity comparison with variants of Greedy in Section 3.

3 Stochastic Submodular Optimization

We follow the general framework of [Vondrák(2008)] whereby the problem is lifted into the continuous domain, a continuous optimization algorithm is designed to maximize the transferred objective, and the resulting solution is rounded. Maximizing subject to a matroid constraint can then be done by first maximizing its multilinear extension over the matroid base polytope and then rounding the solution. Methods such as the projected stochastic gradient ascent can be used to maximize over this polytope.

Critically, we have to assure that the computed local optima are good in expectation. Unfortunately, the multilinear extension lacks concavity and therefore may have bad local optima. Hence, we consider concave continuous extensions of that are efficiently computable, and at most a constant factor away from to ensure solution quality. As a result, such a concave extension could then be efficiently maximized over a polytope using projected stochastic gradient ascent which would enable the application of modern continuous optimization techniques. One class of important functions for which such an extension can be efficiently computed is the class of weighted coverage functions.

The class of weighted coverage functions (WCF).

Let be a set and let be a nonnegative modular function on , i.e. , . Let be a collection of subsets of . The weighted coverage function defined as

is monotone submodular. For all , let us denote by and by the indicator function. The multilinear extension of can be expressed in a more compact way:

(3)

where we used the fact that each element was chosen with probability .

Concave upper bound for weighted coverage functions.

To efficiently compute a concave upper bound on the multilinear extension we use the framework of \citetYaron-cover. Given that all the weights , in (3) are non-negative, we can construct a concave upper bound for the multilinear extension using the following Lemma. Proofs can be found in the Appendix A.

Lemma 1.

For define . Then the Fenchel concave biconjugate of is . Also

Furthermore, is an extension of , i.e. : .

Consequently, given a weighted coverage function with represented as in (3), we can define

(4)

and conclude using Lemma 1 that , as desired. Furthermore, has three interesting properties: (1) It is a concave function over , (2) it is equal to on vertices of the hypercube, i.e. for one has , and (3) it can be computed efficiently and deterministically given access to the sets , . In other words, we can compute the value of using at most operations. Note that is not the tightest concave upper bound of , even though we use the tightest concave upper bounds for each term of .

Optimizing the concave upper bound by stochastic gradient ascent.

Instead of maximizing over a polytope , one can now attempt to maximize over . Critically, this task can be done efficiently, as is concave, by using projected stochastic gradient ascent. In particular, one can control the convergence speed by choosing from the toolbox of modern continuous optimization algorithms, such as Sgd, AdaGrad and Adam. Let us denote a maximizer of over by , and also a maximizer of over by . We can thus write

which is the exact guarantee that previous methods give, and in general is the best near-optimality ratio that one can give in poly-time. Finally, to round the continuous solution we may apply Randomized-Pipage-Rounding [Calinescu et al.(2011)Calinescu, Chekuri, Pál, and Vondrák] as the quality of the approximation is preserved in expectation.

0:  matroid with base polytope , (step size), (maximum # of iterations)
1:   starting point in
2:  for  to  do
3:     Choose at random from a distribution such that
4:     
5:     
6:  end for
7:  
8:   Randomized-Pipage-Round
9:  return such that , .
Algorithm 1 Stochastic Submodular Maximization via concave relaxation

Matroid constraints.

Constrained optimization can be efficiently performed by projected gradient ascent whereby after each step of the stochastic ascent, we need to project the solution back onto the feasible set. For the case of matroid constraints, it is sufficient to consider projection onto the matroid base polytope. This problem of projecting on the base polytope has been widely studied and fast algorithms exist in many cases [Bach et al.(2013), Brucker(1984), Pardalos and Kovoor(1990)]. While these projection algorithms were used as a key subprocedure in constrained submodular minimization, here we consider them for submodular maximization. Details of a fast projection algorithm for the problems considered in this work are presented the Appendix D. Algorithm 1 summarizes all steps required to maximize subject to matroid constraints.

Convergence rate.

Since we are maximizing a concave function over a matroid base polytope , convergence rate (and hence running time) depends on , as well as maximum gradient norm (i.e. with probability ). 1 In the case of the base polytope for a matroid of rank , is , since each vertex of the polytope has exactly ones. Also, from (4), one can build a rough upper bound for the norm of the gradient:

which depends on the weights as well as and is hence problem-dependent. We will provide tighter upper bounds for gradient norm in our specific examples in the later sections. With , and classic results for SGD [Shalev-Shwartz and Ben-David(2014)], we have that

where is the total number of SGD iterations and is the final outcome of SGD (see Algorithm 1). Therefore, for a given , after iterations, we have

Summing up, we will have the following theorem:

Theorem 2.

Let be a weighted coverage function, be the base polytope of a matroid , and and be as above. Then for each , Algorithm 1 after iterations, produces a set such that .

Remark.

Indeed this approximation ratio is the best ratio one can achieve, unless PNP [Feige(1998)]. A key point to make here is that our approach also works for more general constraints (in particular is efficient for simple matroids such as partition matroids). In the latter case, Greedy only gives -approximation and fast discrete methods like Stochastic-Greedy [Mirzasoleiman et al.(2015)Mirzasoleiman, Badanidiyuru, Karbasi, Vondrak, and Krause] do not apply, whereas our method still yields an -optimal solution.

Time Complexity.

One can compute an upper bound for the running time of Algorithm 1 by estimating the time required to perform gradient computations, projection on , and rounding. For the case of uniform matroids, projection and rounding take and time, respectively (see Appendix D). Furthermore, for the applications considered in this work, namely expected influence maximization and exemplar-based clustering, we provide linear time algorithms to compute the gradients. Also when our matroid is the -uniform matroid (i.e. -cardinality constraint), we have . By Theorem 2, the total computational complexity of our algorithm is .

Comparison to Greedy.

Let us relate our results to the classical approach. When running the Greedy algorithm in the stochastic setting, one estimates where are i.i.d. samples from . The following proposition bounds the sample and computational complexity of Greedy. The proof is detailed in the Appendix B.

Proposition 3.

Let be a submodular function defined as (2). Suppose for all and all . Assume denotes the optimal solution for subject to -cardinality constraint and denotes the solution computed by the greedy algorithm on after steps. Then, in order to guarantee

it is enough to have

i.i.d. samples from . The running time of Greedy is then bounded by

where is an upper bound on the computation time for a single evaluation of .

As an example, let us compare the worst-case complexity bound obtained for SGD (i.e. ) with that of Greedy for the influence maximization problem. Each single function evaluation for Greedy amounts to computing the total influence of a set in a sample graph, which makes (here we assume our sample graphs satisfy ). Also, a crude upper bound for the size of the gradient for each sample function is (see Appendix E.1). Hence, we can deduce that SGD can have a factor speedup w.r.t. to Greedy.

4 Applications

We will now show how to instantiate the stochastic submodular maximization framework using several prototypical discrete optimization problems.

Influence maximization.

As discussed in Section 2, the Independent Cascade [Kempe et al.(2003)Kempe, Kleinberg, and Tardos] model defines a distribution over instances that share a set of nodes. The influence of a set of nodes in instance is the fraction of nodes reachable from using the edges . The following Lemma shows that the influence belongs to the class of WCF.

Lemma 4.

The influence function is a WCF. Moreover,

(5)
(6)

where is the set of all nodes having a (directed) path to .

We return to the problem of maximizing given a distribution over graphs sharing nodes . Since is a weighted sum of submodular functions, it is submodular. Moreover,

Let be the uniform distribution over vertices. Then,

(7)

and the corresponding upper bound would be

(8)

This formulation proves to be helpful in efficient calculation of subgradients, as one can obtain a random subgradient in linear time. For more details see Appendix E.1. We also provide a more efficient, biased estimator of the expectation in the Appendix.

Facility location.

Let be a complete weighted bipartite graph with parts and and nonnegative weights . The weights can be considered as utilities or some similarity metric. We select a subset and each selects with the highest weight . Our goal is to maximize the average weight of these selected edges, i.e. to maximize

(9)

given some constraints on . This problem is indeed the Facility Location problem, if one takes to be the set of facilities and to be the set of customers and to be the utility of facility for customer . Another interesting instance is the Exemplar-based Clustering problem, in which is a set of objects and is the similarity (or inverted distance) between objects and , and one tries to find a subset of exemplars (i.e. centroids) for these objects.

The stochastic nature of this problem is revealed when one writes (9) as the expectation , where is the uniform distribution over and . One can also consider this more general case, where ’s are drawn from an unknown distribution, and one tries to maximize the aforementioned expectation.

First, we claim that for each is again a weighted coverage function. For simplicity, let and set , with and .

Lemma 5.

The utility function is a WCF. Moreover,

(10)
(11)

We remark that the gradient of both and can be computed in linear time using a recursive procedure. We refer to Appendix E.2 for more details.

5 Experimental Results

We demonstrate the practical utility of the proposed framework and compare it to standard baselines. We compare the performance of the algorithms in terms of their wall-clock running time and the obtained utility. We consider the following problems:

Figure 1: In the case of Facility location for Blog selection as well as on influence maximization on Epinions, the proposed approach reaches the same utility significantly faster. On the exemplar-based clustering of CIFAR, the proposed approach is outperformed by Stochastic-Greedy, but nevertheless reaches of the Greedy utility in a few seconds (after less than iterations). On Influence Maximization over partition matroids, the proposed approach significantly outperforms Greedy.
  • Influence Maximization for the Epinions network2. The network consists of 75 879 nodes and 508 837 directed edges. We consider the subgraph induced by the top 10 000 nodes with the largest out-degree and use the independent cascade model [Kempe et al.(2003)Kempe, Kleinberg, and Tardos]. The diffusion model is specified by a fixed probability for each node to influence its neighbors in the underlying graph. We set this probability to be , and chose the number of seeds .

  • Facility Location for Blog Selection. We use the data set used in [Glance et al.(2005)Glance, Hurst, Nigam, Siegler, Stockton, and Tomokiyo], consisting of 45 193 blogs, and 16 551 cascades. The goal is to detect information cascades/stories spreading over the blogosphere. This dataset is heavy-tailed, hence a small random sample of the events has high variance in terms of the cascade sizes. We set .

  • Exemplar-based Clustering on CIFAR-10. The data set contains 60 000 color images with resolution . We use a single batch of 10 000 images and compare our algorithms to variants of Greedy over the full data set. We use the Euclidean norm as the distance function and set . Further details about preprocessing of the data as well as formulation of the submodular function can be found in Appendix E.3.

Baselines.

In the case of cardinality constraints, we compare our stochastic continuous optimization approach against the most efficient discrete approaches (Lazy-)Greedy and (Lazy-)Stochastic-Greedy, which both provide optimal approximation guarantees. For Stochastic-Greedy, we vary the parameter in order to explore the running time/utility tradeoff. We also report the performance of randomly selected sets. For the two facility location problems, when applying the greedy variants we can evaluate the exact objective (true expectation). In the Influence Maximization application, computing the exact expectation is intractable. Hence, we use an empirical average of samples (cascades) from the model. We note that the number of samples suggested by Proposition 3 is overly conservative, and instead we make a practical choice of samples.

Results.

The results are summarized in Figure 1. On the blog selection and influence maximization applications, the proposed continuous optimization approach outperforms Stochastic-Greedy in terms of the running time/utility tradeoff. In particular, for blog selection we can compute a solution with the same utility faster than Stochastic-Greedy with . Similarly, for influence maximization on Epinions we the solution faster than Stochastic-Greedy with . On the exemplar-based clustering application Stochastic-Greedy outperforms the proposed approach. We note that the proposed approach is still competitive as it recovers of the value after less than thousand iterations.

We also include an experiment on Influence Maximization over partition matroids for the Epinions network. In this case, Greedy only provides a approximation guarantee and Stochastic-Greedy does not apply. To create the partition, we first sorted all the vertices by their out-degree. Using this order on the vertices, we divided the vertices into two partitions, one containing vertices with even positions, other containing the rest. Figure 1 clearly demonstrates that the proposed approach outperforms Greedy in terms of utility (as well as running time).

Acknowledgments

The research was partially supported by ERC StG 307036. We would like to thank Yaron Singer for helpful comments and suggestions.

Appendix A Proof of Lemma 1

Here we prove the inequality mentioned in the lemma. Proof of the fact of being Fenchel biconjugate is in Appendix F.

We prove the left-hand-side inequality, since the right-hand-side inequality is a consequence of Fenchel biconjugate-ness.

Let . We note from the inequality that . We thus obtain

Now, if then the result is clear. Also, if , then we note that the function is decreasing for , and hence, . The left-hand-side inequality thus follows immediately.

Appendix B Proof of Proposition 3

Note that the total number of subsets of cardinality less than is bounded from above by . For each such set we want the estimate to be at most away from . Also, note that the function is itself a submodular function and maximizing it would give a -approximation to its optimum. Hence, it is enough to have enough samples such that for all subsets of cardinality at most the two values and differ by at most epsilon. By using Hoeffding’s inequality and a union bound over all the subsets of cardinality at most (note that ) we get the result.

Appendix C Proof of Lemmas 4 and 5

c.1 Lemma 4

Proof.

Let , where is the set of vertices reachable from . By construction, there is a one-to-one correspondence between elements of and , namely . For , let be its corresponding subset in , i.e. . It’s obvious that . Setting , makes a WCF. But , so is also a WCF.

Moreover, for each , the set is the set of all elements of that contain , which are precisely those vertices from which there is a (directed) path to . We also relax our notation, and replace any element of by its correspondent in . Hence,

which are poly-time computable since one can find with a simple BFS algorithm in for each .

c.2 Lemma 5

Proof.

Write instead of .

Let , where , and let (set ). Note that there is a natural bijection between and , namely . Let be the modular function with weights , defined on , and define the WCF as

(12)

Since ’s are forming a decreasing chain, and (12) becomes

which is exactly .

Furthermore, is simply the set . Hence, we can write the multilinear extension and the corresponding upper bound as

Appendix D Fast Algorithms for Projection and Rounding

In this section, we show how projection (w.r.t. Mahalanobis norm) can be done in time and rounding in time for the uniform matroid. This projection algorithm also proves to be useful in case of partition matroid polytope. We also discuss a projection method on general matroid base polytopes, based on the method of \citetactive-set-bach, which needs to solve a total number of submodular function minimization (SFM) tasks (details below).

d.1 Efficient projection on the uniform matroid

Let be a diagonal matrix with positive entries, . Our aim is to project a vector on the uniform matroid base polytope defined as

The polytope is the convex hull of all the vectors that have precisely ones and zeros. Projecting onto entails finding a point in , such that

where is the Mahalanobis norm (i.e. the Mahalanobis distance to ). Note that in the special case of , this problem boils down to orthogonal projection of onto . We first transform this problem into an orthogonal projection, and solve that projection in .

(13)

where (13) suggests an orthogonal projection on the polytope . By defining the vector , one has . Theorem 6 shows that this projection can be done in , and Algorithm 2 depicts the algorithm achieving the solution.

Theorem 6.

Let , where is given. Then for any given point one can find the solution to in time. Moreover this solution is unique.

Proof.

Let us begin by writing the KKT optimality conditions for the projected vector . The Lagrangian is defined by

where and . Minimizing the Lagrangian w.r.t. gives for each :

(14)

and also considering complementary slackness, we should have and . If one provides suitable and that satisfy the equations above, then would be the optimal solution. In what follows, we construct and provide suitable .

For each , define , where and are applied element-wise. By definition, one has . Let . We claim that if for a value of , , we are done, since , and it satisfies the KKT conditions: If , by definition of it means that , so we can set and . If , it means , so we can set and . Otherwise, , which in that case we set .

So it suffices to provide an such that . For each , define and . It’s obvious that if then , if then , and otherwise . So is a continuous decreasing function, and so will be . Note that if , then and if , then . So by continuity, there is some such that . Now let be the set of all distinct values among and . It’s clear that for all , is a linear function. By exploiting this fact, we can find by searching through these endpoints. Detailed procedure is explained in Algorithm 2. ∎

1:  Input: vectors and , s.t. .
2:  
3:  
4:  Sort elements in , so that
5:  
6:  for  do
7:      {calculate function value at the new point using current slope }{check if is between and }
8:     if  then
9:        
10:        return the projected vector as follows:
11:     end if
12:      {for these , ’s slope is changing from to }
13:      {for these , ’s slope is changing from to }
14:     
15:  end for
Algorithm 2 Projection on the Scaled Uniform Matroid Polytope

d.2 Efficient projection on Partition matroid base polytope

Let be a ground set and be a partition of . A partition matroid, includes all sets such that for all we have . It’s easy to see that the base polytope would be

In order to project onto , we first note that it becomes a separable objective, partitioned over . This means that it is sufficient to project onto the uniform matroid of , for all . Since each projection takes time, the total process would be .

d.3 Projection on general matroid base polytopes

Let us now ask whether there is an efficient projection algorithm for general matroid polytopes. Here, we argue that the method proposed by \citetactive-set-bach would be a reasonable candidate in the case of general matroid polytopes.

Let be a submodular function, such that , and let be its Lovasz extension. We define the base polytope of as the set

It can be shown [Bach et al.(2013)] that the Lovasz extension is the support function of this polytope, i.e.

(15)

For any consider the task of minimizing the following objective with respect to :

(16)

By using (15), we can rewrite (16) in the following dual form

(17)

in which the latter expression is precisely the projection of on . In \citetactive-set-bach, the authors have exploited the structural properties of the Lovasz extension and the faces of the base polytope to create the so-called “Active-set” algorithm. The Active-set algorithm iteratively solves instances of isotonic regression as well as submodular function minimization tasks, whose overall complexity is less than a single submodular function minimization call‌ (recall that by submodular function minimization, we mean the task of solving ). By knowing (17), the algorithm can be viewed as a sequence of iterative projections on outer-approximations of the base polytope.

For any matroid, its associated rank function is a monotone submodular function. Also, the base polytope for a matroid’s rank function is exactly the matroid base polytope. As a result of (17), we can use the Active-set algorithm to perform projections on the matroid base polytope. Interestingly, in the case of uniform matroids, the main parts of our projection scheme has similar counterparts as in the Active-set scheme. However, runtime complexity is significantly different due to several differences such as optimality checks: In our approach, this check is done in , but in Active-set scheme, in each iteration, one should solve approximately submodular minimization tasks. However, the Active-set approach is more general, as explained above.

d.4 The Randomized-Pipage-Rounding procedure

The randomized pipage rounding procedure was first proposed in [Calinescu et al.(2011)Calinescu, Chekuri, Pál, and Vondrák] for any matroid . Here, we show how this procedure can be efficiently done (in linear time) for the uniform matroid. Suppose we have a matroid and a point in its corresponding base polytope. We want to round to a vertex of the base polytope. In each step of the algorithm, one has a fractional solution and a tight set containing at least two fractional variables (recall that if the matroid rank function is , a set is tight if ; Tight sets are exactly those constraints in the base polytope who are tight at ). It modifies two fractional variables in such a way that their sum remains constant, until some variable becomes integral or a new constraint becomes tight. Note that since the sum of all of elements of is an integer (rank of the matroid), there exist at least two fractional variables in the case that the point is fractional.

For our purpose, we are faced with uniform matroid, which we argue that finding tight constraints is easy, i.e., we can compute the HitConstraint subroutine in a very fast way. This subroutine is given a fractional point and two variables and , and tries to increase and decrease simultaneously, and find a new tight constraint . For sure, one should search for this new tight set through the sets having inside them but not . So let denote the family of all subsets containing and not containing . So we are interested in , the maximum increase in (and decrease in ) that does not violate any polytope condition, but produces a new tight constraint. We claim that is trivial in case of the uniform matroid: . Also the new tight set is either or .

This simple form of the HitConstraint gives an efficient algorithm for Randomized-Pipage-Rounding, which we describe in Algorithm 3. Moreover, one has the following Theorem:

Theorem 7.

Let be the uniform matroid and be a fractional point inside . Then Randomized-Pipage-Rounding returns an integral point in time, such that

The proof of this algorithm’s correctness is similar to the original one given in [Calinescu et al.(2011)Calinescu, Chekuri, Pál, and Vondrák]. It is also noteworthy that our algorithm runs in time compared to , as described in [Calinescu et al.(2011)Calinescu, Chekuri, Pál, and Vondrák].

1:  Input: fractional ; defining the matroid rank
2:  while  fractional do
3:     Select and among fractional variables
4:     if   then
5:        Let
6:        With probability , set and , and with probability , set and .
7:     else
8:        Let
9:        With probability , set