Solving Soft Clustering Ensemble via k-Sparse Discrete Wasserstein Barycenter
Year: 2022

Abstract:
Clustering ensemble is one of the most important problems in ensemble learning. Though it has been extensively studied in the past decades, the existing methods often suffer from the issues like high computational complexity and the difficulty on understanding the consensus. In this paper, we study the more general soft clustering ensemble problem where each individual solution is a soft clustering. We connect it to the well-known discrete Wasserstein barycenter problem in geometry. Based on some novel geometric insights in high dimensions, we propose the sampling-based algorithms with provable quality guarantees. We also provide the systematical analysis on the consensus of our model. Finally, we conduct the experiments to evaluate our proposed algorithms.

heading: 1 Introduction
text: Clustering is a fundamental topic that has important applications in various areas, such as data mining, networking, and bioinformatics [34]. In the past decades, a number of different clustering objectives and algorithms have been proposed. For example, the popular k-means aims to partition the given data set into k clusters and minimize the average squared distance from the input data to the set of cluster centers; the well known k-means clustering algorithms include the Lloyd’s algorithm [47], k-means++ [5], and local search [36]. Other clustering objectives, like hierarchical clustering [49] and density-based clustering [53], are also widely used in practice.
Obviously, different clustering algorithms can obtain different results. Moreover, even for the same clustering algorithm (e.g., the Lloyd’s algorithm), the initialization and data preprocessing (e.g., random projection [12]) steps may yield different clustering results. Therefore, a natural idea is to aggregate these different clustering results so as to achieve a more reliable result. The problem is called clustering ensemble (also termed clustering aggregation or consensus clustering) [27].
However, the current methods still suffer from several issues in theory and practice. Most of existing clustering ensemble methods rely on complicated optimization models, such as the correlation clustering [28], graph partition [25], semi-definite programming [54], matrix completion [63], and spectral clustering [56]; these optimization problems usually have super-linear complexities and thus cannot be efficiently solved for large-scale datasets. Though several heuristic ideas have been proposed for speeding up the computation (e.g., the sampling idea proposed in [61]), they are in lack of rigorous theoretical analysis on their quality guarantees.
Another issue is about the interpretability of consensus. A large number of clustering ensemble models are based on utility function [57, 61, 60] or co-association matrix [26, 56]. From the theoretical perspective, a fundamental question is why these models can yield the final clusterings close to
∗Corresponding author.
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
the ground-truth clustering. The similar consensus question was also studied for the classification problem before [13]. However, the analysis for clustering ensemble is more challenging, because we need to take into account the matchings between different clustering solutions.

heading: 1.1 Our Contributions
text: In this paper, we focus on the more general soft clustering ensemble problem, where each given individual clustering is a soft clustering (also referred to as fuzzy clustering) [61]. In a soft clustering, data points can potentially belong to multiple clusters. For example, a point may be assigned to three clusters with the probabilities 10%, 20%, and 70%, respectively. Compared with hard clustering, soft clustering can provide more realistic and accurate clustering results in many realworld applications [9].
We adopt the geometric model that was studied in [18, 20, 19]. They showed that clustering ensemble can be naturally formulated as a “geometric prototype” problem. But their results are in lack of systematically studies on the efficiency of this model, especially from the theoretical perspective. In this paper, we illustrate that the geometric prototype actually is equivalent with an instance of Discrete Wasserstein Barycenter (DWB) [1] in high dimensions (the formal definition will be given in Section 3). This approach falls under the umbrella of utility function based model, where it uses discrete Wasserstein distance to measure the difference between two clusterings. Compared with other utility functions (e.g., KL-divergence), it has several attractive properties. For example, the discrete Wasserstein distance is symmetric and more robust to noise [42]. More importantly, the DWB based ensemble model can be easily interpreted from the geometric perspective, and thus we can analyze its performance more conveniently. But when applying the DWB model to the soft clustering ensemble problem, we still have several key issues remaining to be solved.
(i) Though a number of DWB algorithms have been developed (as shown in Section 1.2), the clustering ensemble imposes two unique features to the DWB formulation. First, we require the returned DWB to be k-sparse, that is, it should be supported by at most k points in the space (since there are at most k clusters). Also, the number of different clustering solutions can be large in practical scenarios. For instance, to guarantee the consistency of the final ensemble solution to the ground-truth clustering, we may generate a large number of clustering solutions via random initializations or random projections [24]. So from the algorithmic perspective, a natural question is whether we can develop more efficient algorithm for the DWB problem with such features?
(ii) To the best of our knowledge, only Topchy et al. [58] and Jain [35] discussed the consensus of clustering ensemble in theory. However, both of their analyses rely on the assumption that the ground-truth clustering should be the optimal solution of the ensemble model, which is too strong and may not be realistic in practice. Also, the analysis of [58] only considered hard clustering. It is worth noting that the number of hard clusterings on a fixed set of items is finite, but the number of soft clusterings is infinity.
In this paper, we focus on these two issues. First, based on some novel geometric insights, we show that it is possible to achieve a fixed-parameter algorithm for the soft clustering ensemble problem if k is a constant, where the obtained approximation factor is 1 + with being any small number in (0, 1); though this is more a theoretical result, we believe that it is of independent interest for such a combinatorial optimization problem in high dimensions [17]. Moreover, the proposed sampling idea inspires our following speedup for the existing DWB algorithms with provable quality guarantee, even if k is large. Second, we prove that the obtained DWB should be close to the ground-truth clustering if the number of given clustering solutions is large enough. Our idea is quite different from [58, 35]; in particular, our analysis yields a detailed quantitive result for the consensus.

heading: 1.2 Related Works
text: Clustering ensemble. Clustering ensemble was introduced by Strehl and Ghosh [55], where they formulated it as two different graph partitioning problems; one is “Instance-Based Graph Formulation (IBGF)”, and the other is “Cluster-Based Graph Formulation (CBGF)”. Fern and Brodley [25] later proposed a hybrid graph partitioning model for clustering ensemble. Most existing clustering ensemble models can be roughly divided into utility function based [57, 61, 60] and co-association matrix based [26, 46, 56]. The utility function based model is to find the final clustering result via maximizing the total similarities to the set of given clusterings. A co-association matrix is actually
a new representation of the data items, where each entry of the matrix indicates the similarity of a pair of data items based on the information from the given clusterings. The soft clustering ensemble problem was particularly studied in [51, 61]. For more detailed discussion about ensemble clustering, the reader is referred to the survey [27].
Wasserstein distance. The Wasserstein distance is defined for measuring the difference between two probability distributions; if their supports are both discrete sets, the distance is called discrete Wasserstein distance (or Earth Mover’s distance) [59, 52]. Computing the discrete Wasserstein distance actually is equivalent to solving a min-cost max flow problem [2, 37]. Several more efficient discrete Wasserstein distance algorithms were proposed, such as [45, 50]. In the community of machine learning, Cuturi [15] proposed a new objective called “Sinkhorn Distance” that smoothes the transportation problem with an entropic regularization term, and it can be solved much faster than computing the exact discrete Wasserstein distance. Following Cuturi’s work, several improved Sinkhorn algorithms have been proposed in recent years [23, 44, 3, 48].
Wasserstein barycenter. If there are m ≥ 2 different weighted point sets, the problem of discrete Wasserstein barycenter is to compute the average pattern that minimizes the total Wasserstein distances to them [1]. The recent algorithms for computing Wasserstein barycenter include [16, 7, 29, 62, 8, 4, 43]. The computational complexity of Wasserstein barycenter was studied in [39, 11, 43]. Recently, Dognin et al. [21] applied the Wasserstein barycenter based ensemble method to solve several supervised learning problems.

heading: 2 Preliminaries
text: We always use A = {a1, a2, · · · , an} to denote the set of n data items which we want to cluster. Definition 1 (Soft Clustering). Let k ∈ Z+. A k-soft clustering of A can be represented by a set of k vectors C = {S1, · · · , Sk} ⊂ [0, 1]n, where each vector Sj represents an individual cluster and∑k j=1 Sj is equal to the row vector [1, 1, · · · , 1]. Suppose Sjl is the l-th entry of Sj for 1 ≤ l ≤ n. Then, for a fixed l, the set {S1l, · · · , Skl} indicates the degrees of membership of al to the k clusters.
For example, suppose n = 3 and k = 2; the following is a k-soft clustering of A: C = {S1 = [0.1, 0.3, 0.8], S2 = [0.9, 0.7, 0.2]} (e.g., the item a2 belongs to the first and second clusters with probability 30% and 70%, respectively).
In Definition 1, if we restrict each Sj to be binary vector, the clustering C will be a hard clustering. Also, if a clustering C has k′ < k clusters, we can simply add k − k′ dummy clusters where each dummy cluster is just a zero vector [0, 0, · · · , 0]. Following the idea of [18, 6, 20], we define the function ∆ to measure the difference between two clusterings of A. Suppose C = {S1, · · · , Sk} and C′ = {S′1, · · · , S′k} are two different soft clusterings. We build the bipartite graph G from C and C′ as follows: each of the two columns of G contains k vertices corresponding to the k clusters in C and C′ respectively; for any pair of clusters (Sj , S ′ j′) with Sj ∈ C and S′j′ ∈ C′, there is an edge connecting their corresponding vertices in G with a weight equal to their squared Euclidean distance ||Sj − S′j′ ||2. Then, the difference of C and C′, i.e., ∆(C, C′), is the cost of the minimum weight bipartite matching of G. The function ∆ can be computed through the Hungarian algorithm [14]. Assume the minimum weight bipartite matching between C and C′ yields the permutation π of {1, 2, · · · , k}. Namely,
∆(C, C′) = k∑ j=1 ||Sj − S′π(j)|| 2. (1)
Obviously, if ∆(C, C′) = 0, they should be the same clustering solution (the k clusters of C is just reordered by π in C′). To see the rationale behind (1), we can imagine the case that C and C′ are both hard clusterings, i.e., each Sj ∈ C (resp., S′j′ ∈ C′) is a binary vector; so each such vector can be viewed as a subset of A (e.g., Sj represents the set {al | al ∈ A and Sjl = 1}). It is easy to know the symmetric difference of Sj and S′j′ , |Sj \ S′j′ |+ |S′j′ \ Sj |, is equal to ||Sj − Sj′ ||2. Therefore, computing the function ∆(C, C′) in fact is to find the matching of C and C′ that minimizes their total symmetric differences.
Following Definition 1, we then introduce “soft clustering ensemble” below.
Definition 2 (Soft Clustering Ensemble (SCE)). Given m different soft clusterings C1, · · · , Cm of A, the problem of soft clustering ensemble (SCE) is to find the final soft clustering C̃ that minimizes the objective function
1
m m∑ i=1 ∆(C̃, Ci). (2)
For any soft clustering C̃′ and λ ≥ 1, if it achieves an objective value no larger than λ times the minimum value of (2), we say it is a “λ-approximation” for the SCE problem.
Let πi be the permutation between C̃ and Ci for 1 ≤ i ≤ m. To minimize the objective function (2), the major challenge is to find thesem permutations simultaneously. Suppose each Ci = {Si1, · · · , Sik}. Once these m permutations are obtained, the set ∪mi=1Ci is divided into k parts:
{S1π1(j), · · · , S m πm(j) }, 1 ≤ j ≤ k. (3)
If we let the optimal solution C̃ be {S̃1, · · · , S̃k}, from the objective function (2) it is easy to know these k soft clusters should be the means of these k parts, i.e.,
S̃j = 1
m m∑ i=1 Siπi(j), 1 ≤ j ≤ k. (4)
This simple fact will be used in the analysis in our paper. We also have the following hardness result for SCE through the reduction from the NP-hard three-dimensional assignment problem [41]. Theorem 1 (The hardness). When m ≥ 3, optimizing the SCE objective (2) is NP-hard.
The rest of the paper is organized as follows. In Section 3, we discuss the relation between the SCE problem and discrete Wasserstein barycenter. In Section 4, we present our approximation algorithms based on random sampling. In Section 5, we analyze the consensus of the SCE problem under Definition 2. Finally, we illustrate our experimental results in Section 6. Due to the space limit, we leave some proofs and the detailed experimental results to the full version of this paper.

heading: 3 Relation With Discrete Wasserstein Barycenter
text: We introduce the relation between SCE and discrete Wasserstein barycenter in this section. Definition 3 (Discrete Wasserstein Distance [52]). Let P = {p1, p2, · · · , pnP } and Q = {q1, q2, · · · , qnQ} be two sets of weighted points in Rd with nonnegative weights αi and βj for each pi ∈ P and qj ∈ Q, and ∑nP i=1 αi = ∑nQ j=1 βj = 1. Their discrete Wasserstein distance is
Wss (P,Q) = min F nP∑ i=1 nQ∑ j=1 fij ||pi − qj ||ss, (5)
where || · ||s indicates the ls-norm and F = {fij | 1 ≤ i ≤ nP , 1 ≤ j ≤ nQ} is a feasible flow from P to Q, i.e., each fij ≥ 0, ∑nP i=1 fij = βj , and ∑nQ j=1 fij = αi.
In this paper, we only focus on the l2-discrete Wasserstein distance. The l2-discrete Wasserstein barycenter (DWB2) considers the following objective. Given the nonnegative weighted point sets {P1, P2, · · · , Pm} in Rd, where each Pi has the total weights equal to 1, the goal is to find a set of centroid points P̃ , so as to minimize
1
m m∑ i=1 W22 (Pi, P̃ ). (6)
Further, if we require P̃ to have at most k points with some k ∈ Z+, the problem is called “k-sparse DWB2” [11]. It is easy to observe that the k-sparse DWB2 problem is very similar to SCE as described in Definition 2. Intuitively, the obtained optimal ensemble clustering C̃ can be viewed as the k-sparse DWB2, if each point of Ci is assigned the weight 1/k for 1 ≤ i ≤ m.
Theorem 2. Given a set of soft clusterings {C1, · · · , Cm}, the optimal solution of the SCE problem (2) is exactly their optimal k-sparse DWB2. Moreover, for any λ ≥ 1, a λ-approximation of the k-sparse DWB2 can yield a λ-approximate solution of (2).
To prove the above theorem, we need to consider two issues. First, the discrete Wasserstein distance of Definition 3 may result in a many-to-many matching between C̃ and Ci, but the difference function ∆(C̃, Ci) requires a one-to-one matching. Since the discrete Wasserstein distance is actually an instance of the min-cost max flow problem, we can convert the obtained many-to-many matching to a one-to-one matching without increasing the complexity. Second, we should further prove that the obtained (optimal or approximate) k-sparse DWB2 is a feasible soft clustering (i.e., it should satisfy the conditions described in Definition 1).

heading: 4 Approximation Algorithms
text: Due to Theorem 1, we only focus on the approximation algorithms for SCE. It is also worth noting that when m = 2, the optimal k-sparse DWB2 (and the SCE solution) can be easily obtained by computing the matching between C1 and C2; the optimal barycenter should just be the midpoints of the k matched pairs. In this section, we propose the approximation algorithms for optimizing the SCE objective (2) with m ≥ 3. First, we present the following lemma which is the key to our algorithms. Given a point set Q ⊂ Rd, we use µ(Q) and Var(Q) to denote the mean and variance respectively, i.e., µ(Q) = 1|Q| ∑ q∈Q q
and Var(Q) = 1|Q| ∑ q∈Q ||q − µ(Q)||2. Lemma 1. [33] Let δ ∈ (0, 1). Given a point set Q, we suppose that Q′ is a set of t points sampled from Q uniformly at random. Then with probability 1− δ, ||µ(Q)− µ(Q′)||2 ≤ 1δtVar(Q). Remark 1. Lemma 1 reveals that we can estimate µ(Q) by just simple random sampling. For example, if we require the error no larger than Var(Q) with some small > 0, we only need to sample 1δ points from Q and the success probability is 1 − δ. Obviously, the smaller and δ, the larger the required sample size. Moreover, a highlight of Lemma 1 is that the sample size is independent of the dimensionality.
We also need the following lemma from [40]. Lemma 2 indicates that if a point p is close to µ(Q), the total squared distances to the points of Q should be close to Var(Q) as well. Lemma 2. Let Q be a set of points in Rd. For any point p ∈ Rd, 1|Q| ∑ q∈Q ||q − p||2 = Var(Q) + ||µ(Q)− p||2.
In Section 4.1, we propose a fixed-parameter algorithm that returns a (1 + )-approximate solution, if k is assumed to be small. Though it is more a theoretical result, the proposed sampling idea inspires our following improvement in Section 4.2 on the existing alternating minimization algorithms for the case that k is not a constant.
4.1 A Fixed-parameter Algorithm
When k is small, we can achieve a (1 + )- approximation for the SCE problem. We briefly illustrate our idea below.
For ease of presentation, we “temporarily” assume the permutations πi between C̃ and Ci, 1 ≤ i ≤ m, which yield the optimal solution of (2), are given at this moment. For each Ci, we concatenate its k vectors {Si1, · · · , Sik} to be a long vector
vi = (S i πi(1) · · ·Siπi(k)) (7)
in Rkn. See Figure 1 for an illustration. Meanwhile, the ensemble clustering C̃ = {S̃1, · · · , S̃k} can be also represented as a vector
ṽ = (S̃1 · · · S̃k) (8)
Algorithm 1 (1 + )-Approximate SCE Algorithm Input: m k-soft clusterings {C1, C2, · · · , Cm} on a set A of n data items, and two parameters 0 < δ, ≤ 1.
1. Select t = 1δ clusterings uniformly at random, and denote them by Ci1 , Ci2 , · · · , Cit . For 1 ≤ l ≤ t, Cil contains k vectors (i.e., soft clusters) {S il 1 , · · · , S il k } in [0, 1]n. 2. Initiate a candidate set G = ∅. 3. For each Cil , 1 ≤ l ≤ t, enumerate all the k! possible permutations for πil (so there are
(k!)t cases in total). For each case: (a) Let C = {S1, S2, · · · , Sk}, where each
Sj = 1
t t∑ l=1 Silπil (j) , j = 1, 2, · · · , k.
(b) Update G = G ∪ {C}. 4. For each candidate C ∈ G, compute the objective value 1m m∑ i=1 ∆(C, Ci). Let C̄ be the one
with the smallest objective value among G. Output: C̄ as the final solution.
in Rkn. We show that ṽ in fact is the mean of {v1, · · · , vm}. Then we can apply Lemma 1 to estimate the position of ṽ. Finally, since k is assumed to be a constant, we can enumerate all the possible permutations of the sampled clusterings and choose the best one as the final solution. Though Theorem 3 is more a theoretical result, we believe that it is of independent interest as a fixed-parameter solution for such a combinatorial optimization problem in high dimensions. Theorem 3. (i) With probability 1− δ, Algorithm 1 yields a (1 + )-approximate solution of the SCE problem. (ii) The runtime is O(exp( 1δ k log k) · k 2 ·mn)
Proof. First, we need to show that each candidate C ∈ G is a feasible soft clustering, that is, its k vectors {S1, · · · , Sk} generated in Step 3(a) should satisfy the conditions of Definition 1. Since each Sj = 1 t ∑t l=1 S il πil (j) , it must belong to [0, 1]n. Also,
k∑ j=1 Sj = k∑ j=1 1 t t∑ l=1 Silπil (j) = 1 t t∑ l=1 k∑ j=1 Silπil (j) = [1, 1 · · · , 1], (9)
where the final equality comes from the fact that each sampled Cil is a feasible soft clustering. Consequently, each candidate C from G is also a feasible soft clustering. Now, we consider the induced objective value of the best candidate selected from G. As discussed before (see Figure 1), each soft clustering can be converted to a long vector in Rkn. Denote by vi the vector of Ci, for 1 ≤ i ≤ m; similarly, denote by ṽ the vector of the optimal solution C̃. We use V to denote the set {v1, · · · , vm}. Obviously, ṽ = µ(V ) (from the fact (4)). Since {Ci1 , Ci2 , · · · , Cit} are the randomly selected t = 1δ clusterings from the input, together with Lemma 1, we have
||ṽ − ṽ′||2 ≤ Var(V ) (10) with probability 1− δ, where ṽ′ = 1t ∑t l=1 vil . Since 1 m ∑m i=1 ||vi − ṽ′||2 = Var(V ) + ||ṽ − ṽ′||2 (by Lemma 2), the inequality (10) implies
1
m m∑ i=1 ||vi − ṽ′||2 ≤ (1 + )Var(V ). (11)
Hence, ṽ′ is a (1 + )-approximate solution of the objective (2). Because we do not know those permutations πil , 1 ≤ l ≤ t, we cannot directly obtain ṽ′. Through enumerating all the (k!)t cases, we can claim that there must exist one candidate in G that yields a (1 + )-approximation. The time complexity contains two parts, i.e., constructing |G| and selecting the best candidate from G. Since |G| = (k!)t = O(exp( 1δ k log k)), the first part takes O(exp( 1 δ k log k) · t · kn) time.
For the second part, we use the Hungarian algorithm to compute the one-to-one matching from Ci, 1 ≤ i ≤ m, to each candidate of G. Note that the complexity of the Hungarian algorithm is O(k3+k2·n) = O(k2n) (n k), where the term k2n is due to the time for building the k×k bipartite graph in Rn. Therefore the second part takesO(exp( 1δ k log k)·k
2·mn) time. Overall, the second part dominates the whole complexity, and the runtime of Algorithm 1 is O(exp( 1δ k log k) · k 2 ·mn).

heading: 4.2 When k Is Not a Constant
text: Now we consider the case that k is not a constant. Due to the discussion in Section 3, we can directly apply any existing k-sparse DWB2 algorithm to compute the final ensemble clustering. An interesting observation is that a couple of widely used DWB2 algorithms follow the alternating minimization framework [16, 7, 62] (several previous ensemble clustering algorithms also used this alternating minimization idea [31, 20]): in each iteration, the algorithm alternatively performs the following two steps:
1. update the matchings (i.e., the Wasserstein flows) from those Cis to the temporary barycenter; 2. update the k points of the temporary barycenter based on the new matchings.
Different algorithms may adopt different strategies for implementing these two steps. Eventually, the result will converge to a local optimum (note it is NP-hard to achieve a global optimum according to Theorem 1). A bottleneck of this framework, especially when m is large, is the computation for the Wasserstein distances over all the m clusterings C1, · · · , Cm. In fact, similar to Algorithm 1, we can still apply Lemma 1 to reduce the time complexity for this bottleneck. Let C̃⊥ be the temporary barycenter at the beginning of the current iteration. Suppose C̃> is the updated barycenter if we compute all the m Wasserstein distances. Let 0 ≤ , δ ≤ 1. If we randomly select 1 δ clusterings and only compute these 1 δ Wasserstein distances to update C̃⊥ to be C̃ ′ >, we have the following result via the same idea of Theorem 3. Lemma 3. With probability 1− δ, 1m ∑m i=1 ∆(C̃′>, Ci) ≤ (1 + ) · 1 m ∑m i=1 ∆(C̃>, Ci).
Lemma 3 indicates that we can avoid computing all the m Wasserstein distances and still achieve a barycenter close to C̃> via random sampling. Note that a key difference with Algorithm 1 is that we do not need to enumerate all the permutations since the 1 δ matchings can be determined by C̃⊥.

heading: 5 Analysis on the Consensus
text: Denote by Cgt the ground-truth clustering over the n data items of A. In this section, we analyze the consensus of the objective (2). Specifically, when the number m increases, we are wondering whether the optimal solution C̃ will converge to the ground-truth clustering Cgt. Let Ω be the set of all the possible soft clusterings over A, following a probability measure ρ. Namely, for any soft clustering C ∈ Ω, the probability of obtaining it is ρ(C), and ∫ C∈Ω ρ(C)dC = 1.
First, we propose the following assumption that there is an upper bound for the difference between Cgt and any C ∈ Ω. Assumption 1. There exists a value L > 0, such that ∆(C, Cgt) ≤ L for any C ∈ Ω.
Actually, this assumption is easy to understand in the context of clustering ensemble. Since each clustering solution C is obtained by some reasonable clustering algorithm, it makes sense to assume that these solutions are not arbitrarily far from Cgt. Assumption 1 directly implies the following Lemma which can be proved by the fact ∆(C, C′) ≤ 2∆(C, Cgt) + 2∆(C′, Cgt). Lemma 4. For any C and C′ ∈ Ω, ∆(C, C′) ≤ 4L.
As mentioned in Section 1.1, to the best of our knowledge, only Topchy et al. [58] and Jain [35] discussed the consensus of clustering ensemble in theory. However, both of their analyses need the assumption that the ground-truth clustering should be exactly the one achieving the smallest expected objective value over all the possible clustering solutions:
Cgt = arg min Ĉ∈Ω ∫ Ω ρ(C)∆(Ĉ, C)dC, (12)
which may be too strong in reality. Also, they did not provide the quantitive analysis, e.g., the numerical relation between the convergence and the value m. Here, we relax the assumption (12). In particular, we allow Cgt 6= arg minĈ∈Ω ∫ Ω ρ(C)∆(Ĉ, C)dC. We instead assume any clustering solution that achieves sufficiently low expected objective value, should be close to Cgt. Let Opt := minĈ∈Ω ∫ Ω ρ(C)∆(Ĉ, C)dC, and Copt be the clustering solution that achieves Opt. Note that Copt and Cgt are not necessary to be the same one, and our ultimate goal is to find a solution close to Cgt.
Assumption 2. There exist two numbers c, ξ ≥ 0, such that for any Ĉ ∈ Ω, if its expected objective value ∫ Ω ρ(C)∆(Ĉ, C)dC ≤ (1 + c)Opt, we have ∆(Ĉ, Cgt) ≤ ξ.
Obviously, when c = ξ = 0, Assumption 2 will be as same as the assumption (12) from [58]. So our assumption is more relaxed. We also need to point out that the objective (2) (and the k-sparse Wasserstein barycenter) may have multiple isolated global optimums. But under Assumption 1, we restrict Ω to a local region and thus it is reasonable to assume Assumption 2 to be true.
To begin analyzing the consensus, we fix a clustering Ĉ ∈ Ω first. We let xi = ∆(Ĉ, Ci) for i = 1, 2, · · · ,m, and view each xi as an independent random variable. Then, we denote their mean 1 m ∑m i=1 xi as x̄. Obviously, E[x̄] = ∫ Ω ρ(C)∆(Ĉ, C)dC. We then study the difference between x̄ and E[x̄]. From Lemma 4, we know xi ∈ [0, 4L] for 1 ≤ i ≤ m. Through the Hoeffding’s inequality [30], for any η > 0 we have
Prob [ |x̄− E [x̄] | > ηE [x̄] ] < 2 exp(−mη 2(E[x̄])2
2L2 ).
Consequently, we have the following Lemma.
Lemma 5. Fix Ĉ ∈ Ω. For any η, δ ∈ (0, 1), if m ≥ 8L 2
η2(E[x̄])2 log 2 δ , with probability 1 − δ,
1 m ∑m i=1 ∆(Ĉ, Ci) ∈ (1± η) ∫ Ω ρ(C)∆(Ĉ, C)dC.
But Lemma 5 is only for a fixed Ĉ. We need to extend the Lemma to any Ĉ ∈ Ω. To realize this goal, we discretize the set Ω first.
Discretization. We use B(p, r) to denote the ball centered at a point p with radius r ≥ 0 in Rn. We suppose the k soft clusters (i.e., the k vectors in Rn) of Cgt are Sgt,1, · · · , Sgt,k. From Assumption 1, we know for each Ĉ ∈ Ω, its k vectors should be covered by the region R = ∪kj=1B(Sgt,j , √ L). Imagine that we draw a uniform grid insideR with the grid side length being equal to ϑ√ n
(the value of ϑ will be determine later). Thus, for any two points inside the same cell of the grid, their distance is no larger than √ n · (ϑ/ √ n)2 = ϑ. Moreover, by using the formula for ball volume in Rn, we know the size of Γj , which denotes the set of the grid points inside B(Sgt,j , √ L), is O ((√ πeL ϑ )n) .
So the set Ωgrid := Γ1 × Γ2 × · · · × Γk contains N = O ((√ πeL ϑ )kn) different k-tuple points (i.e.,
soft clusterings). If we replace δ by δ/N in Lemma 5 and take the union bound over all the soft clusterings of Ωgrid, we can obtain the following result.
Lemma 6. For any η, δ ∈ (0, 1), if m ≥ O ( knL2
η2(Opt)2 log L ϑδ
) , with probability 1 − δ,
1 m ∑m i=1 ∆(Ĉ, Ci) ∈ (1± η) ∫ Ω ρ(C)∆(Ĉ, C)dC for any Ĉ ∈ Ωgrid.
Remark 2. We replace E [x̄] by Opt in the lower bound of m in Lemma 6, since E [x̄] is always no smaller than Opt.
Now, we consider the clusterings in Ω \ Ωgrid. For each point p ∈ R, we use N (p) to denote its nearest grid point. Similarly, for any Ĉ = {Ŝ1, · · · , Ŝk} ∈ Ω \ Ωgrid, we denote its “nearest clustering” as N (Ĉ), which contains the k vectors {N (Ŝ1), · · · ,N (Ŝk)}. We should prove that∣∣∆(N (Ĉ), C)−∆(Ĉ, C)∣∣ is small for any C ∈ Ω, as long as ϑ is sufficiently small. Consequently, we can extend the result of Lemma 6 from Ωgrid to all the clusterings in Ω.
Lemma 7. For any η, δ ∈ (0, 1), if m ≥ O ( knL2
η2(Opt)2 log kL ηδOpt
) , with probability 1 − δ,
1 m ∑m i=1 ∆(Ĉ, Ci) ∈ (1± 7η) ∫ Ω ρ(C)∆(Ĉ, C)dC for any Ĉ ∈ Ω.
The value of ϑ is set to be η·Opt√ 8kL
in the proof of Lemma 7 (the detailed proof is placed to our full paper). Finally, we achieve the consensus theorem under Assumption 1 and 2. Theorem 4 (Consensus). Let δ ∈ (0, 1) and η < c7(2+c) . Suppose {C1, · · · , Cm} are drawn i.i.d
from Ω with the probability ρ(·). If m ≥ O ( knL2
η2(Opt)2 log kL ηδOpt ) , with probability 1− δ, ∆(C̃, Cgt) is
no larger than ξ, where C̃ is the optimal solution for the objective (2). Proof. Suppose ∆(C̃, Cgt) > ξ. From Assumption 2, we know ∫ Ω ρ(C)∆(C̃, C)dC > (1 + c)Opt, which implies
1
m m∑ i=1 ∆(C̃, Ci) ≥ (1− 7η)(1 + c)Opt (13)
via Lemma 7. Moreover, by using Lemma 7 again, we have
1
m m∑ i=1 ∆(Copt, Ci) ≤ (1 + 7η)Opt. (14)
Since we let η < c7(2+c) , (13) and (14) together imply 1 m ∑m i=1 ∆(Copt, Ci) < 1 m ∑m i=1 ∆(C̃, Ci), which is contradict with the fact that C̃ is the optimal solution for the objective (2). Hence the inequality ∆(C̃, Cgt) ≤ ξ is true.

heading: 6 Experimental Results
text: We evaluate the practical performance of our proposed algorithm in this section. All the experimental results were obtained on a server equipped with 2.8GHz Intel CPU, 8GB main memory, and Matlab 2019a. We consider three real datasets: USPS has 11000 data items in R256 with k = 10 [32]; IRIS [22] has 150 data items in R4 with k = 3; CIFAR-10 [38] has 10000 data items in R3072 with k = 10. Similar with [24, 13], we apply random projections to generate the clustering solutions (in each random subspace, we use k-means to cluster the data). We consider two representative baselines: the bipartite graph partition method BGP [25]; the top-down method FURTHEST [28]. Our sampling idea of Section 4.2 is incorporated into the alternating minimization Wasserstein barycenter algorithm [62], which is denoted as AM-r with r representing the sample rate (e.g., AM-1 means we directly run the algorithm on the original data without sampling).
We set m = 1000 (i.e., the number of generated clustering solutions for ensemble) and show the results in Figure 2 (a)-(i). We can see our Wasserstein barycenter based algorithm significantly outperforms other baselines in terms of the objective value (2), the Wasserstein distance to ground truth, and the runtime on the datasets USPS and CIFAR-10. Only for the smallest dataset IRIS, the baselines are faster (actually all the runtimes are very close for this small dataset). Also, we study the convergence of our obtained result, i.e., its Wasserstein distance to the ground-truth clustering as m increases, in Figure 2 (j)-(l). We can see in general the convergence performs better when the sample rate is larger. Due to the space limit, we leave the detailed experimental results to our full paper.

heading: 7 Conclusion and Future Work
text: In this paper, we connect the soft clustering ensemble problem to k-sparse discrete Wasserstein barycenter. There are several interesting problems deserved to study in future. For example, the robustness of the Wasserstein barycenter based clustering ensemble is in lack of discussions so far. In particular, we can consider its robustness under adversarial attacks, e.g., the poisoning and evasion attacks [10]. Also, we believe it is important to study some other relevant issues, e.g., the privacy-preserving problem and the fairness problem, for clustering ensemble.

heading: Acknowledgment
text: We would like to thank the anonymous reviewers for their helpful suggestions and comments.

Review summaries:
This manuscript shows the connection between soft clustering ensemble and the discrete Wasserstein barycenter problem. Based on the connection, a sampling-based algorithm is proposed with some provable quality guarantees. The manuscript provides some theoretical analysis about the quality of the obtained consensus clustering results.

This paper proposes a new approach for soft clustering ensemble, using a new problem formulation (a “geometric prototype”), which is equivalent to a Discrete Wasserstein Barycenter problem.  The proposed algorithm achieves, with a fixed probability 1 – δ, a (1+ ε)-approximate solution to the soft clustering ensemble problem, when the number of clusters k is fixed and small. As the authors state, this is a theoretical result since the runtime is more than exponential in k, 1/ δ, and 1/ε.  In addition, the authors propose an improvement on the existing alternating minimization algorithms in the case of k not being a constant. Then, they prove that the obtained solution converges to the ground-truth clustering as the number m of clustering increases.  Finally, they evaluate the performance of the alternating minimization Wasserstein barycenter algorithm with the proposed improvement, and they compare it to three baselines from the literature. The results show that their algorithm outperforms the baselines in terms of the objective function value and, for the largest datasets, also in terms of runtime. 

This paper studies the ensemble clustering problem with soft clustering indicators. They reveal an interesting equivalence between the soft clustering ensemble problem and the $k$-sparse discrete Wasserstein barycenter. The authors proved the hardness result of SCE and reported approximation algorithms for two different regimes (fixed $k$ or not). On the modeling side, they provided consensus results for the SCE problem under a weaker assumption on the ground-truth clustering set than prior results. They also conducted numerical experiments to support their results.

Authors connected the soft clustering ensemble and discrete Wasserstein barycenters. They proposed a (1+eps) approximation algorithm and adopted a uniform sampling scheme to reduce the time complexity of computing Wasserstein barycenters. Convergence is proved and verified via experiments.