Year: 2021

Abstract:
Cardinality constrained combinatorial optimization requires selecting an optimal subset of k elements, and it will be appealing to design data-driven algorithms that perform TopK selection over a probability distribution predicted by a neural network. However, the existing differentiable TopK operator suffers from an unbounded gap between the soft prediction and the discrete solution, leading to inaccurate estimation of the combinatorial objective score. In this paper, we present a self-supervised learning pipeline for cardinality constrained combinatorial optimization, which incorporates with Gumbel-Sinkhorn-TopK (GS-TopK) for neardiscrete TopK predictions and the re-parameterization trick resolving the nondifferentiable challenge. Theoretically, we characterize a bounded gap between the Maximum-A-Posteriori (MAP) inference and our proposed method, resolving the divergence issue in the previous differentiable TopK operator and also providing a more accurate estimation of the objective score given a provable tightened bound to the discrete decision variables. Experiments on max covering and discrete clustering problems show that under comparative time budget, our method outperforms state-of-the-art Gurobi solver and the novel one-shot learning method Erdos Goes Neural.

heading: 1 INTRODUCTION
text: In this paper, we aim to solve a family of cardinality constrained combinatorial optimization: given the cardinality budget k, our aim is to find the optimal set with size k which minimizes the objective function f(x). Denote ||x||0 ≤ k as the cardinality constraint, we have:
min x f(x) s.t. ||x||0 ≤ k (1)
Extensive applications of these problems can be found in machine learning and operations research, for example, discovering coresets for dataset summarization and continual learning (Borsos et al., 2020), discovering the most influential seed users in social networks (Chen et al., 2021), and planning facility locations in operation research (Liu, 2009). Although there has been a success in discovering the important submodularity (Fujishige, 1991) for a family of these problems with a bounded approximation ratio of greedy algorithms, the ubiquitous NP-hard challenge is always motivating us to develop improved algorithms with better efficiency-accuracy trade-off. The recent success of machine learning has inspired researchers to develop learning-based algorithms for combinatorial optimizations (Vinyals et al., 2015; Dai et al., 2016; Karalias & Loukas, 2020; Wang et al., 2021b), and this paper falls in line with these works.
According to the taxonomy by a recent survey (Peng et al., 2021), existing machine learning methods for combinatorial optimization can be mainly categorized into two groups: multi-step methods (Vinyals et al., 2015; Dai et al., 2016; Lu et al., 2019; Wang et al., 2021a) and one-shot methods (Wang et al., 2019; Karalias & Loukas, 2020; Li et al., 2019). The former line of works usually formulates the combinatorial optimization solving procedure as a Markov Decision Process (MDP) and refers to novel reinforcement learning methods (Mnih et al., 2013; Schulman et al., 2017) to tackle the MDP. Perhaps the largest drawback of such a multi-step pipeline is that the training and
inference procedures can be inefficient. Based on some previous efforts (Wang et al., 2019; Li et al., 2019), there is a recent seminal work (Karalias & Loukas, 2020) aiming to develop a one-shot pipeline for general combinatorial problems. Although being a fast and general pipeline, Karalias & Loukas (2020) do not encode constraint in the model prediction and thus the estimation of the combinatorial objective score may be inaccurate.
In this paper, we aim to facilitate the one-shot learning pipeline with cardinality constraints, such that the constraints can be softly encoded by the output of a neural network. Regarding the input to TopK as a probability distribution, the maximum-a-posterior (MAP) inference is done by sorting all probabilities and selecting the K-largest ones, however, this process is non-differentiable. Though appealing, it is non-trivial to encode the constraint in a differentiable manner with bounded gap w.r.t. the MAP inference and the bound of the existing SOFT-TopK algorithm (Xie et al., 2020) diverges if the k-th and (k+1)-th probabilities are equal. Unfortunately, such a dilemma is possible because the probabilities are outputs from a neural network. Besides, concerning the discrete decision variables in combinatorial optimization, it is also non-trivial to build an accurate estimation of the objective score with neural network output (which has to be continuous in order to preserve the gradient).
To tackle the non-differentiable challenge and the ill-estimated objective score, we refer to the seminal reparameterization technique with Gumbel noise, which can mimic the sampling procedure from the continuous probability distributions to near-discrete decision variables, while still preserving the gradient. Such a technique has been found successful for Top-1 (Jang et al., 2017), permutation (Mena et al., 2018), and sorting (Grover et al., 2019) distributions, and in this paper, we improve the SOFT-TopK operator (Xie et al., 2020) by introducing the Gumbel re-parameterization trick, and we theoretically prove that its gap can be bounded with the existence of the Gumbel noise, in contrast to the gap of SOFT-TopK (Xie et al., 2020) that may turn diverged.
To this end, we present Gumbel-Sinkhorn-TopK (GS-TopK) to facilitate the one-shot selfsupervised learning pipeline (Karalias & Loukas, 2020) with constrained neural network output. The proposed GS-TopK operator improves the SOFT-TopK (Xie et al., 2020), where we theoretically characterize a bounded gap by introducing non-zero noise factors. Besides, we show how to integrate the proposed GS-TopK method into the combinatorial optimization pipeline, with an accurate estimation of the objective score, formulating a self-supervised learning pipeline following (Karalias & Loukas, 2020). The contributions of this paper are summarized as follows:
1) We propose a self-supervised learning approach to solve cardinality constrained combinatorial optimizations in one shot. The cardinality constraints are encoded in the model output via our proposed GS-TopK, and our pipeline is differentiable via the re-parameterization trick. 2) We theoretically characterize a bounded gap between the MAP inference and the approximate solution by Gumbel-Sinkhorn-TopK without requiring the large gap of k-th and (k+ 1)-th probabilities, which improves the result of the state-of-the-art SOFT-TopK algorithm (Xie et al., 2020). This property is crucial for an accurate estimation of the combinatorial objective score in learning. 3) Experiments on max covering and discrete clustering problems show that under a comparative or smaller inference time budget, our method surpasses the state-of-the-art one-shot learning methods e.g. (Karalias & Loukas, 2020), and also surpasses the state-of-the-art Gurobi and SCIP solvers.

heading: 2 RELATED WORK

heading: 2.1 LEARNING OF COMBINATORIAL OPTIMIZATION
text: Recently it is a trending topic of developing machine learning based combinatorial optimization methods to achieve a better trade-off between time and accuracy by modern data-driven approaches. Here we follow the taxonomy by Peng et al. (2021) and discuss some representative works.
Multi-step methods predict the solution by a multi-step decision process. Since the fundamental work (Vinyals et al., 2015) that proposes an encoding-decoding architecture, researchers have explored this general pipeline to various real-world combinatorial optimization problems. The most popular framework follows the reinforcement learning pipeline introduced by Khalil et al. (2017), where the state is encoded by a graph neural network, and multi-step actions are taken to sequentially construct the solution. Applications have been explored on job scheduling (Zhang et al., 2020), computer resource allocation (Mao et al., 2019), quadratic assignment programming (Liu et al., 2020), bin-packing (Duan et al., 2019), just to name a few. The advantage of multi-step methods is that they can handle combinatorial constraints that are non-trivial to be encoded in a differentiable manner, by
adding constraints on the feasible actions in the multi-step roll-out procedure. However, the training of reinforcement learning may be inefficient and vulnerable, and there are also some recent efforts on developing one-shot combinatorial optimization learning methods.
One-shot methods predict the solution by a single forward pass. Despite the challenges in encoding general constraints, one-shot learning methods have been successful on certain problems where the constraint can be encoded by the neural network. For example, in graph matching (Wang et al., 2019; 2021b) the constraint can be encoded by Sinkhorn algorithm (Sinkhorn & Rangarajan, 1964), and also for graph similarity regression problems (Bai et al., 2019; Li et al., 2019) that are unconstrained. The seminal work (Karalias & Loukas, 2020) aims to develop a general pipeline for one-shot learning of combinatorial optimization, by encoding the violation of the constraint as part of its objective score. However, our experimental results show that compared to Karalias & Loukas (2020), it will be more appealing if we can explicitly encode the constraint in the neural network output.

heading: 2.2 OPTIMAL TRANSPORT AND GUMBEL RE-PARAMETERIZATION
text: Optimal Transport. The optimal transport (OT) is a linear programming problem that aims to find the optimal transportation plan with minimal cost. Though the development of the Sinkhorn algorithm can be date back to 1960s (Sinkhorn & Rangarajan, 1964), it has drawn increasing attention from the machine learning community since the pioneering work (Cuturi, 2013) which introduces the entropic-regularization for the optimal transport problem. Since the Sinkhorn algorithm (Sinkhorn & Rangarajan, 1964) incorporates row- and column-wise normalization which is differentiable, OT can be incorporated within any end-to-end deep learning pipelines. Thus, OT has been found successful for deep generative models (Arjovsky et al., 2017; Patrini et al., 2019), solving jigsaw puzzles (Santa Cruz et al., 2018), and deep graph matching (Fey et al., 2020; Yu et al., 2020). In complement to the recent line of differentiable sorting papers (Petersen et al., 2021; Cuturi et al., 2019), Xie et al. (2020) develop the differentiable TopK algorithm namely SOFT-TopK, by formulating the TopK problem as an OT problem. However, the theoretical gap between the prediction of SOFT-TopK and the result obtained by sorting diverges if the k-th and (k + 1)-th largest items are equal. In this paper, we prove that this gap can be bounded by introducing Gumbel noise.
Gumbel Re-parameterization. The Gumbel distribution with the re-parameterization trick is shown effective to mimic the sampling procedure from a given distribution in a differentiable manner. Jang et al. (2017); Maddison et al. (2017) take the initiative to develop Gumbel-Softmax, where one-hot vectors are sampled from a softmax probability distribution. The non-differentiable sampling procedure turned to be differentiable by the namely re-parameterization trick (Kingma & Welling, 2013). Inspired by Jang et al. (2017), Mena et al. (2018) propose a Gumbel-Sinkhorn method to learn the latent permutations, which is relevant to the GS-TopK method developed in this paper. Besides, Paulus et al. (2020) develop a general framework for generating combinatorial distributions via random distributions. However, these existing methods do not study the application in combinatorial optimization learning, and the gap between the MAP inference and the result after Gumbel re-parameterization is not well characterized in Mena et al. (2018). In this paper, we explore the application of the re-parameterization trick with Gumbel noise in combinatorial optimization, also we theoretically characterize the gap between the MAP inference and the proposed Gumbel method.

heading: 3 THE PROPOSED GUMBEL-SINKHORN-TOPK

heading: 3.1 PRELIMINARY: SOFT-TOPK ALGORITHM
text: The TopK problems are ubiquitous in machine learning, yet the MAP inference approach (performing sorting at the time cost of O(m logm) and selecting the first k items) is unfriendly for current gradient-based deep neural networks. To introduce differentiable TopK for machine learning, Xie et al. (2020) formulate the TopK problem as an optimal transport (OT) problem. Specifically, for selecting the k largest items from the probability vector s of size m (k < m), the OT problem refers to moving k items to one destination (selected), and the other m− k elements to another destination (not selected). Thus, the marginal distributions of the OT are:
c = [1 1 ... 1]︸ ︷︷ ︸ m items
r = [ m− k k ]} 2 destinations (2)
Algorithm 1: Gumbel-Sinkhorn-TopK (GS-TopK) for Solving Cardinality Constrained Combinatorial Optimization Input: List s with m items; TopK size k; Sinkhorn factor τ ; noise factor σ; sample size #G. 1 for a ∈ {1, 2, ..., #G} do 2 for all si, s̃i = si − σ log(− log(ui)), where ui is sampled from (0, 1) uniform distribution;
3 D̃ = [ s̃1 −min(s) ... s̃m −min(s) max(s)− s̃1 ... max(s)− s̃m ] ; c = [1 1 ... 1]︸ ︷︷ ︸
m items
; r = [ m− k k ] ;
4 T̃a = exp(−D̃/τ); 5 while not converged do 6 D̃r = diag(T̃a1 r); T̃a = D̃−1r T̃a; 7 D̃c = diag(T̃ > a 1 c); T̃a = T̃aD̃−1c ;
Output: A list of transport matrices [T̃1, T̃2, ..., T̃#G].
We can design the destinations to be the min/max values of s, such that the TopK items are moved to max(s), and the other items are moved to min(s). Then, the distance matrix of OT is given as:
D = [ s1 −min(s) s2 −min(s) ... sm −min(s) max(s)− s1 max(s)− s2 ... max(s)− sm ] . (3)
While the OT is formulated as an integer linear programming problem:
min T
tr(T>D) s.t. T ∈ {0, 1}m×2,T1 = r,T>1 = c, (4)
where T is the transportation plan which also corresponds to the TopK solution of the problem, and 1 is a column vector whose all elements are 1s. The optimal solution T∗ should be equivalent to the solution by firstly sorting all items and then selecting the TopK items. It is also equivalent to the MAP inference by regarding s1, ..., sm as probabilities. Following the differentiable Sinkhorn algorithm for OT (Cuturi, 2013; Sinkhorn & Rangarajan, 1964), the binary constraint on T is relaxed to continuous values [0, 1], and Eq. (4) is modified with entropic regularization:
min Tτ
tr(Tτ>D) + τh(Tτ ) s.t. Tτ ∈ [0, 1]m×2,Tτ1 = r,Tτ>1 = c, (5)
where h(Tτ ) = ∑ i,j T τ ij log T τ ij is the entropic regularizer (Cuturi, 2013). Given any real-valued matrix D, Eq. (5) is solved by firstly enforcing the regularization factor τ : Tτ = exp(−D/τ). Then Tτ is row- and column-wise normalized alternatively:
Dr = diag(T τ1 r), Tτ = D−1r Tτ , Dc = diag(Tτ>1 c), Tτ = TτD−1c , (6)
where means element-wise division. We denote Tτ∗ as the converged solution, which is the optimal solution to Eq. (5).
Theorem 2 of (Xie et al., 2020). Denote xk, xk+1 as the k-th (k + 1)-th largest items, we have
||T∗ −Tτ∗||F ≤ 2mτ log 2
|xk − xk+1| , (7)
where the gap is controlled by τ if |xk − xk+1| > 0. However, the above bound diverges if xk = xk+1, which is unavoidable because xk, xk+1 are outputs by a neural network. It motivates us to develop our improved version Gumbel-Sinkhorn-TopK by introducing the Gumbel noise.

heading: 3.2 GUMBEL-SINKHORN-TOPK (GS-TOPK) ALGORITHM
text: In this section, we present our proposed Gumbel-Sinkhorn-TopK (GS-TopK) as summarized in Alg. 1 and we theoretically characterize the gap between the MAP inference and the prediction by GS-TopK. As the re-parameterization trick (Jang et al., 2017), instead of sampling from a distribution that is non-differentiable, we add random variables to probabilities predicted by neural networks. The Gumbel distribution is characterized as:
gσ(u) = −σ log(− log(u)), (8)
where σ controls the variance and u is from (0, 1) uniform distribution. We can update s and D as:
s̃i = si + gσ(ui) (9)
D̃ = [ s̃1 −min(s) s̃2 −min(s) ... s̃m −min(s) max(s)− s̃1 max(s)− s̃2 ... max(s)− s̃m ] . (10)
Again we formulate the integer linear programming version of the OT with Gumbel noise:
min Tσ
tr(Tσ>D̃) s.t. Tσ ∈ {0, 1}m×2,Tσ1 = r,Tσ>1 = c, (11)
where the optimal solution to Eq. (11) is denoted as Tσ∗. Without loss of generality, in the following we sort s1, s2, s3, ..., sm in descending order, into X = x1, x2, x3, ..., xm, and xk, xk+1 are the k-th and (k + 1)-th largest items, respectively. We characterize the gap between T∗ and Tσ∗ under the expectation over Gumbel distribution.
Now we present the following two lemmas to help establish our theoretical results, in terms of a bounded gap between the MAP inference and the approximate solution by Gumbel-Sinkhorn-TopK, Lemma 1 (Bias by the Gumbel Noise). For the integer linear programming solutions with and without Gumbel noises, we have
Eu [||T∗ −Tσ∗||F ] ≤ √
4m
1 + exp |xk−xk+1|σ
(12)
Proof sketch: The gap between T∗,Tσ∗ is derived by counting the expected number of times that, after adding the Gumbel noise, an element may become larger than xk or smaller than xk+1. For the top k elements, the probability that it will become smaller than xk+1 is upper bounded by 1/(1 + exp |xk−xk+1|σ ). For the last (m−k) elements, the probability that it will become larger than xk is also upper bounded by 1/(1 + exp |xk−xk+1| σ ). Since ||T
∗ −Tσ∗||2F changes at most by 4 if an item becomes larger/smaller than xk/xk+1, and we have m items, then we can proof the upper bound of this gap. The detailed proof is referred to Appendix A.1.
To make the integer linear programming problem feasible for gradient-based deep learning methods, we also relax the integer constraint and add the entropic regularization term:
min T̃
tr(T̃>D̃) + h(T̃) s.t. T̃ ∈ [0, 1]m×2, T̃1 = r, T̃>1 = c, (13)
which is solved by Sinkhorn algorithm on D̃: firstly T̃ = exp(−D̃/τ), then
D̃r = diag(T̃1 r), T̃ = D̃−1r T̃, D̃c = diag(T̃>1 c), T̃ = T̃D̃−1c . (14)
Here we denote the optimal solution to Eq. (13) as T̃∗. We theoretically characterize the gap between Tσ∗ and T̃∗, which corresponds to the optimal solutions of Eq. (11) and Eq. (13), respectively. Lemma 2 (Gap to a Discrete Solution). Under probability (1− ), we have
Eu [ ||Tσ∗ − T̃∗||F ] ≤ (log 2)mτ ∑ i 6=j Ω(xi, xj , σ, ), (15)
where
Ω(xi, xj , σ, ) = 2σ log
( σ − |xi−xj |+2σlog(1− ) ) + |xi − xj | ( π 2 + arctan xi−xj 2σ ) (1− )((xi − xj)2 + 4σ2)(1 + exp xi−xkσ )(1 + exp xk+1−xj σ ) , (16)
Proof sketch: The gap between Tσ∗ and T̃∗ is a generalization from the Theorem 2 of Xie et al. (2020). We denote xπk , xπk+1 as the k-th and (k + 1)-th largest items after disturbed by the Gumbel noise, and our aim becomes to prove the upper bound of Eu [ 1/(|xπk + gσ(uπk)− xπk+1 − gσ(uπk+1)|) ] , where the probability density function of gσ(uπk) − gσ(uπk+1) can be bounded by f(y) = 1/(y2 + 4). Thus we can compute the upper bound under probability (1− ) by integration. The detailed proof is referred to Appendix A.2.
Theorem 1. For the proposed Gumbel TopK solver, denote T∗ as the solution by deterministic TopK algorithm, T̃∗ as the solution of Gumbel-Sinkhorn-TopK. With probability (1− ), we have
Eu [∣∣∣∣∣∣T∗ − T̃∗∣∣∣∣∣∣
F
] ≤ √
4m
1 + exp |xk−xk+1|σ + (log 2)mτ ∑ i6=j Ω(xi, xj , σ, ), (17)
where σ is the noise factor and τ is the temperature of Sinkhorn algorithm. m is the number of candidates to be selected. Ω(xi, xj , σ, ) is defined in Eq. (16).
Proof: Based on Lemma 1 and Lemma 2, Theorem 1 is proved by triangle inequality: Eu [∣∣∣∣∣∣T∗ − T̃∗∣∣∣∣∣∣
F
] ≤Eu [||T∗ −Tσ∗||F ] + Eu [ ||Tσ∗ − T̃∗||F ]
Remarks. As derived above, the gap between T∗ and T̃ is composed of two parts: the bias by the Gumbel noise (in Lemma 1) and the gap to a discrete solution (in Lemma 2). Interestingly, the first term becomes smaller given a smaller σ, and the second term is smaller given a larger σ or a smaller τ . See the toy example in Fig. 1, a larger σ can tighten the gap of GSTopK to a discrete solution, which is welcomed in our combinatorial optimization learning scenario where T̃∗ is used to compute the objective score in a self-supervised learning pipeline.
If |xk − xk+1| > 0, with σ → 0+, Eq. (17) degenerates to the bound derived by (Xie et al., 2020) and only differs by a constant factor (see Appendix A.3 for details):
lim σ→0+
Eu [ ||T∗ − T̃∗||F ] ≤ (π log 2)mτ
(1− )|xk−xk+1|
It makes a strong assumption that |xk−xk+1| > 0, and the bound diverges if xk = xk+1. Un-
fortunately, xk, xk+1 are predictions by a neural network, whereby such an assumption may not be satisfied. In comparison, given σ > 0, our conclusion in Eq. (17) is bounded for any xk, xk+1.

heading: 3.3 APPLICATION TO SELF-SUPERVISED COMBINATORIAL OPTIMIZATION LEARNING
text: In the following, we discuss the application of GS-TopK to self-supervised combinatorial optimization learning, by taking the max covering problem as an example. We consider the max covering problem with m sets and n objects. Each set may cover any number of objects, and each object is
Algorithm 2: Gumbel-Sinkhorn-TopK Learning for Solving the Max Covering Problem Input: bipartite adjacency A; values v; learning rate α; GS-TopK parameters k, τ, σ, #G.
1 if Training then 2 Randomly initialize neural network weights θ; 3 if Inference then 4 Load a pretrained neural network weights θ; Jbest = 0; 5 while not converged do 6 s = GraphSageθ(A); [T̃1, T̃2, ..., T̃#G] = GS-TopK(s, k, τ, σ, #G); 7 for all i, J̃i = min(T̃i[2, :]A, 1) · v; J = mean([J̃1, J̃2, ..., J̃#G]); 8 if Training then 9 update θ with respect to the gradient ∂J∂θ and learning rate α by gradient ascent;
10 if Inference then 11 update s with respect to the gradient ∂J∂s and learning rate α by gradient ascent; 12 for all i, J̃i = TopK(T̃i[2, :])A · v; Jbest = max([J̃1, J̃2, ..., J̃#G], Jbest);
13 if Homotopy Inference then 14 Shrink the values of τ, σ and jump to line 5;
Output: Learned network weights θ (if training)/The best objective Jbest (if inference).
associated with a value. We aim to find k sets (where k < m) such that the covered objects have the maximum sum of values. The general illustration of our pipeline can be found in Fig. 2.
Firstly, we build a bipartite graph whose two disjoint sets of vertices are the sets and the objects. An edge is defined if an object is covered by a set. Denote v ∈ Rn as the value of each object, and A ∈ {0, 1}m×n as the adjacency matrix, the problem is formulated as
max x n∑ j=1
( I ( m∑ i=1 xiAij ) · vj ) s.t. x ∈ {0, 1}m, ||x||0 ≤ k (18)
where I(x) is an indicator I(x)i = 1 if xi ≥ 1 else I(x)i = 0. To encode the bipartite graph, we exploit three layers of GraphSage (Hamilton et al., 2017) followed by a fully-connected layer with sigmoid to predict the probability of selecting each set, which is denoted as s ∈ [0, 1]m. The probabilities s are fed into the GS-TopK algorithm, which outputs a batch of near-discrete transportation matrices [T̃1, T̃2, ..., T̃#G], whose second row is regarded as a continuous approximation of the decision variables. The objective value is estimated as:
J̃i = min(T̃i[2, :]A, 1) · v, J = mean([J̃1, J̃2, ..., J̃#G]) (19)
For the discrete clustering problem, we aim to select k objects from the full set of m objects, by minimizing the sum of distance for each object between itself and its nearest selected object. It differs from the k-means clustering by adding the “discrete” constraint that cluster centers have to be a subset of objects, which is also known as facility location problem in operations research (Liu, 2009). Denote ∆ ∈ R+m×m as the distance matrix for each pair of objects, we formulate the discrete clustering problem as follows (please refer to Appendix B for details):
min x m∑ j=1 min({∆i,j | ∀xi = 1}) s.t. x ∈ {0, 1}m, ||x||0 ≤ k (20)
As illustrated in Alg. 2, during training, the Adam optimizer (Kingma & Ba, 2014) is applied for learning the neural network weights. During inference, the neural network prediction is regarded as initialization, and we again optimize the distribution w.r.t. the combinatorial objective score by Adam, during which a search procedure is also performed among all Gumbel samples.
Homotopy GS-TopK. Based on our Theorem 1, the gap between T∗ and T̃∗ can be tightened by shrinking τ and σ, which motivates us to develop a homotopy variant (Xiao & Zhang, 2012; Xu et al., 2016) of GS-TopK with gradually decreased τ, σ, such that T̃∗ and J̃ can provide more accurate estimation during inference. See line 13 in Alg. 2 for details.

heading: 4 EXPERIMENTS AND DISCUSSION

heading: 4.1 PROTOCOLS AND BASELINES
text: We apply the proposed GS-TopK algorithm to learning two representative problems of cardinality constrained NP-hard combinatorial optimization: i) max covering and ii) discrete clustering problem. For each problem, both a smaller-sized case and a larger-sized case are considered, and we build separate training/testing datasets for learning methods. In this paper, we referred to randomly generated datasets following previous machine learning for combinatorial optimization papers (Khalil et al., 2017; Kool et al., 2019) concerning the limited large-scale open-source dataset. For the max covering problem, our data distribution follows the distribution in ORLIB1. For the discrete clustering problem, the data points are uniformly generated in a 2D plane.
We compare the following representative and diverse baselines with our proposed method:
Classic Algorithms. Classic algorithms are well-developed strong baselines with theoretical guarantees. In our experiment, we consider the greedy algorithm as the baseline for the max covering problem, which has the worst-case approximation ratio of (1 − 1/e) due to the submodularity of the max covering problem (Fujishige, 1991). For the discrete clustering problem, we resort to the classic kmeans (Macqueen, 1967) and kmeans++ (Arthur & Vassilvitskii, 2007) with improved initialization techniques. These algorithms are fast but effective, even surpassing commercial solvers with only 1/1000 inference time on certain test cases.
General Purpose Solvers. We also compare two general-purpose solvers Gurobi 9.0 (under educational license) (Gurobi Optimization, LLC, 2021) and SCIP 7.0 (Gamrath et al., 2020). Gurobi is the state-of-the-art commercial solver for general combinatorial optimization tasks, and it is known heavily engineered by various speedup techniques. SCIP is the state-of-the-art open-source solver, which is also regarded as a baseline. We formulate the combinatorial problems as integer programming and then call the solver to solve the problems, with a limited time budget. Our implementation of the solvers is with the Google ORTools API 2.
One-shot Learning Baselines. We compare with the general Erdos Goes Neural (EGN) (Karalias & Loukas, 2020) framework that models each decision variable as a Bernoulli distribution and encodes the constraint as part of the learning objective. Since the authors of (Karalias & Loukas, 2020) do not consider cardinality constrained problems in their experiment, we re-implement and tune the EGN model based on their official implementation3. We also compare with SOFT-TopK (Xie et al., 2020), which is theoretically characterized as a special case of our proposed GS-TopK when σ = 0. Also to ensure the effectiveness of the re-parameterization trick, we facilitate SOFT-TopK with Gumbel sampling during inference. All one-shot learning methods are with the same model structure for a fair comparison and are trained and tested on separate datasets. Besides, we empirically find the learning process of one-shot methods converges within tens of minutes. Since the reinforcement learning methods (Khalil et al., 2017; Chen & Tian, 2019) require significantly more training time, we only include one-shot methods for a fair comparison.
1http://people.brunel.ac.uk/˜mastjjb/jeb/orlib/scpinfo.html 2https://developers.google.com/optimization 3https://github.com/Stalence/erdos_neu

heading: 4.2 EXPERIMENTAL RESULTS
text: Table 1 and Table 2 report results on max covering and discrete clustering problems, respectively. An interesting finding is that the classic algorithms are efficient but also very effective. However, if we can afford a higher time budget, there seems no universal approach to improve the results of the classic algorithms. In comparison, Gurobi and SCIP solvers can trade-off time for better accuracy, but they are less effective than our learning methods when considering the objective scores under the same time budget. Among one-shot learning methods, Erdos Goes Neural (Karalias & Loukas, 2020) is the only method that does not encode the constraint in the model output, and the experimental results suggest that it is appealing if we can encode the constraint in the model output. Directly applying Gumbel sampling to SOFT-TopK seems either an effective improvement, and satisfying results are achieved when exploiting the Gumbel re-parameterization trick and enabling gradientbased optimization over the combinatorial objective score. We further validate the effectiveness of Homotopy GS-TopK, where the objective scores are improved at the cost of more inference time.
Further Discussions. Our insight for the effectiveness of GS-TopK over other state-of-the-art oneshot learning methods (Karalias & Loukas, 2020; Xie et al., 2020) is that it can provide a more accurate estimation of the combinatorial objective score with discrete decision variables. Fig. 3 shows that our GS-TopK estimates a lower bound of the true objective score at the beginning 200 iterations while providing an accurate estimation afterward. Besides, the Homotopy version of our algorithm can gradually reduce the variance and help converge to a better solution.
REPRODUCIBILITY STATEMENT
The following efforts are made to ensure the reproducibility of this paper:
• We provide a complete proof about our Lemmas and Theorems in the Appendix, see Appendix A;
• We discuss the hyper-parameter configurations and the experiment setup in Sec 4.1; • We provide the implementation details on max covering problem in Sec. 3.3, and the im-
plementation details on discrete problem problem in Appendix B;
• We discuss the way of generating random dataset in Sec 4.1; • The code will be made publicly available once this paper is accepted.

heading: A PROOF OF THEOREMS
text: Before starting the detailed proof of the lemmas and theorems, firstly we recall the notations used in this paper:
• T∗ = TopK(D) is the optimal solution of the integer linear programming form of the OT problem Eq. (4), which is equivalent to the solution by firstly sorting all items and then select the TopK items;
• Tτ∗ = Sinkhorn(D) is the optimal solution of the entropic regularized form of the OT problem Eq. (5) solved by Sinkhorn algorithm. It is also the output by SOFT-TopK (Xie et al., 2020);
• Tσ∗ = TopK(D̃) is the optimal solution of the integer linear programming form of the OT problem after disturbed by the Gumbel noise Eq. (11), which is equivalent to the solution by firstly adding the Gumbel noise, then sorting all items and finally select the TopK items;
• T̃∗ = Sinkhorn(D̃) is the optimal solution of the entropic regularized form of the OT problem after disturbed by the Gumbel noise Eq. (13) solved by Sinkhorn algorithm. It is also the output by our proposed GS-TopK.
Lemma 3. Given real numbers xi, xj , and ui, uj are from i.i.d. (0, 1) uniform distribution. After Gumbel perturbation, the probability that xi + gσ(ui) > xj + gσ(uj) is:
P (xi + gσ(ui) > xj + gσ(uj)) = 1
1 + exp−xi−xjσ . (21)
Proof. Since gσ(ui) = −σ log(− log(ui)), P (xi + gσ(ui) > xj + gσ(uj)) is equivalent to the probability that the following inequality holds:
xi − σ log(− log(ui)) > xj − σ log(− log(uj)) (22)
And we have
xi − xj > σ log(− log(ui))− σ log(− log(uj)) (23) xi − xj σ > log ( log(ui)
log(uj)
) (24)
e xi−xj σ > log(ui)
log(uj) (25)
Since uj ∈ (0, 1), log(uj) < 0. Then we have
log(uj) < log(ui)e − xi−xj σ (26)
log (uj) < log
( u exp− xi−xj σ
i
) (27)
uj < u exp−
xi−xj σ
i (28)
Since ui, uj are i.i.d. uniform distributions, the probability when the above formula holds is
∫ 1 0 ∫ uexp− xi−xjσi 0 duj dui = ∫ 1 0 u exp− xi−xj σ i dui = 1 1 + exp−xi−xjσ (29)
Thus the probability that xi + gσ(ui) > xj + gσ(uj) after Gumbel perturbation is:
P (xi + gσ(ui) > xj + gσ(uj)) = 1
1 + exp−xi−xjσ (30)
A.1 PROOF OF LEMMA 1
Proof of Lemma 1. If one of the top-k items becomes smaller than xk+1 after distorted by the Gumbel noise, then ||T∗ −Tσ∗||F increases at most by 2. Without loss of generality, here we assume xi as one of the top-k items i.e. i ≤ k. According to Lemma 3, the probability that xi becomes smaller than xk+1 after perturbation is
P (xi + gσ(ui) < xk+1 + gσ(uk+1)) = 1
1 + exp−xk+1−xiσ (31)
Similarly, if one of the last-(m − k) items becomes larger than xk after distorted by the Gumbel noise, then ||T∗ −Tσ∗||F also increases at most by 2. Without loss of generality, here we assume xj as one of the last-(m− k) items i.e. j ≥ k + 1. According to Lemma 3, the probability that xj becomes larger than xk after perturbation is
P (xj + gσ(uj) > xk + gσ(uk)) = 1
1 + exp−xj−xkσ (32)
To compute the expectation of ||T∗ −Tσ∗||, we sum all probabilities and we have
Eu [ ||T∗ −Tσ∗||2F ] ≤ k∑ i=1
4
1 + exp−xk+1−xiσ + m∑ j=k+1
4
1 + exp−xj−xkσ (33)
where the constant 4 appears because ||T∗ −Tσ∗||2F will change at most by 4 if one more item crosses the boundary. Since for all i ≤ k and j ≥ k + 1, we have
xk+1 − xi ≤ xk+1 − xk (34) xj − xk ≤ xk+1 − xk (35)
Thus we have
Eq. (33) ≤ k∑ i=1
4
1 + exp−xk+1−xkσ + m∑ j=k+1
4
1 + exp−xk+1−xkσ (36)
= m∑ i=1
4
1 + exp−xk+1−xkσ (37)
= 4m
1 + exp |xk−xk+1|σ
(38)
Since Eu [||T∗ −Tσ∗||F ] ≤ √ Eu [ ||T∗ −Tσ∗||2F ] , we have
Eu [||T∗ −Tσ∗||F ] ≤ √
4m
1 + exp |xk−xk+1|σ
(39)
A.2 PROOF OF LEMMA 2
Proof of Lemma 2. By disturbing X with i.i.d. Gumbel noise, we have
X ′ = x1 + gσ(u1), x2 + gσ(u2), x3 + gσ(u3), ..., xm + gσ(um) (40)
where gσ(u) = −σ log(− log(u)) is the Gumbel noise modulated by noise factor σ, and u1, u2, u3, ..., um are i.i.d. uniform distribution. We define π as the permutation of sorting X ′ in descending order, i.e. xπ1 + gσ(uπ1), xπ2 + gσ(uπ2), xπ3 + gσ(uπ3), ..., xπm + gσ(uπm) are in descending order.
Recall Theorem 2 in (Xie et al., 2020), for x1, x2, x3, ..., xm we have
||T∗ −Tτ∗||F ≤ 2mτ log 2
|xk − xk+1| (41)
By substituting X with X ′ and taking the expected value, we have
Eu [ ||Tσ∗ − T̃∗||F ] ≤ Eu [ 2mτ log 2
|xπk + gσ(uπk)− xπk+1 − gσ(uπk+1)|
] (42)
Based on Lemma 3, the probability that πk = i, πk+1 = j is
P (πk = i, πk+1 = j) = 1
1 + exp−xi−xjσ ∑ ∀π ( k−1∏ a=1
1
1 + exp−xπa−xiσ m∏ b=k+2
1
1 + exp−xj−xπbσ ) (43)
where the first term denotes xi + gσ(ui) > xj + gσ(uj), the second term denotes all conditions that there are (k − 1) items larger than xi + gσ(ui) and the rest items are smaller than xj + gσ(uj).
In the following we derive the upper bound of Eu [
1
|xπk+gσ(uπk )−xπk+1−gσ(uπk+1 )|
] . We denote
Ai,j as ui, uj ∈ Ai,j , s.t. xi + gσ(ui)− xj − gσ(uj) > (44)
where is a sufficiently small number. Then we have
Eu [ 1∣∣xπk + gσ(uπk)− xπk+1 − gσ(uπk+1)∣∣ ]
= ∑ i6=j P (πk = i, πk+1 = j) Eui,uj∈Ai,j [
1
|xi + gσ(ui)− xj − gσ(uj)|
] (45)
= ∑ i6=j
( 1
1 + exp−xi−xjσ ∑ ∀π ( k−1∏ a=1
1
1 + exp−xπa−xiσ m∏ b=k+2
1
1 + exp−xj−xπbσ
)
Eui,uj∈Ai,j [
1
|xi + gσ(ui)− xj − gσ(uj)|
]) (46)
= ∑ i6=j
( 1
1 + exp−xi−xjσ ∑ ∀π ( k−1∏ a=1
1
1 + exp−xπa−xiσ m∏ b=k+2
1
1 + exp−xj−xπbσ
)
Eui,uj∈Ai,j [
1
|xi − σ log(− log(ui))− xj + σ log(− log(uj))|
]) (47)
= ∑ i6=j
( f(xi − xj , σ, z)
∑ ∀π ( k−1∏ a=1
1
1 + exp−xπa−xiσ m∏ b=k+2
1
1 + exp−xj−xπbσ
)) (48)
We denote f (δ, σ, z) as:
f (δ, σ, z) = 1
1 + exp− δσ Eui,uj
[ 1
|δ − σ log(− log(ui)) + σ log(− log(uj))|
] (49)
s.t. δ − σ log(− log(ui)) + σ log(− log(uj)) > z > 0 (50)
For the probability terms in Eq. (48), for all permutations π, there must exist πa, πb, such that
1 1 + exp−xπa−xiσ ≤ 1 1 + exp−xk−xiσ (51)
1 1 + exp−xj−xπbσ ≤ 1 1 + exp−xj−xk+1σ (52)
Thus we have
Eq. (48) ≤ ∑ i6=j
( f(xi − xj , σ, z)
1
1 + exp−xk−xiσ
1
1 + exp−xj−xk+1σ
) (53)
≤ ∑ i6=j f(xi − xj , σ, z) (1 + exp xi−xkσ )(1 + exp xk+1−xj σ )
(54)
By transforming Eq. (21) in Lemma 3, and substituting xj − xi by y, we have
Eq. (21)⇒P (gσ(ui)− gσ(uj) > xj − xi) = 1
1 + exp−xi−xjσ (55)
⇒P (gσ(ui)− gσ(uj) > y) = 1
1 + exp yσ (56)
⇒P (gσ(ui)− gσ(uj) < y) = 1− 1
1 + exp yσ =
1
1 + exp− yσ (57)
where the right side is the form of the cumulative distribution function (CDF) of standard Logistic distribution by setting σ = 1:
CDF(y) = 1
1 + exp (−y) (58)
Thus− log(− log(ui)) + log(− log(uj)) is equivalent to the Logistic distribution whose probability density function (PDF) is
PDF(y) = dCDF(y)
dy =
1
exp (−y) + exp y + 2 (59)
and in this proof we exploit an upper bound of PDF(y):
PDF(y) = 1 exp (−y) + exp y + 2 ≤ 1 y2 + 4
(60)
Based on the Logistic distribution, we can replace −σ log(− log(ui)) + σ log(− log(uj)) by σy where y is from the Logistic distribution. Thus we can derive the upper bound of f(δ, σ, z) as
follows
f (δ, σ, z) = 1
1 + exp− δσ ·
∫∞ −δ/σ+z 1 δ+σy
1 exp (−y)+exp y+2dy∫∞
−δ/σ+z 1 exp (−y)+exp y+2dy (61)
= 1
1 + exp− δσ ·
∫∞ −δ/σ+z 1 δ+σy 1 exp (−y)+exp y+2dy
1− 11+exp (δ/σ−z) (62)
= 1
1 + exp− δσ ·
∫∞ −δ/σ+z 1 δ+σy 1 exp (−y)+exp y+2dy
exp (δ/σ−z) 1+exp (δ/σ−z)
(63)
= 1
1 + exp− δσ ·
∫∞ −δ/σ+z 1 δ+σy 1 exp (−y)+exp y+2dy
1 1+exp (−δ/σ+z)
(64)
= 1 + exp (− δσ + z)
1 + exp− δσ
∫ ∞ −δ/σ+z
1
δ + σy
1
exp (−y) + exp y + 2 dy (65)
≤ 1 + exp (− δσ + z)
1 + exp− δσ
∫ ∞ −δ/σ+z
1
δ + σy
1
y2 + 4 dy (66)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
2σ log ( (zσ − δ)2 + 4σ2 ) − 2δ arctan ( z−δ/σ
2
) − 4σ log z + πδ
4δ2 + 16σ2
(67)
≤ 1 + exp (− δσ + z)
1 + exp− δσ ·
2σ log ( (zσ + |δ|)2 + 4σ2 ) − 2δ arctan ( z−δ/σ
2
) − 4σ log z + πδ
4δ2 + 16σ2
(68)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
2σ log ( (zσ + |δ|)2 + 4σ2 ) − 2δ arctan ( z−δ/σ
2
) − 2σ log z2 + πδ
4δ2 + 16σ2
(69)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
2σ log (
(zσ+|δ|)2+4σ2 z2
) − 2δ arctan ( z−δ/σ
2
) + πδ
4δ2 + 16σ2 (70)
≤ 1 + exp (− δσ + z)
1 + exp− δσ ·
2σ log (
(zσ+|δ|+2σ)2 z2
) − 2δ arctan ( z−δ/σ
2
) + πδ
4δ2 + 16σ2 (71)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
4σ log ( zσ+|δ|+2σ
z
) − 2δ arctan ( z−δ/σ
2
) + πδ
4δ2 + 16σ2 (72)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
4σ log ( zσ+|δ|+2σ
z
) + δ ( π − 2 arctan ( z−δ/σ
2 )) 4δ2 + 16σ2
(73)
≤ 1 + exp (− δσ + z)
1 + exp− δσ ·
4σ log ( zσ+|δ|+2σ
z
) + |δ| ( π − 2 arctan ( z−δ/σ
2 )) 4δ2 + 16σ2
(74)
≤ 1 + exp (− δσ + z)
1 + exp− δσ ·
4σ log ( zσ+|δ|+2σ
z
) + |δ| ( π − 2 arctan ( − δ2σ )) 4δ2 + 16σ2
(75)
= 1 + exp (− δσ + z)
1 + exp− δσ ·
4σ log ( zσ+|δ|+2σ
z
) + |δ| ( π + 2 arctan ( δ 2σ )) 4δ2 + 16σ2
(76)
where Eq. (66) is because 1exp (−y)+exp y+2 ≤ 1 y2+4 , and Eq. (74) is because π−2 arctan( z−δ/σ
2 ) ≥ 0. With probability (1− ), we have
z = log 1 + exp δσ
1− ≥ − log (1− ) (77)
1 + exp (− δσ + z) 1 + exp− δσ = 1 1− (78)
Thus
f(δ, σ, z) ≤ Eq. (76) = 1 1−
4σ log ( zσ+|δ|+2σ
z
) + |δ| ( π + 2 arctan ( δ 2σ )) 4δ2 + 16σ2
(79)
≤ 1 1−
4σ log ( σ − |δ|+2σlog(1− ) ) + |δ| ( π + 2 arctan ( δ 2σ )) 4δ2 + 16σ2
(80)
Thus we have
Eq. (54) ≤ ∑ i6=j
 4σ log ( σ − |xi−xj |+2σlog(1− ) ) + |xi − xj | ( π + 2 arctan ( xi−xj 2σ )) (1− )(4(xi − xj)2 + 16σ2)(1 + exp xi−xkσ )(1 + exp xk+1−xj σ )  (81) In conclusion, with the probability (1− ), we have Eu [ ||Tσ∗ − T̃∗||F
] ≤ ∑ i6=j
(2 log 2)mτ ( 4σ log ( σ − |xi−xj |+2σlog(1− ) ) + |xi − xj | ( π + 2 arctan xi−xj 2σ )) (1− )(4(xi − xj)2 + 16σ2)(1 + exp xi−xkσ )(1 + exp xk+1−xj σ )
(82)
= ∑ i6=j
(log 2)mτ ( 2σ log ( σ − |xi−xj |+2σlog(1− ) ) + |xi − xj | ( π 2 + arctan xi−xj 2σ )) (1− )((xi − xj)2 + 4σ2)(1 + exp xi−xkσ )(1 + exp xk+1−xj σ )
(83) =(log 2)mτ ∑ i 6=j Ω(xi, xj , σ, ) (84)
And we denote Ω(xi, xj , σ, ) as
Ω(xi, xj , σ, ) = 2σ log
( σ − |xi−xj |+2σlog(1− ) ) + |xi − xj | ( π 2 + arctan xi−xj 2σ ) (1− )((xi − xj)2 + 4σ2)(1 + exp xi−xkσ )(1 + exp xk+1−xj σ )
(85)
Finally, Theorem 1 is proved by triangle inequality by jointly considering Lemma 1 and Lemma 2.
A.3 REMARKS W.R.T. SOFT-TOPK
In the following, we add some remarks about the relationship between our conclusion of GS-TopK and the conclusion derived by the authors of SOFT-TopK (Xie et al., 2020): the SOFT-TopK is a special case of our proposed algorithm when we set σ = 0. We have the following conclusion: Theorem 2. Assume the values of xk, xk+1 are unique4, under probability (1− ), we have
lim σ→0+
Eu [ ||T∗ − T̃∗||F ] ≤ (π log 2)mτ
(1− )|xk − xk+1| (86)
4For the ease of a compact proof, we make this assumption that the values of xk, xk+1 are unique. If there are duplicate values of xk, xk+1, the bound only differs by a constant multipler therefore does not affect our conclusion: SOFT-TopK (Xie et al., 2020) is a special case of our proposed approach when σ = 0.
which differs from the conclusion in (Xie et al., 2020) by only a constant factor.
Proof. Since σ → 0+, the first term becomes 0. For the second term, we discuss four conditions as shown in Fig. 4, except for the following condition: xi = xk, xj = xk+1.
Condition 1. If xi ≥ xk, xj ≤ xk+1 (equalities do not hold at the same time), we have at least xi − xk > 0 or xk+1 − xj > 0. Then we have
lim σ→0+
1
(1 + exp xi−xkσ )(1 + exp xk+1−xj σ ) = 0 (87)
⇒ lim σ→0+ Ω(xi, xj , σ, ) = 0 (88)
Condition 2. For any case that xi < xj , we have xi − xj < 0, thus
lim σ→0+ arctan xi − xj σ = −π 2
(89)
⇒ lim σ→0+
π 2 + arctan xi − xj σ = 0 (90)
⇒ lim σ→0+ Ω(xi, xj , σ, ) = 0 (91)
Condition 3. If xi ≥ xj ≥ xk (equalities do not hold at the same time), we have xi−xk > 0. Then we have
lim σ→0+
1
1 + exp xi−xkσ = 0 (92)
⇒ lim σ→0+ Ω(xi, xj , σ, ) = 0 (93)
Condition 4. If xk+1 ≥ xi ≥ xj (equalities do not hold at the same time), we have xk+1 − xj > 0. Then we have
lim σ→0+
1
1 + exp xk+1−xj
σ
= 0 (94)
⇒ lim σ→0+ Ω(xi, xj , σ, ) = 0 (95)
Therefore, if xi 6= xk and xj 6= xk+1, the second term Ω(xi, xj , σ, ) degenerates to 0 when σ → 0+. Thus we have the following conclusion by only considering xi = xk, xj = xk+1:
lim σ→0+
Eu [ ||T∗ − T̃∗||F ] ≤
(log 2)mτ ( |xk − xk+1| ( π 2 + arctan xk−xk+1 2σ )) (1− )(xk − xk+1)2
(96)
≤ (π log 2)mτ (1− )|xk − xk+1|
(97)
B IMPLEMENTATION DETAILS ON DISCRETE CLUSTERING PROBLEM
For the discrete clustering problem, recall that we aim to select a subset of objects as cluster centers. It differs from the k-means clustering problem whose cluster centers may not have to be a subset of the objects, and the constraint for discrete clustering is more welcomed for real-world facility location planning applications. Denote ∆ ∈ R+m×m as the distance matrix for each pair of objects. Following our implementation of the max covering problem, we can also derive a GPU-friendly formulation of the discrete clustering problem using PyTorch syntax:
min x
sum(min(∆[x, :], dim = 0)) s.t. x ∈ {0, 1}m, ||x||0 ≤ k (98)
We also model the discrete clustering problem as a graph and encode the problem structure by a graph neural network. Specifically, for points that lie on the 2D plane, we define there is an edge between two points if their distance is smaller than a predefined threshold. Then the geometric information is encoded by 3-layer SplineCNN (Fey et al., 2018) which was found successful for geometric feature learning. Similar to the max covering problem, we rewrite the GPU-friendly form of the objective score in a differentiable manner, being estimated by the output of GS-TopK algorithm:
J̃i = sum(softmax(−T∆[x, :], dim = 0, keepdim = True) ∆[x, :]) J = mean([J̃1, J̃2, ..., J̃#G])
(99)
where we use a softmax operator with temperature T to approximate the min operator in a differentiable way. means the tensor-product operator in PyTorch. The GS-TopK learning algorithm is summarized in Alg. 3.

heading: C EXPERIMENT SETUP
text: Our algorithms are implemented by PyTorch and the graph neural network modules are based on (Fey & Lenssen, 2019). In our paper, we optimize the hyperparameters by greedy search on a small subset of problem instances (∼5) and set the best configuration of hyperparameters for both GS-TopK, SOFT-TopK. For the max covering problem, we empirically set learning rate = 0.1, τ = 0.05, σ = 0.15 for GS-TopK, τ = 0.05 for SOFT-TopK, and (τ, σ) = (0.05, 0.15), (0.04, 0.10), (0.03, 0.05) for Homotopy GS-TopK. For the discrete clustering problem, we set learning rate = 0.1, τ = 0.05, σ = 0.25 for GS-TopK, τ = 0.05 for SOFT-TopK, and we set (τ, σ) = (0.05, 0.25), (0.04, 0.15), (0.03, 0.05) for Homotopy GS-TopK. We set #G = 1000 for max covering, #G = 500 for discrete clustering. All experiments are done on a workstation with i7-9700K@3.60GHz, 16GB memory, and 2080Ti GPU.

heading: D ABLATION STUDY ON HYPERPARAMETERS
text: Firstly, we want to add some remarks about the selection of hyperparameters:
• #G (number of Gumbel samples): #G affects how many samples are taken during training and inference for GS-TopK. A larger #G (i.e. more samples) will be more appealing,
Algorithm 3: Gumbel-Sinkhorn-TopK Learning for Solving the Discrete Clustering Problem Input: the distance matrix ∆; learning rate α; softmax temperature T ; GS-TopK parameters
k, τ, σ, #G. 1 if Training then 2 Randomly initialize neural network weights θ; 3 if Inference then 4 Load a pretrained neural network weights θ; Jbest = +∞; 5 while not converged do 6 s = SplineCNNθ(∆); [T̃1, T̃2, ..., T̃#G] = GS-TopK(s, k, τ, σ, #G); 7 for all i,
J̃i = sum(softmax(−T∆[T̃i[2, :], :], dim = 0, keepdim = True) ∆[T̃i[2, :], :]); 8 J = mean([J̃1, J̃2, ..., J̃#G]); 9 if Training then
10 update θ with respect to the gradient ∂J∂θ and learning rate α by gradient descend; 11 if Inference then 12 update s with respect to the gradient ∂J∂s and learning rate α by gradient descend; 13 for all i, J̃i = TopK(T̃i[2, :])A · v; Jbest = min([J̃1, J̃2, ..., J̃#G], Jbest);
14 if Homotopy Inference then 15 Decline the values of τ, σ and jump to line 5;
Output: Learned network weights θ (if training)/The best objective Jbest (if inference).
because GS-TopK will have a more accurate estimation of the objective score, and it will have a higher probability of discovering better solutions. However, #G cannot be arbitrarily large because we are running GS-TopK in parallel on GPU, and the GPU has limited memory. Therefore, in experiments, we set a large enough #G (e.g. #G = 1000) and ensure that it can fit into the GPU memory of our workstation (2080Ti, 11G).
• τ (entropic regularization factor of Sinkhorn): Theoretically, τ controls the gap of the continuous Sinkhorn solution to the discrete solution, and a smaller τ will lead to a tightened gap. This property is validated by our theoretical findings in Lemma 2. Unfortunately, τ cannot be arbitrarily small, because a smaller τ requires more Sinkhorn iterations to converge. Therefore, given a fixed number of Sinkhorn iterations (100) to ensure the efficiency of our algorithm, we need trial-and-error to discover the suitable τ for both SOFT-TopK and GS-TopK. The grid search results below show that our selection of τ fairly balances the performances of both SOFT-TopK and GS-TopK.
• σ (Gumbel noise factor): As derived in Theorem 1, we need to balance between the two gaps (Lemma 1 and Lemma 2) by a suitable σ. Since |xk − xk+1| is unknown, we cannot directly find the theoretically optimal σ given τ . In the experiments, we firstly determine a τ , and then find a suitable σ by greedy search on a small subset (∼5) of problem instances.
We conduct an ablation study about the sensitivity of hyperparameters by performing an extensive grid search near the configuration used in our max covering experiments (τ = 0.05, σ = 0.15,#G = 1000). We choose the k=50, m=500, n=1000 max covering problem, and we have the following results for GS-TopK (ours) and SOFT-TopK (Xie et al., 2020) (higher is better):
Under the configuration used in our paper, both SOFT-TopK and GS-TopK have relatively good results. By grid search, we also discover better hyperparameters for both GS-TopK and SOFTTopK, but it is worth noting that SOFT-TopK is still inferior to GS-TopK with the best-discovered hyperparameters. Our grid search result shows that our GS-TopK is not very sensitive to σ if we have τ = 0.05 or 0.1, and the result of τ = 0.01 is inferior because the Sinkhorn algorithm may not converge. The results of #G = 1000 are all better than #G = 800, suggesting that a larger #G is appealing if we have enough GPU memory. It is also discovered that SOFT-TopK seems to be able to accept a smaller value of τ compared to GS-TopK, probably because adding the Gumbel noise will increase the divergence of elements thus performs in a sense similar to decreasing τ when considering the convergence of Sinkhorn.

Review summaries:
The paper proposes an unsupervised method for cardinality-constrained combinatorial problems. By building on previous work on Soft TopK through optimal transport (OT), the authors provide an enhanced Soft TopK framework that includes Gumbel noise in the OT problem. This provides a more complete theoretical bound on the original OT problem but also strong experimental performance on two combinatorial optimization problems.

This paper advances one-shot solution generation methods for cardinality constrained optimization problems. The recently developed SOFT-TopK algorithm suffers from an optimality bound discontinuity. This paper introduces a Gumbel randomization trick to the algorithm, prove that this removes the discontinuity in the bound, and implement their algorithm in two combinatorial problems, where they demonstrating achieving better solutions in less time than baselines. 

This submission studies the learning-based formulation of cardinality constrained combinatorial optimization using probabilistic methods, specifically using Gumble-Sinkhorn method. By introducing Gumble noise into the existing Soft-TopK (Xie et al., 2020) formulation, the authors claimed that the bound of the difference between approximate objective function value and the ground-truth value is tightened, compared to the results in Xie et al., 2020. The authors presented Gumbel-Sinkhorn-TopK (GS-TopK) and its application in TopK combinatorial optimization problems. Experimental results for max covering and k-means clustering were presented to demonstrate the claim.