Optimal clustering by merge-based branch-and-bound

1.   Introduction

Clustering is a fundamental problem that must often be solved as a part of more complicated tasks in pattern recognition, image analysis and other fields of science and engineering ^[1,2,3,4]. Clustering is also needed for designing a codebook in vector quantization ^[4]. The clustering problem is defined here as follows. Given a set of N data points X = {x₁, x₂, …, x_N}, partition the data set into M clusters so that a given objective function f is minimized.

The clustering problem in its combinatorial form has been shown to be NP-hard ^[5]. No polynomial time algorithm is known to find the globally optimal solution, and sub-optimal solutions are usually obtained by heuristic algorithms. Despite the known limitations implicated by the NP-completeness, solving the optimal clustering problem has theoretical interest that can provide insight to the problem itself. It might also have practical implications in the case of problem instances of limited size.

Agglomerative clustering is an approach for generating the clustering hierarchically. The clustering starts by initializing each data point as its own cluster. Two clusters are merged at each step and the process is repeated until the desired number of clusters is obtained. Ward's method ^[6] is a method that selects the cluster pair to be merged so that it increases mean square error least. In the vector quantization context, this is known as the pairwise nearest neighbor (PNN) method due to ^[7]. In the rest of this paper, we denote it as the PNN method.

The PNN method is interesting here because of its conceptual simplicity and because of the optimality of the single merge step. This step reduces a given clustering from m clusters to m-1 clusters by minimizing the optimization function value. Even though the step is optimal, there is no guarantee of optimality of the final clustering resulting from a series of locally optimal merge steps. The main idea of the PNN method, however, can be generalized so that we do not optimize only a single merge but over multiple merge steps.

In this paper, we present an optimal clustering algorithm derived from the PNN method. It is easy to see that any clustering can be produced by a series of merge operations. Every merge reduces the number of clusters by one. It therefore takes exactly N-M steps to generate a clustering with M groups from the set of N points. Optimal clustering can be found by considering all the possible merge sequences and finding the one that minimizes the objective function. The idea can be implemented as a branch-and-bound technique that uses a search tree for finding the optimal clustering, and a suitable bounding criterion to cut out non-optimal branches of the tree.

The relation of the proposed method to the PNN method is demonstrated in Fig. 1. At the first step, all possible merges are generated. At the second step, the PNN method would continue from the locally optimal result whereas branch-and-bound technique will study all branches. The root of the search tree represents the case where all data points are assigned to their own clusters. At the level N-m, there are all possible solutions to m clusters. The final clustering with M groups is located at the level N-M. All branches of the tree must be generated to find the optimal clustering.

We consider also two sub-optimal variants that compromise the optimality but work in polynomial time. The first algorithm generates the search tree only down to a fixed limit. The best result of the fixed depth is then taken as the new starting point, and the same procedure is repeated until the desired number of clusters is reached. The second variant proceeds only one level in the tree along the path towards to the direction of the best solution within the search range. After this, a completely new search tree is generated for one level further, and the process is then repeated.

The rest of the paper is organized as follows. In Section 2, we give formal description of the clustering problem, and review previous literature of optimal clustering. In Section 3 we recall the PNN method. The branch-and-bound technique is introduced in Section 4. We first study the redundancy of the search tree in Section 4.1, and then present a method to generate minimum redundancy tree in Sections 4.2 and 4.3. A bounding criterion is formalized in Section 4.4 to cut irrelevant branches of the tree. Two suboptimal polynomial time variants are introduced in Sections 5. Experimental tests are made in Section 6, and conclusions drawn in Section 7. The proposed branch-and-bound method was first reported in a conference ^[8], and the polynomial time variants later in ^[9].

2.   Optimal clustering

Given a set of N data points X = {x₁, x₂, …, x_N}, clustering aims at solving the partition P = {p₁, p₂, …, p_N}, which defines for each data point the index of the cluster where it belongs to. Cluster s_a is defined as the set of data points that belong to the same partition a:

Clustering is represented as the set S = {s₁, s₂, ..., s_M}, and a set of cluster centroids C = {c₁, c₂, …, c_M}.

The most important choice in clustering is the objective function f for evaluating the goodness of clustering. When the data points belong to the Euclidean vector space, a commonly used function is the mean square error between the data points and their cluster centroids. Given a partition P and the cluster centroid C, it is calculated as:

The choice of the function depends on the application and there is no general solution of which measure should be used. However, once the objective function is decided the clustering problem can be formulated as a combinatorial optimization problem.

Optimal solution can be solved by constructing all M^N/M! possible groupings of N data points into M groups and selecting the optimal one. This can be implemented by brute force by permuting all possible partitions of the data points. We refer this as partition-based approach.

This approach has been used for developing branch-and-bound technique by Koontz et al. ^[10]. They construct the partition by assigning data points to the clusters one by one, and by using a partial error for bounding non-optimal solutions. They have reported to solve problems of size N = 40, and of size N = 120 by compromising the optimality with a stronger bouncing criterion.

Cheng ^[11] finds clusters in a binary matrix where the values correspond to files/transactions relationships, and the goal is to organize the database for minimizing disc access. Problem size of 45 × 20 matrix has been solved by the method.

In ^[12], the clustering problem was formulated as a partitioning on edge-weighted graph, in which the goal is to minimize the sum of weights of the edges within the clusters ^[12]. Problem size of N = 145 has been considered.

Iyer and Aronson ^[13] proposed a parallel implementation with a linear speed-up with respect to the number of processors available. This increases the size of problems that is possible to solve by the algorithm, but only by a logarithmic factor.

Brusco ^[14] added data points to the solution one by one until all points have been added. Each point has M possible clusters to join, and branch-and-bound technique is used to prune the search tree. They managed to provide optimal solutions for problems with up to N = 240 objects in case of M = 8 well-separated clusters. For randomly distributed data optimal solution was found up to N = 60.

Since then, optimality of clustering has received very little attention. Guns et al. extended Brusco's approach using constrained programming ^[15]. Optimal clustering was studied in cutting and packing problem in 2-d space ^[16]. Optimal graph structure was presented in ^[17] for graph-based clustering. Missing-value problem was shown to fit nicely into the optimal clustering framework in ^[18]. All these results have mainly theoretical interest and do not generalize beyond small-scale instances.

Approximation algorithms have also been considered but with a limited success. Feder and Greene have shown that one cannot approximate optimal cluster size for fixed number of clusters in polynomial time size within a factor close to 2, unless P = NP ^[19]. A solution with the time complexity of O(N log M) was then proposed. Mettu and Plaxton have proposed a randomized O(1)-approximation algorithm that works in O(NM) time ^[20].

Polynomial time solutions are known but only for some special cases. For example, the clustering has been formulated as a graph problem in ^[21] using the assumption that the data set can be adequately represented by an adjacency graph. The method finds the partition that minimizes the maximum flow between the sub graphs (clusters) in the given partition. Efficient polynomial time solutions are known for this problem. The clustering problem involved in parallel scheduling ^[22] has been formulated as finding exact minimum makespan constrained by a tree structure.

It is also noted that efficient algorithm exists for the 1-dimensional special case. An example is the scalar quantization problem, for which O(NM) solution is known ^[23].

3.   Pairwise nearest neighbor method

The pairwise nearest neighbor (PNN) method ^[6,7] generates the clustering hierarchically using a sequence of merge operations as described in Fig. 2. In each step of the algorithm, the number of the clusters is reduced by merging two nearby clusters:

The cost of merging two clusters s_a and s_b is the increase in the MSE-value caused by the merge. It can be calculated using the following formula ^[7]:

where n_a and n_b denote to the sizes of the corresponding clusters s_a and s_b.

PNN applies local optimization strategy: all possible cluster pairs are considered and the one increasing the objective function least is chosen:

A single merge step of PNN is optimal but there is no guarantee of optimality of the final clustering resulting from a series of locally optimal merge steps. The time complexity of PNN varies from O(N²) to O(N³) depending on the implementation and data set ^[24].

4.   Merge-based branch-and-bound technique

We described next a branch-and-bound technique that generates optimal clustering by a sequence of merge operations. It is easy to see that any clustering can be produced by merging the data points into the groups one by one. Every merge operation reduces the number of clusters by one. It therefore takes exactly N-M steps to generate a clustering with M clusters, independent of the order of the merge operations.

For example, consider the example shown in Fig. 1, in which we have five data points {A, B, C, D, E}. The resulting clustering can be generated by the following three merge operations:

Initial: {A} {B} {C} {D} {E}

Step 1: {AE} {B} {C} {D}

Step 2: {AE} {BC} {D}

Step 3: {AE} {BCD}

All alternative merge sequences can be represented as a search tree. The root of the tree represents the starting point in which every data point is assigned to its own cluster (N clusters), and its descendants represent all possible clusterings of N-1 clusters. In general, every node in the tree represents a single clustering with m clusters, and its children represent the clusterings that have been produced by merging any two of m existing clusters.

4.1.   Redundancy of the search tree

The search tree includes a lot of redundant solutions (referred to redundant search tree) as the same partition can be constructed with different orders of the merge operations. The partitions are generated from the search tree as follows. At the first step, there are N·(N-1)/2 alternatives for the merge operation, and therefore, equally many partitions at the level (m = N-1). Similarly, at the second level (m = N-2) there are (N-1)·(N-2)/2 alternatives for the merge operation, and they are independent of the merge operation made at the previous level.

In general, every node has (m)·(m-1)/2 children at the level with m clusters. We can therefore derive the total number of merge sequences, which lead to a partition of M clusters as follows:

At the same time, we know from ^[25] that the total number of different partitions equals to Stirling's number of second kind ^[26]:

We can calculate the average number of sequences per different partitions as:

where N and M are nonnegative, and N$\geqslant$M. In other words, the search tree contains significant amount of redundant solutions.

4.2.   Permutating non-redundant clusters

We consider next a single cluster represented as a list of the data points, and merge operation as the catenation of the two lists. For example, the clustering in Fig. 1 is represented as the pair of lists (AE) (BCD), and their merge as (AEBCD). Using this representation, we can see that the same cluster has several different representations. The cluster (BCD), for example, has the following representations:

(BCD)

(BDC)

(CBD)

(CDB)

(DBC)

(DCB)

The data points x_i can be ordered by their index i in the data set. We therefore use the following condition to prevent redundant representations of the same cluster.

Condition 1: The only valid representation for a cluster s_j = {x₁, x₂, ..., x_nj} is the ordered sequence

(x₁ x₂ ... x_nj).

For example, the only valid representation for {BCD} in the previous example is then (BCD). When using this condition, we can still represent all possible clusterings but without redundant representations for a single cluster.

The condition 1 can be applied in the Branch-and-bound algorithm using the following merge condition.

Condition 2: Clusters s_a and s_b can be merged iff: $\;\;i < j\;\;\forall \;\;{x_i} \in {s_a}, {x_j} \in {s_b}$.

Because of this merge rule, the order of the data points will be automatically retained. Yet, every possible cluster can be obtained by merging the data points starting to merge from the smallest index one by one.

For example, the cluster (ABE) is possible to obtain but only using the sequence that merges first (A) + (B), and then (AB) + (E). On the other hand, the cluster pair (AE) (B) cannot be merged because the resulting cluster (AE) + (B) = (AEB) is not a valid representation as the data points are not sorted. Furthermore, if the current clustering were (AE) (B) (C) (D), we could not merge the cluster (AE) with any other cluster anymore and (AE) would inevitably remain as such in the final clustering.

4.3.   Permutating minimum redundancy search tree

The condition 2 removes the redundancy in the case of representing a single cluster but it is still possible to construct the same cluster via different paths (merge sequences) in the search tree. For example, the cluster (BCD) can be constructed using two different paths:

Sequence 1: (B) (C) (D) → (BC) (D) → (BCD)

Sequence 2: (B) (C) (D) → (B) (CD) → (BCD)

Furthermore, the clustering (AE) (BCD) can be reached by six different merge sequences, of which two are shown in Table 1.

It is therefore not enough to limit only the intra cluster representation, but we must also limit the permutation of the search paths in the tree. We do this by applying the following permutation criterion:

Condition 3: Clusters s_a and s_b can be merged iff $\;\;a \geqslant {a_0}\; \wedge \;b > a$,

where a₀ is the index of the first cluster in the previous merge. In other words, we force the algorithm to permutate the cluster pairs in a predefined order so that the index of the first cluster is always monotonically non-decreasing during the process. So, if we have merged clusters s_a0 and s_b0 at the previous level, we can consider only cluster pairs s_a and s_b such that a ≥ a₀. The second term (b > a) is induced by condition 2.

We can see that any individual cluster s_j = (s_j1 s_j2 s_j3... s_jn) is constructed by the following sequence of merge operations: (s_j1) + (s_j2) → (s_j1 s_j2) + (s_j3) → (s_j1 s_j2 s_j3) + (s_j4), and so on. This fulfills the constraint b > a. Any clustering {s₁, s₂, ..., s_m} can then be generated by constructing the clusters one by one in the order from s₁ to s_m. This fulfills the constraint a ≥ a₀.

Furthermore, there is no other merge sequence that could construct the same clustering without contradicting the conditions 2 and 3. These conditions together guarantees that every cluster can be constructed in only one order, and the condition 2 (a ≥ a₀) that the clusters are constructed in a unique sequence from smallest index to largest. Thus, the use of the conditions 2 and 3 produces non-redundant search tree.

In practice, the condition 2 can be too complicated to be implemented in practice. We therefore introduce the following condition that simplifies it.

Condition 4: Clusters s_a and s_b can be merged iff a = a₀ $\Rightarrow$ b ≥ b₀.

This implication says that if we merge the same cluster s_a as previously, the index of the second cluster must be greater than that of the previously merged cluster s_b0. In other cases (a > a₀), the inequality b > a from condition 3 is sufficient to satisfy also the condition 2.

The non-redundant search tree for the previous example is illustrated in Fig. 3. At the first level, the permutation creates the merges: (AB) (AC) (AD) (AE) (BC) (BD) (BE). We can see that after the merge (B)+(C), the merge (A)+(C) do not appear anymore because of condition 3, and it already exists in the branch where (AC) was constructed before (BD). Furthermore, if the previous merges created clusters (AC) and (BD), the cluster (ACBD) does not appear anymore as it has already been created by the sequence (A)+(B), (AB)+(C), (ABC)+(D).

The condition 3 removes redundant clusterings but there still exist partial branches that cannot be completed (resulting in too few clusters). For example, after the merge sequence (AB) and (CE), there would be no more valid merges left because the permutation criterion does not allow us to add new point in the cluster (AB) anymore, and because the merge (CE)+(D) would break the intra cluster order. Such branches can be eliminated using rather simple bounds for the permutation loop.

When we permutate new cluster pairs for merge, we always start to permute from the previous cluster s_a0 and consider all potential pairs in a loop. Clusters pairs (s_a, s_b) to be considered are such that 1 ≤ a < b ≤ m, where m is the current number of clusters. The index a also indicates how many clusters have been completed by now. This is because they are not allowed to be included in the merge operation any more due to condition 3.

The same applies also to the clusters that have index between a and b. Since we have m-M more merges to be performed, we know that there must be equally many valid cluster pairs. Concluding from this, we can derive the upper bound for the index a at any stage of the process.

Condition 5: The first cluster index in the merge is upper limited to a ≤ M.

Greater values than this would lead to situation that we cannot complete the clustering with M clusters. Using the conditions 1, 3, 4 and 5, we can now present the algorithm for generating non-redundant search tree as the pseudo code shown in Fig. 4.

4.4.   Bounding criterion

The search can be terminated earlier when we know that the current branch cannot lead to a better solution than the best solution found so far. The termination is based on the fact that every merge operation increases the MSE-value of the solution, see Fig. 5. Thus, when we have generated the first solution in the search tree, we can use its MSE-value as the upper bound for the optimal solution. Other branches of the tree can then be terminated if the following condition is true:

Bounding criterion 1: $MS{E_t} \geqslant MS{E_{\min }}$,

where MSE_t is the value of the current solution after the t^th merge, and MSE_min is the value of the best solution found so far. We call this as a simple bounding.

It has been shown in ^[27] that the merge costs of the PNN method are monotonically increasing if the cluster pair with minimum cost is always merged. We denote the series of merge costs by d₁, d₂, ..., d_N-M, where d_t is the merge cost at the t^th merge. The monotony property implicates:

The criterion has been shown to apply to the PNN method where we always select the merge with minimum cost. From this property, we could derive a stronger termination criterion for the branch-and-bound algorithm:

Bounding criterion 2: $MS{E_{\text{t}}} + \left( {N - M - t} \right) \cdot {d_t} \geqslant MS{E_{\min }}$.

Here (N-M-t) indicates the number of forth coming merges in the algorithm, and d_t is the previous merge cost. This bounding criterion is based on the assumption that all forth-coming merge operations increase the MSE no less than the previous merge operation. This would allow termination of sub-optimal solution earlier than the criterion 1.

The only problem of using the bounding criterion 2 is that the monotony property does not necessarily hold true in the algorithm of Fig. 4. In the branch-and-bound method, we can also perform sub-optimal merges, which can result in a non-monotonic series of merge costs. Consequently, we could terminate a path to the optimal solution because of using the stronger criterion.

The optimal solution can be reached via several different search paths. It is expected that at least one of these paths fulfills the monotony property, and therefore, termination of other redundant paths would not be a big problem. The algorithm in Section 3.3, however, generates non-redundant search tree and there is no guarantee that only path to the optimal solution would meet the monotony property. The consequence of this is that optimality cannot be guaranteed if the stronger termination criterion was used with the non-redundant search tree.

5.   Polynomial time variants

The time complexity of the branch-and-bound technique is exponential regardless the bounding criterion used. The practical usability of the method is therefore limited to small special cases only. We propose next two sub-optimal variants that compromise the optimality but work in polynomial time.

5.1.   Piecewise optimization

The first method, called Piecewise optimization, divides the original problem into a series of smaller sub problems that are solved independently. The input at each stage of the algorithm is N clusters (whole data set in the beginning), and the output is the optimal clustering to N-Z clusters, where Z is a parameter of the algorithm. The result is then input to the same procedure, and the process is repeated until the desired number of M clusters is reached. The method is illustrated in Fig. 6, and its pseudo code given in Fig. 7.

If the size of the search tree of a single sub-problem is Z, we need to repeat the algorithm $\left\lceil {(N - M)/Z} \right\rceil$ times. The quality of the result depends on the parameter Z; greater values will give better clustering result at the cost of longer run time. The extreme case is when we set Z = N-M, which would result to the same algorithm as the branch-and-bound technique in Section 3. By setting Z = 1, on the other hand, the method would be the same as the PNN method.

At the starting point of the algorithm when we have N clusters as the input, we have O(N²) possible merge operations to be considered. Two subsequent merge operations result in O(N⁴) different alternatives. In general, the time complexity of z subsequent merge operations is O(N^2Z), and the overall Piecewise algorithm (N/Z)·O(N^2Z) = O(N^2Z+1/Z). The algorithm works in polynomial time if Z is small enough to be considered as a constant. On the other hand, the time complexity increases exponentially as a function of Z. The algorithm is therefore useful only with very small values such of Z. For example, the time complexities for Z = 2 and Z = 3 are O(N⁵) and O(N⁷), respectively.

5.2.   Look-ahead optimization

The Piecewise optimization traverses the search tree through the local optima. At each stage of the algorithm, the last merge is the most critical because the algorithm considers only one level further. Slightly better (but Z times slower) variant, denoted as Look-ahead optimization, can be designed by the following modifications.

As in the Piecewise optimization, we generate complete search tree to the level Z and search for the optimal clustering with N-Z clusters. Instead of moving to this local optimum at the level Z, we proceed only one level in the tree along the path towards the direction of the local optimum. After this, we regenerate a completely new search tree starting from the level N-1, and then repeat the procedure N-M times. The process is illustrated in Fig. 8.

This variant is less critical for the local minima in practice. As a drawback, the time complexity of the algorithm is Z times that of the previous variant; that is O(N^2Z+1). With large Z-values, we would also do unnecessary work as a part of the search tree would be re-generated several times. In practice, the algorithm can be realized only with very small Z-values.

6.   Experiments

We consider the data sets that are summarized in Table 2. The first set (B) includes a randomly selected pixel blocks from a 256 × 256 size image Bridge. The second set (S) is artificially generated two-dimensional data set with varying complexity in terms of spatial data distributions with M = 15 predefined clusters. The third set (SS2) is a standard clustering test problem of ^[25], pp. 103−104. The data set contains 89 postal zones in Bavaria (Germany) and their attributes are the number of self-employed people, civil servants, clerks and manual workers in these areas. The attributes are normalized to the scale [0, 1] according to their minimum and maximum values.

The results with the B sets are summarized in Figs. 9−11 with different number of clusters (M = 2, 5, 9) using Pentium 450 MHz computer. Comparative results are given for the following methods:

● Partition-based full search (implementation by Juha Kivijärvi).

● Partial partition-based branch-and-bound ^[10]

● Merge-based full search

● Merge-based branch-and-bound

6.1.   Results

Proposed methods work faster when the number of clusters is small (M = 2), or large (M = 9). These correspond to situations, in which the number of possible solutions is smallest. The use of bounding criterion improves the full search remarkably when the number of clusters is large (M = 5 and M = 9). In these cases, the proposed merge-based branch-and-bound is also faster than the partition-based counterpart. In the case of M = 2, on the other hand, the partition-based approach is superior. This is reasoned by the fact that the more merge steps are required the less there are clusters in the solution. We can also conclude that the practical usability of any of the optimal algorithms tested here are limited to very small size problem instances only.

Similar results are also reported in Fig. 12 for the polynomial time variants in the case of set S. The polynomial time variants are faster and can process larger data sets than the optimal branch-and-bound technique although only small sub tree sizes (Z = 2, 3) were applied.

The qualities of the sub-optimal variants are compared in Table 3 with other heuristic methods. Random clustering is obtained by randomly selecting M data points and using them as cluster representatives, and then partition the data set optimally by minimizing Euclidean distance from the data points to the cluster representatives. K-means is an iterative clustering algorithm, which is also known as GLA in vector quantization context due to ^[28]. PNN is implemented as in ^[24], and GAIS refers to a genetic algorithm which is the best clustering algorithm that we experimented in ^[29].

In comparison to PNN, the polynomial time variants (Piecewise and Look-ahead) manage to produce better clustering in two cases (N = 70, 80) but the difference is marginal. For the rest of the data there are no differences between PNN and branch-and-bound variants. GAIS, on the other hand, gives somewhat better results than PNN and BB. It is noted that the optimality of these results is not known but it is expected that GAIS results are very close to optimum.

The results for SS2 are summarized in Table 4. They indicate that the Piecewise and Look-ahead algorithms can improve over PNN and K-means but at the cost of significant increase in run time. The results are still worse than that of GAIS. The additional time spent for the improvement is simply not worth as better results can be obtained sophisticated heuristics like genetic algorithm ^[29] much faster.

6.2.   Discussion

Here we discuss the pros and cons of the studied merge-based branch-and-bound approach and compare it against alternative splitting approach. This approach puts all the points in a single cluster and then splits one cluster at a time until M clusters are reached. These two approaches are compared in Fig. 13.

The main benefit of the merge-based approach is that the search tree is narrower. There are only O(N²) possible merges at each step compared to O(2^N) of possible splits at every step. The split approach would also require minimum redundancy search tree to avoid the same clustering result to be found via different order of splits. The merge-based search tree is also deeper as there are O(N) merges in total compared to O(M) of split. Theoretically this provides more effective branch-and-bound as each successful bound would cut bigger proportion of the tree while the opposite splitting approach reaches candidate solutions faster in M steps.

However, a bigger problem is that the last merge is the costliest of all merges. For example, MSE-values for M = (30, 20, 10, 7, 2) are (0.036, 0.070, 0.191, 0.314, 4.086). The early merges are therefore insignificant compared to the ones at the end, and sub-optimal merges at the earlier stages rarely trigger the bounding criterion. The situation of the splitting approach would be opposite.

When compared to the results in literature, we see that Brusco's method ^[14] found optimal solutions up to N = 60 with M = 8 clusters in about 6½ hours. This is significantly larger than any of our results which remained smaller than N = 20. The bounding criterion of the merge-based approach is too poor because of the increasing nature of the merge errors. The opposite splitting-based approach is expected to work better due to this reason.

One factor for better efficiency is the nearest neighbor separation in ^[14], which places the points with smallest pairwise distance at the opposite ending in the ordering. The idea is to have diverse dataset early where disruptive points appear early in the ordering. This makes even smaller mistakes to trigger the bounding criterion. Kaminka ^[15] enhanced this further by applying furthest point heuristic by Gonzalez ^[30]. This heuristic (also known as minmax) has been shown to provide much better initialization for k-means than random centroids in ^[31].

Another important factor is to choose the initial solution well. In the context of travelling salesperson problem, an undocumented rule of thumb states that the initial solution should be as good as 1% for the branch-and-bound to be significantly more efficient to exhaustive search according to ^[32]. In our branch-and-bound we have used the PNN algorithm but the results for SS2 (0.336 with M = 7) are actually 2.6% worse than obtained by GAIS (0.314).

Brusco ^[14] also reached optimal solutions for problems size up to N = 240 with M = 8 clusters in case of well-separated clusters. However, such data is significantly easier, and we can expect to reach optimal results by a simple but effective heuristic algorithm called random swap ^[33] even for higher values of N. The only purpose of branch-and-bound would then be merely to verify the optimality of the result.

Finally, we note that our results have been ran already in 2004 using 450 MHz computer, around the same time that those of Brusco and others. At first sight one might expect these results being already quite old. However, our current machine runs in 3.7 GHz which is merely 8:1 increase in speed from 2004 to 2022 (see also Fig. 14). Since the problem has exponential in nature, time complexity of 2^N would mean that we could solve problems of size N+3 in 2022 compared to problem of sizes N in 2004. This is hardly significant.

Our machine itself has 20 cores available for single test research experiments. While branch-and-bound indeed would parallelize nicely, this would mean another 20:1 increase, so we could hope to solve problems of size N+7.

7.   Conclusions

We have introduced a merge-based branch-and-bound technique to solve optimal clustering. The proposed algorithm works better than the partition-based counterparts when the number of clusters is high. In the case of a small number of clusters, however, the partition-based approach works faster. Nevertheless, all variants have exponential time complexity, and therefore, the results are mainly of theoretical interest only.

Two polynomial time algorithms were also introduced inspired by the proposed branch-and-bound technique but with limited success. Further improvement could be achieved by using stronger bounding criteria in the polynomial time variants.