Clustering problems arise from variety of applications, such as documents/web pages clustering 
, pattern recognition, image linking, image segmentation 
, data compression via vector quantization and nearest neighbor search [5, 6, 7]. In the last three decades, various clustering algorithms have been proposed. Among these algorithms, k-means  remains a popular choice for its simplicity, efficiency and moderate but stable performance across different problems. It was known as one of top ten most popular algorithms in data mining . On one hand, k-means has been widely adopted in different applications. On the other hand, continuous efforts have been devoted to enhance the performance k-means as well.
Despite its popularity, it actually suffers from several latent issues. Although the time complexity is linear to data size, traditional k-means is still not sufficiently efficient to handle the web-scale data. In some specific scenarios, the running time of k-means could be even exponential in the worst case [10, 11]. Moreover, k-means usually only converges to local optima. As a consequence, recent research has been working on either improving its clustering quality [12, 13] or efficiency [14, 15, 16, 13, 17, 18]. K-means has been also tailored to perform web-scale image clustering [2, 19].
There are in general three steps involved in the clustering procedure. Namely, 1. initialize k cluster centroids; 2. assign each sample to its closest centroid; 3. recompute cluster centroids with assignments produced in Step 2 and go back to Step 2 until convergence. This is known as Lloyd iteration procedure . The iteration repeats Step 2 and Step 3 until the centroids do not change between two consecutive rounds. Given are cluster centroids, are samples to be clustered, above procedure essentially minimizes the following objective function:
In Eqn. 1, function q() returns the closest centroid for sample . Unfortunately, searching an optimal solution for the above objective function is NP-hard. In general k-means only converges to local minimum . The reason that k-means maintains its popularity is mainly due to its linear complexity in terms of the number of samples to be clustered. The complexity is , given as the number of iterations to converge. Compared with other well-known clustering algorithms such as DBSCAN  and Mean shift , this complexity is considerably low. However, the efficiency of traditional k-means cannot cope with the massive growth of data in Internet. In particular, in the case that the size of data (m), the number of clusters (k) and the dimension (d) are all very large, k-means becomes unbearably slow. The existing efforts [16, 18] in enhancing the scalability of k-means for web-scale tasks often come with price of lower clustering quality. On the other hand, k-means++ proposed in [12, 17] focuses on enhancing the clustering quality by a careful design of the initialization procedure. However, k-means is slow down as a few rounds of scanning over the dataset is still necessary in the initialization.
In this paper, a novel variant of k-means is proposed, which aims to make a better trade-off between clustering quality and efficiency. Inspired by the work in , a novel objective function is derived from Eqn. 1. With the development of this objective function, the traditional k-means iteration procedure has been revised to a simpler form, in which the costly initial assignment becomes unnecessary. In addition, driven by the objective function, sample is moved from one cluster to another cluster when we find this movement leads to higher objective function score, which is known as incremental clustering [1, 23]. These modifications lead to several advantages.
K-means clustering without initial assignment results in better quality as well as higher speed efficiency.
K-means iteration driven by an explicit objective function converges to considerably lower clustering distortion in faster pace.
Different from traditional k-means, it is not necessary to assign a sample to its closest centroid in each iteration, which also leads to higher speed.
In addition, when clustering in hierarchical bisecting fashion, the proposed method achieves the highest scalability among all top-down hierarchical clustering methods. Extensive experiments are conducted to contrast the performance of proposed method with k-means and its variants including tasks document clustering, nearest neighbor search (NNS) with product quantization  and image clustering.
The remainder of this paper is organized as follows. The reviews about representative works on improving the performance of traditional k-means are presented in Section 2. In Section 3, the clustering objective functions are derived based on Eqn. 1. Based on the objective function, Section 4 presents the clustering method. Extensive experiement studies over proposed clustering method are presented in Section 5. Section 6 concludes the paper.
2 Related Works
Clustering is a process of partitioning a set of samples into a number of groups without any supervised training. Due to its versatility in different contexts, it has been studied in the last three decades . As the introduction of Web 2.0, millions of data in Internet has been generated on a daily basis. Clustering becomes one of the basic tools to process such big volume of data. As a consequence, traditional clustering methods have been shed with new light. People are searching for clustering methods that are scalable [16, 17, 18] to web-scale data. In general, boosting the performance of traditional k-means becomes the major trend due to its simplicity and relative higher efficiency over other clustering methods.
In general, there are two major ways to enhance the performance of k-means. For the first kind, the aim is to improve the clustering quality. One of the important work comes from Ostrovsky et al. [12, 17]. The motivation is based on the observation that k-means converges to a better local optima if the initial cluster centroids are carefully selected. According to , k-means iteration also converges faster due to the careful selection on the initial cluster centroids. However, in order to adapt the initial centroids to the data distribution, k rounds of scanning over the data are necessary. Although the number of scanning rounds has been reduced to a few in , the extra computational cost is still inevitable.
In each k-means iteration, the processing bottleneck is the operation of assigning each sample to its closest centroid. The iteration becomes unbearably slow when both the size and the dimension of the data are very large. Noticed that this is a nearest neighbor search problem, Kanungo et al.  proposed to index dataset in a KD Tree  to speed-up the sample-to-centroid nearest neighbor search. However, this is only feasible when the dimension of data is in few tens. Similar scheme has been adopted by Dan et al. 
. Unfortunately, due to the curse of dimensionality, this method becomes ineffective when the dimension of data grows to a few hundreds. A recent work takes similar way to speed-up the nearest neighbor search by indexing dataset with inverted file structure. During the iteration, each centroid is queried against all the indexed data. Attributing to the efficiency of inverted file structure, one to two orders of magnitude speed-up is observed. However, inverted file indexing structure is only effective for sparse vectors.
only pick a small portion of the whole dataset to update the cluster centroids each time. For the sake of speed efficiency, the number of iterations is empirically set to small value. It is therefore possible that the clustering terminates without a single pass over the whole dataset, which leads to higher speed but also higher clustering distortion. Even though, when coping with high dimensional data in big size, the speed-up achieved by these methods is still limited.
Apart from above methods, there is another easy way to reduce the number of comparisons between samples and centroids, namely performing clustering in a top-down hierarchical manner [1, 28, 29]. Specifically, the clustering solution is obtained via a sequence of repeated bisections. The clustering complexity of k-means is reduced from to . This is particularly significant when n, d and k are all very large. In addition to that, another interesting idea from [1, 29] is that cluster centroids are updated incrementally [1, 23]. Moreover, the update process is explicitly driven by an objective function (called as criterion function in [1, 29]). Unfortunately, objective functions proposed in [1, 28, 29] are based on the assumption that input data are in unit length. The clustering method is solely based on Cosine distance, which makes the clustering results unpredictable when dealing with data in the general -space.
In this paper, a new objective function is derived directly from Eqn. 1, which makes it suitable for the whole -space. In other word, objective function proposed in  is the special case of our proposed form. Based on the proposed objective function, conventional egg-chicken k-means iteration is revised to a simpler form. On one hand, when applying the revised iteration procedure in direct k-way clustering, k-means is able to reach to considerably lower clustering distortion within only a few rounds. On the other hand, as the iteration procedure is undertaken in top-down hierarchical clustering manner (specifically bisecting), it shows faster speed while maintaining relatively lower clustering distortion in comparison to traditional k-means and most of its variants.
3 Clustering Objective Functions
In this section, the clustering objective functions upon which our k-means method is built are presented. Basically, two objective functions that aim to optimize the clustering results from different aspects are derived. Furthermore, we also show that these two objective functions can be reduced to a single form.
In order to facilitate the discussions that are followed, several variables are defined. Throughout the paper, the size of input data is given as n, while the number of clusters to be produced is given as k. The partition formed by a clustering method is represented as . Accordingly, the sizes of clusters are given as . The composite vector of a cluster is defined as . The cluster centroid 111We refer to as column vector across the paper. is defined by its members,
The inner-product of is given by , which is expanded as following form.
Re-arrange the above equation, we have
The sum of pairwise -distance within one cluster is given as
3.2 Objective Functions
In this section, two objective functions (also known as criterion functions ) are developed. In addition, with the support of the results obtained in Section 3.1, these objective functions will be reduced to simple forms, which enable them to be carried out efficiently in the incremental optimization procedure.
According to , objective functions are categorized into two groups. One group of the functions considers the tightness of clusters, while another focuses on alienating different clusters. In this paper, we focus on producing a clustering solution defined over the elements within each cluster. It therefore does not consider the relationship between the elements assigned to different clusters.
The first objective function we consider is to minimize the distance of each element to its cluster centroid, which is nothing more than the objective function of k-means.
The above equation is simplified as
Since the input data are fixed, is a constant. As a result, minimizing Eqn. 8 is equivalent to maximizing following function
Although objective function in Eqn. 9 is in the same form as the first objective function in , they are derived from different initial objectives. More importantly, in our case, there is no constraint that input data should be in unit length.
The second internal objective function that we will study minimizes the sum of the average pairwise distance between the elements assigned to each cluster, weighted according to the size of each cluster.
Plug Eqn. 6 in, we have
In Eqn. 11, is close to 1, the above criterion function can be approximated as
Noticed that similar optimization objectives have been discussed under Cosine similarity measure in . As it is shown that above two objective functions are the same in the paper. This is different from the result obtained in our case (general -space). As shown above, in -space, the objective functions for and are only approximately the same. The advantage that two objective functions are reduced to the same form is that, when we try to optimize one objective function, we optimize another in the mean time. Specifically, when we minimize the distances from elements to their cluster centroid, the average intra-cluster distance is minimized in the meantime. Since these two objective functions can be simplified to the same form, only objective function is discussed in the rest of paper.
Although objective function in Eqn. 9 is derived from Eqn. 1, the former is much easier to operate in the incremental k-means procedure. As it will be shown in the next section, it is quite convenient to evaluate whether Eqn. 9 attains a higher score (implies lower distortion in terms of Eqn. 1) when a sample is moved from one cluster to another.
4 K-means Driven by Objective Function
In this section, with the objective function developed in Section 3, two iterative clustering procedures are presented. Namely, one produces k clusters directly (called as direct k-way k-means), while another produces k clusters by bisecting input data sequentially k-1 times (called as bisecting k-means). Both clustering strategies are built upon incremental clustering [1, 23] and driven by objective function (Eqn. 9).
4.1 Clustering Algorithm
The basic idea of incremental clustering is that one sample is moved from cluster to as soon as this movement leads to higher score of objective function . To facilitate our discussion, the new function value as sample is moved from to is formulated as following.
In each iteration of the clustering, sample is randomly selected. The algorithm checks whether moving from its current cluster to any other cluster will lead to higher (i.e., ). If it is the case, is moved to another cluster. The clustering procedure is detailed in Alg. 1.
As seen from Step 3 of Alg. 1, the initialization of our method is different from most of the current practice of k-means, there is no assignment of each sample to its closest initial centroid. On the contrary, each sample is assigned with a random cluster label (ranges from 1 to k). This allows to calculate an initial score of and the composite vector D of each cluster. It is possible to do the initial assignment following the way of k-means or k-means++ . However, as will be revealed in Section 5, initialization under either k-means manner or k-means++ manner improves the clustering quality slightly. However, extra computation is required in such kind of initial assignment.
During each iteration, each sample is checked in random order. The optimization in Step 8-10 seeks the movement of that leads to highest increase of function score. From the optimization point of view, the algorithm reduces the clustering distortion greedily. From another point of view, the seeking process is comparable to the sample-to-centroid assignment in traditional k-means. They are actually on the same computational complexity level.
Whereas it is not necessary that we must seek the best movement for . As we discover by experiment, it is feasible that moving to another cluster as long as we find is greater than 0. On one hand, this will speed-up the iteration. On the other hand such kind of scheme usually takes more rounds to reach to the same level of distortion. However, we discover that such kind of less greedy scheme results in lower clustering distortion if the iteration loops for sufficient number of times.
Moving from one cluster to another (Step 9) is very convenient to take. It includes the operation that updates the cluster label of and the operation that updates the composite vector for cluster and , viz., , .
Note that this incremental updating scheme is essentially different from online learning vector quantization (LVQ) , in which the cluster centroids are updated incrementally. In the above iteration procedure, no cluster centroids are explicitly produced. As a result, there is no need to update cluster centroid. The clustering iteration is explicitly driven by an objective function rather than driven by the discrepancy between cluster centroids and their cluster members. As revealed later in the experiment, compared to LVQ, Alg 1 is more efficient and leads to considerably lower distortion.
Direct k-way Boost k-means
1: Input: matrix2: Output:3: Assign with a random cluster label;4: Calculate and ;5: while not convergence do6: for each (in random order) do7: Seek that maximizes ;8: if then9: Move from current cluster to ;10: end if11: end for12: end while
Figure 1 illustrates three iterations of Alg. 1 in 2D case. As shown in the figure, the initial clutsering result is random and messy. Samples belonging to different clusters are totally mixed up. However, only after one round of iteration, the clustering result becomes much more compact. The clustering terminates at the 10th round, where Lloyd’s condition is reached. The optimality of this procedure is analyzed in Appendix A and its convergence is proved in Appendix B.
Overall, method presented in Alg. 1 is different from traditional k-means in three major aspects. Firstly, no initial assignment is required. Moreover, the egg-chicken loop in the traditional k-means has been replaced by a simpler stochastic optimization procedure. Furthermore, unlike traditional k-means, it is not necessary to seek the best movement for each sample in the iteration.
The method presented in Alg. 1 is on the same complexity level as traditional k-means (i.e., ), which is unbearably slow when dealing with large-scale data.
The method is revised into a top-down hirarchical clustering version for large-scale clustering. Specifically, at each time, one intermediate cluster is selected and bisected into two smaller clusters by calling Alg. 1. The details of this method are given in Alg. 2.
As shown in Alg. 2, priority queue Q pops out one cluster for bisecting each time. As discussed in , there are basically two ways to organize the priority queue. One can prioritze the cluster with biggest size or the one with highest average intra-cluster distance to split. Similar as , we find splitting the biggest cluster usually demonstrates more stable performance. As a result, the queue is sorted in descending order by the cluster sizes in our practice.
Bisecting Boost k-means
1: Input: matrix2: Output:3: Treat X as one cluster ;4: Push into a priority queue Q;5: i = 1;6: while do7: Pop cluster from queue Q8: Call Alg. 1 to bisect into ;9: Push into queue Q;10: i = i + 1;11: end while
It is possible to partition the intermediate cluster into more than two clusters each time. In the following, we are going to show that this bisecting scheme achieves highest scalability among all alternative top-down secting schemes.
4.2 Scalability Analysis
In this section, the computation complexity of Alg. 2 is studied by considering the total number of comparisons required in the series of bisecting clustering. The number of iterations in each bisecting is assumed to be a constant by taking the average number of iterations.
In order to facilitate the analysis while without loss of generality, we assume that each intermediate cluster in Alg. 2 is partitioned evenly. In addition, we generalize Alg. 2 to an s-secting algorithm. Namely, an intermediate cluster is partitioned to s () clusters. Now we consider the size of series of intermediate clusters that are produced when performing sequential secting. Given q is the depth of splitting, it is easy to see . The sizes of all intermediate clusters are given as following.
As a result, the number of samples to be visited during the clustering procedure is
Considering that one sample has to compare with centroids each time, the total number of comparisons is
Given n and k are fixed, Eqn. 16 increases monotonically with respect to s. As a result, the number of comparisons reaches to the minimum (i.e., ) when . To this end, it is clear that bisecting is the most efficient secting scheme.
Compared with Alg. 1, the complexity of Alg. 2 is reduced to , where is the average number of iterations in each bisecting. Compared with t in traditional k-means, is much smaller given the scale of clustering problem is much smaller in terms of both the size of input data and the number of clusters to be produced. As a result, the complexity of Alg. 1 has been largely reduced since term has been multiplied by a much smaller factor .
Although Alg. 2 is efficient, the clustering result produced by Alg. 2 unfortunately does not satisfy with Lloyd’s condition. This problem is illustrated in Figure 2. As one of the clusters is further partitioned into two (from Figure 2(a) to Figure 2(b)), the partition over 2D space is formed by centroids changes. Cluster C claims bordering points from cluster B. However, points from cluster B cannot be reassigned to cluster C if no further intervention is involved. This is actually an underfitting issue and exists for any hierarchical clustering method. Fortunately, this issue can be alleviated by adopting Alg. 1 as a refinement procedure after Alg. 2 outputs k clusters. To do so, extra time is required. It therefore becomes a problem of balancing between efficiency and quality.
According to our observation, it is possible to further speed-up the proposed boost k-means. After a few iterations, both k-means and boost k-means will be trapped in a local minima. Only samples that bordering between different clusters should be shuffled from one cluster to another. As a result, given a sample, it is no need to search for the best movement among k clusters. Instead, sample only needs to compare with top- () centroids to search the suitable movement. We find that, this simple modification typically leads to times speed-up while without significant performance degradation.
|Probability based ||k-means++ ||BsBKM++||BKM(kpp)||BsBKM(kpp)|
In this section, the effectiveness of proposed clustering method, namely boost k-means (BKM) is studied under different scenarios. In the first experiment, dataset SIFT1M  is adopted to evaluate the clustering quality. In the second experiment, BKM is tested on the nearest neighbor search task based on product quantizer (PQ)  in which this method is adopted for quantizer training. In the third experiment, BKM has been applied to traditional document clustering. Following the practice of [1, 29], 15 document datasets222Available at http://glaros.dtc.umn.edu/gkhome/fetch/sw/cluto/datasets.tar.gz have been adopted. In the last experiment, the scalability of BKM has been tested on large-scale image clustering task, for which the number of images we use is as large as 10 million.
In our study, the performance from traditional k-means is treated as comparison baseline. In addition, representative k-means variants, such as Mini-Batch , Repeated Bisecting k-means (RBK) , online Learning Vector Quantization (LVQ)  and k-means++  are considered in the comparison. For Mini-Batch, our configuration makes sure that the iteration covers 10% of the input data. The configuration is fixed across all the experiments. For RBK, we select the objective function that maximizes the average Cosine similarity between samples within one cluster, which is the special case of ours given the input data is -normalized. LVQ is similar to k-means except that in each round, a cluster centroid is upated as soon as a sample is assigned. The updating rate starts from 0.01 and decreases at a pace of in one iteration.
As shown in Table I
, there are variants of k-means depending on cluster initialization and data partitioning strategies (e.g., direct k-way or bisecting). This is also true for the proposed BKM. In the table, ‘initial assignment’ refers to the operation of assigning each sample to its closest initial centroid. When the initial assignment is based on random seeding like traditional k-means, it is denoted as ‘rnd’. When it is based on probability distribution seeding as k-means++, it is denoted as ‘kpp’. Initialization without initial assignment is denoted as ‘non’. In the experiments, all the variants out of these different configurations on k-means as well as BKM are considered. Performance evaluation is separately conducted for k-way and bisecting clustering method. Noted that BsBKM(rnd) is the same as RBK if the input data is-normalized. The experiment in this section is conducted using 1 million SIFT features . The features are clustered into 10,000 partitions and the average distortion error is calculated for performance evaluation.
In addition, we also study the performance trend of BKM when Steps 7-10 in Alg. 1 are modified to moving the sample as soon as . The variants under this modification are denoted as BKM(xxx)+Fast. 333Note that this is not applicable for bisecting BKM. All the methods considered in the paper are implemented in C++ and the simulations are conducted on a PC with 2.4GHz Xeon CPU and 32G memory setup.
5.1 Evaluation of Clustering Distortion
Since k-means and most of its variants share the same objective function (Eqn. 1), it is straightforward to evaluate the clustering performance by checking to what degree the objective is reached. The average distortion (given in Eqn. 17) is adopted for evaluation , which takes average over Eqn. 1,
For above equation, the lower the distortion value, the better quality of the clustering result is.
The first experiment mainly studies the behavior of the proposed BKM under different initializations. The average distortion curves produced by variants direct k-way BKM are given in Figure 3(a) as a function of numbers of iteration. Traditional k-means is treated as baseline for performance comparison. The result shows that clustering distortion of BKM drops faster than traditional k-means. The average distortion from traditional k-means is around 40,450 after 130 iterations. In contrast, BKM without initial assignment (BKM(non)) is able to reach to the same distortion level after only 7 iterations. Moreover, we find that initializing BKM as traditional k-means way (BKM(rnd)) or as k-means++ (BKM(kpp)) allows the iteration to start from a low distortion level. Nevertheless the advantage over BKM(non) fades away after 15 iterations. In comparison to BKM(non), the extra cost of adopting initial assignment in BKM is relatively high.
The second experiment studies the performance trend of Alg. 1 in case when Steps 7-10 do not seek the best movement (BKM(xxx)+Fast). As shown in Figure 3(b), the distortion drops slower than BKM(non) which seeks the best movement. However, lower distortion is achievable by BKM(rnd)+Fast when reaching to sufficiently number of iterations (e.g., 20 iterations). This indicates that when the optimization scheme is more greedy, it is likely to get trapped in a worse local optima. This observation applies to BKM under different kinds of initialization. Noted that the time cost for BKM(xxx)+Fast is lower than that of BKM that seeks the best movement in each iteration. Whereas, BKM(xxx)+Fast usually needs a few more number of iterations to reach to the similar distortion level. Overall, as investigated in Section 5.4, BKM(xxx)+Fast is faster than BKM(xxx).
Figure 3(c) studies the trend of average distortion among the proposed BKM (specifically BKM(non)), traditional k-means, k-means++, Mini-Batch and LVQ. For all the methods presented, their distortion decreases steadily as the iteration continues. A big performance gap is observed between Mini-Batch and other k-means variants. In addition k-means and k-means++ share similar distortion curve. BKM(non) outperforms k-means and k-means++ by requiring only 7 iterations. Most of the methods including k-means and k-means++ take more than 120 iterations to finally converge. On the other hand, little distortion is observed after 20 iterations, which implies the possibility of terminating the iteration at 20. Although similiar as BKM, LVQ updates the intermediate clusters incrementally, updating cluster centroid directly turns out to be inefficient, which leads to considerably poor performance.
Since k-means and its variants are all sensitive to initialization, the performance fluctuates from one run to another. The candelstick chart shown in Figure 3(d) further confirms the significance of the improvement achieved by BKM. This chart is plotted with 128 clustering runs () on SIFT100K  for each method. As shown in the figure, although the performance flucturates for all the methods, the variations are minor. Similar as previous observation, there is no significant difference between traditional k-means and k-means++. In contrast, the performance gap between BKM and traditional k-means is much more significant than the performance variations across different runs.
Table II shows the average distortion of different k-means variants under bisecting strategy. The result from k-means (after 130 iterations) is presented for the comparison. As shown from the table, the average distortion from all bisecting methods are on the level of . Methods built upon Alg. 1 always perform better. The average distortion from all bisecting clustering methods are much higher than that of k-means. They are actually only close to the distortion level of k-means after one iteration. However, the merit of clustering with bisecting strategy is that it is more than 20 times faster than k-means of a single iteration. The relatively poor clustering quality produced by bisecting strategy is mainly due to the issue of underfitting (as discussed in Section 4.2). The clustering results can be further refined by Alg. 1 as shown on the 3rd row of Table II.
As learned from above experiments, on one hand initial assignment under k-means manner or under k-means++ manner is able to improve the performance of BKM slightly. On the other hand, the initial assignment slows down the method considerably. A trade-off has to be made. In the following experiments, only the results from two representative configurations of BKM, namely BKM(non) and BKM(rnd)+Fast are presented. We leave the other possible configurations to the readers.
5.2 Nearest Neighbor Search by Product Quantizer (PQ)
In this section, BKM is applied for visual vocabulary training using product quantization . Following the practice of , 100K SIFT features are used for product quantizer training, while SIFT1M set  is encoded with the trained product quantizers as the reference set for nearest neighbor search (NNS). The obtained recall@top-k is averaged over 1,000 queries for each method. In the experiment, two different settings are tested for product quantizer. Namely, the 128-dimensional SIFT vector is encoded with 8 and 16 product quantizers respectively. For clarity, the evaluations are separately conducted for direct k-way and bisecting k-means.
Recall@top-100 for direct k-way are presented in Figure 4(a)-(d) under two different settings ( and ). As seen from the figures, the performances from k-means, k-means++ and BKM(non) are all very close to each other under different settings. The product quantizer trained with bisecting clustering methods shows only 0.1-1.3% lower performance than that of direct k-way methods.
This basically indicates that product quantizer itself is tolerant to clustering quality. The performance of Mini-Batch and RBK is around 2-6% lower than the other methods. The poor performance of RBK basically indicates the optimization objective function under Cosine similarity is not directly feasible for general -space.
5.3 Document Clustering
In this section, the performance of proposed method is evaluated under the context of document clustering. Following in , 15 document datasets are used for evaluation. The documents has been represented with TF/IDF model and normalized to unit length. Similar to , entropy is adopted for the evaluation as following
where c is the number of classes. Eqn. 18 evaluates to what degree that elements from the same class are put in one cluster. The lower the value, the better the result is. In the experiment, each method performs clustering for 10 runs, and the run with the lowest entropy is presented in Table III. The presented entropy are averaged over 15 datasets.
|k = 5||k = 10||k = 15||k = 20|
In general, methods based on BKM perform considerably better. Furthermore, methods with bisecting strategy demonstrate slightly better performance than that of direct k-way in the document clustering task, which shares similar observation as . Overall, BsBKM(non) shows the best performance. While the performance of RBK is close to BsBKM(non). These two methods are quite similar except that no initial assignment is involved in BsBKM(non). This indicates the advantage of no initial assignment in this scenario, which allows clustering to converge to a better local optima.
5.4 Scalability Test on Image Clustering
In this section, the scalability of the proposed k-means is tested on image clustering. The experiment is conducted on 10 million Flickr images (Flickr10M), which are a subset of YFCC100M . Hessian-Affine  keypoints are extracted from each image and are described by RootSIFT feature . Finally, the RootSIFT features from each image are pooled by VLAD  with a small visual vocabulary of size 64. The resulting 8,192-dimensional feature is further mapped to 512 dimensions by PCA. Following , the final VLAD vector is normalized to unit length. In the direct k-way clustering case, we set the number of maximum iterations for all methods to 20. While for the bisecting case, there is no threshold on the number of iterations. The results reported in this section have been averaged over 10 runs for each method.
In the first experiment, clustering methods are tested in the way that the scale of input images varies from 10K to 10M. A fixed number of clusters, i.e., 1,024 is used regardless of the size of dataset. The time costs for direct k-way and bisecting methods are presented in Figure 5(a)-(b). Accordingly, the average distortion of all the methods are presented in Figure 6(a).
As shown in the figures, BKM exhibits slightly faster speed over k-means and its variants across different scales of input data under both direct k-way and bisecting cases. The speed-up becomes more significant as the scale of input data increases. The higher efficiency of these methods is mainly attributed to no requirement of initial assignment. Compared with BKM(non), BKM(rnd)+Fast takes extra time. However, the cost of initial assigment is compensated later for not seeking the best movement.
Compared with direct k-way clustering, methods with bisecting strategy achieve much higher scalability. In particular, BsBKM(non) shows the highest scalability. It only takes less than 94 minutes to cluster 10 million vectors (in 512 dimensions) into 1,024 clusters. The efficiency of Mini-Batch is close to BsBKM(non). However, as shown in Figure 6(a), the clustering quality is poor in most of the cases. Overall, BKM(rnd)+Fast achieves the highest speed efficiency and lowest distortion among all direct k-way clustering methods. While in the bisecting case, BsBKM(non) shows the best performance in terms of both speed efficiency and clustering quality. Similar to the experiments in Section 5.1, the average distortion introduced by bisecting clustering is much higher than direct k-way due to the problem of under-fitting.
In addition, the scalability of clustering methods is tested in the way that the number of clusters by varying from 1,024 to 8,192, while the scale of input data is fixed to 1 million. Figure 5(c)-(d) show the time cost of all 9 methods. Accordingly, the average distortion from all these 9 methods are presented in Figure 6(b). As shown in the figures, for all direct k-way clustering methods, the time cost increases linearly as the number of clusters increases. Mini-Batch is no longer efficient as k increases. In contrast, the time cost of all bisecting methods remains steady across different cluster numbers. In terms of clustering quality, as seen from Figure 6(b), in both direct k-way and bisecting cases, clustering driven by the proposed optimization procedure (Alg. 1) performs considerably better. A clear trend is observed from Figure 6(b), methods based on Alg. 1 shows increasingly higher performance than the rest as k grows. Overall, clustering driven by proposed optimization process shows higher speed and better quality. The highest speed is achieved by BsBKM(non), for which only 8 minutes are required to cluster 1 million high dimensional data into 8,192 clusters. Due to extra cost in initial assignment, bisecting with traditional k-means and k-means++ still shows around 35% slower speed than BsBKM(non).
As a summary, clustering based on Alg. 1 shows superior performance in terms of both speed efficiency and quality under different scenarios. This is mainly due to the nature of incremental updating scheme, which allows the cluster structures to be fine-tuned in a more efficient way. When the proposed Alg. 1 is performed under bisecting manner (i.e., BsBKM(non)), it shows two orders of magnitude faster than traditional k-means.
We have presented a novel k-means variant. Firstly, a clustering objective function that is feasible for the whole -space is developed. Supported by the objective function, the traditional k-means clustering has been modified to simpler form. In this novel k-means variant, we interestingly find that neither the costly initial assignment nor the seeking of closest centroid for each sample in the iteration are necessary. This leads to higher speed and considerably lower clustering distortion. Furthermore when the proposed clustering method is undertaken in the ways of top-down bisecting, it achieves the highest scalability and best quality among all hierarchical k-means variants. Extensive experiments have been conducted in different contexts and on various datasets. Superior performance over most of the k-means variants is observed across different scenarios.
This work is supported by National Natural Science Foundation of China under grants 61572408. The authors would like to express their sincere thanks to Prof. George Karypis from University of Minnesota, USA for his detailed explanation about the implementation of repeated bisecting k-means.
Appendix A: Optimality of Incremental K-means
As shown in Eqn. 13 and Eqn. 9, two optimal objectives are quite similar. In this section, we show that optimal solution with respect to objective function (Eqn. 9) can be reached with incremental updating scheme presented in Alg. 1.
: For contradiction, let be an optimal solution and assume that there exists one element d and clusters and such that . Now consider the clustering solution . Let , , and , be the composite and centroid vectors of cluster and , respectively. Let , then
Let’s define , are the average pairwise inner product in cluster and respectively. In addition, and are given as the average inner-products between d and elements in and respectively, viz , and . Above Equation is rewritten as
Given the fact that , we have , which is contradicting.
Appendix B: Convergence of Incremental k-means
and are two clusters. d is initially part of , and is the composite of exclude d, is the centroid of exclude d, is the composite and centroid of cluster , the move condition of d from to should satisfied
This equation can be rewritten as:
Now if we assume that both and are sufficiently large, then and will be close to 1. Under these assumptions, we can get
Now , are defined as the average pairwise inner product in cluster and respectively. and are given as the average inner-products between d and elements in and respectively, viz , and , the following inequation holds.
-  Y. Zhao and G. Karypis, “Empirical and theoretical comparisons of selected criterion functions for document clustering,” Machine Learning, vol. 55, pp. 311–331, 2004.
-  Y. Avrithis, Y. Kalantidis, E. Anagnostopoulos, and I. Z. Emiris, “Web-scale image clustering revisited,” in ICCV, pp. 1502–1510, 2015.
-  J. Shi and J. Malik, “Normalized cuts and image segmentation,” Trans. PAMI, vol. 22, pp. 888–905, Aug. 2000.
-  J. Sivic and A. Zisserman, “Video Google: A text retrieval approach to object matching in videos,” in ICCV, pp. 1470–1477, Oct. 2003.
-  H. Jégou, M. Douze, and C. Schmid, “Product quantization for nearest neighbor search,” Trans. PAMI, vol. 33, pp. 117–128, Jan. 2011.
-  M. Muja and D. G. Lowe, “Scalable nearest neighbor algorithms for high dimensional data,” Pattern Analysis and Machine Intelligence, IEEE Transactions on, vol. 36, pp. 2227–2240, 2014.
-  A. Babenko and V. Lempitsky, “Additive quantization for extreme vector compression,” in CVPR, pp. 931–938, 2014.
-  S. P. Lloyd, “Least squares quantization in PCM,” IEEE Trans. Information Theory, vol. 28, pp. 129–137, 1982.
-  X. Wu, V. Kumar, J. R. Quinlan, J. Ghosh, Q. Yang, H. Motoda, G. J. McLachlan, A. Ng, B. Liu, P. S. Yu, Z.-H. Zhou, M. Steinbach, D. J. Hand, and D. Steinberg, “Top 10 algorithms in data mining,” Knowledge and Information System, vol. 14, pp. 1–37, Dec. 2007.
-  N. Ailon, R. Jaiswal, and C. Monteleoni, “Streaming k-means approximation,” NIPS, pp. 10–18, 2009.
-  A. Vattani, “k-means requires exponentially many iterations even in the plane,” Discrete and Computational Geometry, vol. 45, no. 4, pp. 596–616, 2011.
-  D. Arthur and S. Vassilvitskii, “K-means++: The advantages of careful seeding,” in Proceedings of the Eighteenth Annual ACM-SIAM Symposium on Discrete Algorithms, pp. 1027–1035, 2007.
-  M. Shindler, A. Wong, and A. W. Meyerson, “Fast and accurate k-means for large datasets,” in NIPS, pp. 2375–2383, 2011.
-  T. Kanungo, D. M. Mount, N. S. Neanyahu, C. D. Piatko, R. Silverman, and A. Y. Wu, “An efficient k-means clustering algorithm: Analysis and implementation,” Trans. PAMI, vol. 24, pp. 881–892, Jul. 2002.
-  C. Elkan, “Using the triangle inequality to accelerate,” in ICML, 2003.
-  D. Sculley, “Web-scale k-means clustering,” in In Proceedings of the 19th international conference on World wide web, pp. 1177–1178, 2010.
-  B. Bahmani, B. Moseley, A. Vattani, R. Kumar, and S. Vassilvitskii, “Scalable k-means++,” In Proceedings of the VLDB Endowment, vol. 5, no. 7, pp. 622–633, 2012.
-  A. Broder, L. Garcia-Pueyo, V. Josifovski, S. Vassilvitskii, and S. Venkatesan, “Scalable k-means by ranked retrieval,” in Proceedings of the 7th ACM international conference on Web search and data mining, pp. 233–242, 2014.
-  Y. Gong, M. Pawlowski, F. Yang, L. Brandy, L. Boundev, and R. Fergus, “Web scale photo hash clustering on a single machine,” in CVPR, pp. 19–27, 2015.
-  L. Bottou and Y. Bengio, “Convergence properties of the kmeans algorithm,” Advances in Neural Information Processing Systems, pp. 585–592, 1995.
-  M. Ester, H. peter Kriegel, J. Sander, and X. Xu, “A density-based algorithm for discovering clusters in large spatial databases with noise,” in In Proceedings of Knowledge Discovery and Data Mining, pp. 226–231, 1996.
-  D. Comaniciu and P. Mee, “Mean shift: A robust approach toward feature space analysis,” Trans. PAMI, vol. 24, pp. 603–619, May 2002.
-  R. O. Duda, P. E. Hart, and D. G. Stork, Pattern Classification. Wiley-Interscience, 2 ed., 2001.
-  A. K. Jain, M. N. Murty, and P. J. Flynn, “Data clustering: A review,” ACM Computer Survey, vol. 31, pp. 264–323, Sep. 1999.
-  J. L. Bentley, “Multidimensional binary search trees used for associative searching,” Communications of ACM, vol. 18, pp. 509–517, Sep. 1975.
-  D. Pelleg and A. Moore, “Accelerating exact k-means algorithms with geometric reasoning,” in Proceedings of the Fifth ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, pp. 277–281, ACM, 1999.
-  A. Goswami, R. Jin, and G. Agrawal, “Fast and exact out-of-core k-means clustering,” in Proceedings of the Fourth IEEE International Conference on Data Mining, pp. 83–90, 2004.
-  A. K. Jain and R. C. Dubes, Algorithms for Clustering Data. Prentice-Hall, Inc., 1988.
-  Y. Zhao and G. Karypis, “Hierarchical clustering algorithms for document datasets,” Data Mining and Knowledge Discovery, vol. 10, no. 2, pp. 141–168, 2005.
-  T. Kohonen, M. R. Schroeder, and T. S. Huang, eds., Self-Organizing Maps. Secaucus, NJ, USA: Springer-Verlag New York, Inc., 3rd ed., 2001.
-  D. Lowe, “Distinctive image features from scale-invariant keypoints,” IJCV, vol. 60, pp. 91–110, Nov. 2004.
-  B. Thomee, D. A. Shamma, G. Friedland, B. Elizalde, K. Ni, D. Poland, D. Borth, and L.-J. Li, “YFCC100M: The new data in multimedia research,” Communications of the ACM, vol. 59, no. 2, pp. 64–73, 2016.
-  K. Mikolajczyk and C. Schmid, “Scale and affine invariant interest point detectors,” IJCV, vol. 60, pp. 63–86, Oct. 2004.
-  R. Arandjelovic and A. Zisserman, “Three things everyone should know to improve object retrieval,” in CVPR, pp. 2911–2918, Jun. 2012.
-  H. Jégou, F. Perronnin, M. Douze, J. Sánchez, P. Pérez, and C. Schmid, “Aggregating local descriptors into compact codes,” Trans. PAMI, vol. 34, pp. 1704–1716, Sep. 2012.