Efficient Centroid-Linkage Clustering

Authors: mohammadhossein  bateni
laxman  dhulipala
willem  fletcher
kishen n gowda
d  ellis hershkowitz
jakub  łącki
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 18:39:24

Abstract: We give an algorithm for Centroid-Linkage Hierarchical Agglomerative Clustering (HAC), which computes a $c$-approximate clustering in roughly $n^{1+O(1/c^2)}$ time. We obtain our result by combining a new centroid-linkage HAC algorithm with a novel fully dynamic data structure for nearest neighbor search which works under adaptive updates.

We also evaluate our algorithm empirically. By leveraging a state-of-the-art nearest-neighbor search library, we obtain a fast and accurate centroid-linkage HAC algorithm. Compared to an existing state-of-the-art exact baseline, our implementation maintains the clustering quality while delivering up to a $36\times$ speedup due to performing fewer distance comparisons.

Main paper: 
 introduction: 
 Hierarchical Agglomerative Clustering (HAC) is a widely-used clustering method, which is available in many standard data science libraries such as SciPy [39] , scikit-learn [34] , fastcluster [29] , Julia [20] , R [35] , MATLAB [28] , Mathematica [25] and many more [31, 32, 36] . HAC takes as input a collection of n points in R d . Initially, it puts each point in a separate cluster, and then proceeds in up to n -1 steps. In each step, it merges the two "closest" clusters by replacing them with their union.The formal notion of closeness is given by a linkage function; choosing different linkage functions gives different variants of HAC. In this paper, we study Centroid-Linkage HAC wherein the distance between two clusters is simply the distance between their centroids. This method is available in most of the aforementioned libraries. Other common choices include single-linkage (the distance between two clusters is the minimum distance between a point in one cluster and a point in the other cluster) or average-linkage (the distance between two clusters is the average pairwise distance between them). HAC's applicability has been hindered by its limited scalability. Specifically, running HAC with any popular linkage function requires essentially quadratic time under standard complexity-theory assumptions. This is because (using most linkage functions) the first step of HAC is equivalent to finding the closest pair among the input points. In high-dimensional Euclidean spaces, this problem was shown to (conditionally) require Ω(n 2-α ) time for any constant α > 0 [23] . In contrast, solving HAC in essentially Θ(n 2 ) time is easy for many popular linkage functions (including centroid, single, average, complete, Ward), as it suffices to store all-pairs distances between the clusters in a priority queue and update them after merging each pair of clusters.In this paper, we show how to bypass this hardness for Centroid-Linkage HAC by allowing for approximation. Namely, we give a subquadratic time algorithm which, instead of requiring that the two closest clusters be merged in each step, allows for any two clusters to be merged if their distance is within a factor of c ≥ 1 of that of the two closest clusters.Contribution 1: Meta-Algorithm for Centroid-Linkage. Our first contribution is a simple metaalgorithm for approximate Centroid-Linkage HAC. The algorithm assumes access to a dynamic data structure for approximate nearest-neighbor search (ANNS). A dynamic ANNS data structure maintains a collection of points S ⊂ R d subject to insertions and deletions and given any query point u ∈ S returns an approximate nearest neighbor v ∈ S. Formally, for a c-approximate NNS data structure, v satisfies D(u, v) ≤ c • min x∈S\{u} D(u, x).1 Here, D(•, •) denotes the Euclidean distance. By using different ANNS data structures we can obtain different Centroid-Linkage HAC algorithms (which is why we call our method a meta-algorithm). We make use of an ANNS data structure where the set of points S is the centroids of the current HAC clusters. Roughly, our algorithm runs in time Õ(n) times the time it takes an ANNS data structure to update or query.While finding distances is the core part of running HAC, we note that access to an ANNS data structure for the centroids does not immediately solve HAC. Specifically, an ANNS data structure can efficiently find an (approximate) nearest neighbor of a single point, but running HAC requires us to find an (approximately) closest pair among all clusters. Furthermore, HAC must use the data structure in a dynamic setting, and so an update to the set of points S (caused by two clusters merging) may result in many points changing their nearest neighbors.Our meta-algorithm requires that the dynamic ANNS data structure works under adaptive updates.2 Namely, it has to be capable of handling updates, which are dependent on the prior answers that it returned (as opposed to an ANNS data structure which only handles "oblivious" updates).To illustrate this, consider a randomized ANNS data structure D. Clearly, a query to D often has multiple correct answers, as D can return any near neighbor within the promised approximation bound. As a result, an answer to a query issued to D is dependent on the internal randomness of D. Let us assume that D is used within a centroid-linkage HAC algorithm A and upon some query returns a result u. Then, algorithm A uses u to decide which two clusters to merge, and the centroid p of the newly constructed cluster is inserted into D. Since p depends on u, which in turn is a function of the internal randomness of D, we have that the point that is inserted into D is dependent on the internal randomness of the data structure. Hence, it is not sufficient for A to return a correct answer for each fixed query with high probability (meaning at least 1 -1/poly(n)). Instead, A must be able to handle queries and updates dependent on its internal randomness. We note that a similar issue is prevalent in many cases when a randomized data structure is used as a building block of an algorithm. As a result, the subtle notion of adaptive updates is a major area of study [40, 7, 33, 22, 21] .Contribution 2: Dynamic ANNS Robust to Adaptive Updates. Our second contribution is a dynamic ANNS data structure for high-dimensional Euclidean spaces which, to the best of our knowledge, is the first one to work under adaptive updates.To obtain our ANNS data structure, we show a black-box reduction from an arbitrary randomized dynamic ANNS data structure to one that works against adaptive updates (see Theorem 2) . The reduction increases the query and update times by only a O(log n) factor. We apply this reduction to a (previously known) dynamic ANNS data structure which is based on locality-sensitive hashing [3, 19] and requires non-adaptive updates.By combining our dynamic ANNS data structure for adaptive updates with our meta-algorithm, we obtain an O(c)-approximate Centroid-Linkage HAC algorithm which runs in time roughly n 1+1/c 2 .Furthermore, our data structure for dynamic ANNS will likely have applications beyond our HAC algorithm. Specifically, ANNS is commonly used as a subroutine to speed up algorithms for geometric problems. Oftentimes, and similarly to HAC, formal correctness of these algorithms require an ANNS data structure that works for adaptive updates but prior work often overlooks this issue. Our new data structure can be used to fix the analysis of such works. For instance, an earlier work on subquadratic HAC algorithms for average-and Ward-linkage [1] overlooked this issue and our new data structure fixes the analysis of this work.3 Similarly, a key bottleneck in the dynamic k-centers algorithm in [5] was that it had to run nearest neighbor search over the entire dataset (instead of just the k-centers). By replacing the ANNS data structure used therein with our ANNS algorithm, we believe that that the n ϵ running times from that paper could be improved to k ϵ .Contribution 3: Empirical Evaluation of Centroid-Linkage HAC. Finally, we use our Centroid-Linkage meta algorithm to obtain an efficient HAC implementation. We use a state-of-the-art static library for ANNS that we extend to support efficient dynamic updates when merging centroids.We empirically evaluate our algorithm and find that our approximate algorithm achieves strong fidelity with respect to the exact centroid algorithm, obtaining ARI and NMI scores that are within 7% of that of the exact algorithm, while achieving up to 36× speedup over the exact implementation, even when both implementations are run sequentially. Our implementation can be found at https://github.com/kishen19/CentroidHAC.We note that the ANNS we use in the empirical evaluation of our meta-algorithm is different from our new dynamic ANNS (Contribution 2). This is because modern practical methods for ANNS (e.g., based on graph-based indices [38, 18] ) have far surpassed the efficiency of data structures for which rigorous theoretical bounds are known. Closing this gap is a major open problem.Improving the efficiency of HAC has been an actively researched problem for over 40 years [8, 31, 30, 32] . A major challenge common to multiple linkage functions has been to improve the running time beyond Θ(n 2 ), which is the time it takes to compute all-pairs distances between the input points.For the case of low-dimensional Euclidean spaces, a running time of o(n 2 ) is possible for the case of squared Euclidean distance and Ward linkage, since the ANNS problem becomes much easier in this case [41] . However, breaking the Θ(n 2 ) barrier is (conditionally) impossible in the case of high-dimensional Euclidean spaces (without an exponential dependence on the dimension) [23] .To bypass this hardness result, a recent line of work focused on obtaining approximate HAC algorithms [26] . Most importantly, Abboud, Cohen-Addad and Houdrouge [1] showed an approximate HAC algorithm which runs in subquadratic time for Ward and average linkage. These algorithms rely on the fact in the case of these linkage functions the minimum distance between two clusters can only increase as the algorithm progresses, which is not the case for centroid linkage; e.g. consider 3 equi-distant points as in Figure 1 . We also note that, as mentioned earlier, our dynamic ANNS data structure fixes a gap in the analysis of both algorithms in the paper. In fact, not only can the minimum distances shrink when performing Centroid-Linkage HAC, but when doing O(1)approximate Centroid-Linkage HAC distances can get arbitrarily small; see Figure 2 . This presents an additional issue for dynamic ANNS data structures which typically assume lower bounds on minimum distances.v 1 v 2 v 4 v 3 2 2 2 1 (a) First merge. v 1 v 2 v 4 v 3 2 2 2 (b) Result of merge. v 1 v 2 v 4 v 3 2 2 2 3 2 (c) Second merge. v 1 v 2 v 4 v 3 2 2 2 (d) Result of merge.Compared to the algorithm of Abboud, Cohen-Addad and Houdrouge, our centroid-linkage HAC algorithm introduces two new ideas. First, we show how to handle the case when the minimum 3We contacted the authors of [1] and they confirmed this gap in their analysis. distance between two clusters decreases as a result of two cluster merging. Second, we introduce an optimization that allows us not to consider each cluster at each of the (logarithmically many) distance scales. This optimization improves the running time bound by a logarithmic factor when each merge "makes stale" a small number of stored pairs and in practice the number of such stale merges occur only 60% of the time (on average) compared to the number of actual merges performed.Another line of work considered a different variant of HAC, where the input is a weighted similarity graph [15] [16] [17] 6] . The edge weights specify similarities between input points, and a lack of edge corresponds to a similarity value of 0. By assuming that the graph is sparse, this representation allows bypassing the hardness of finding distances, which leads to near-linear-time approximate algorithms for average-linkage HAC [15] , as well as to efficient parallel algorithms [16, 6] .In terms of the theoretical and empirical part of this paper considering different algorithms, such a gap has also been observed for various other prominent clustering algorithms. For instance, for correlation clustering, while the best known algorithms are obtained by solving a linear program (see [9] ), in practice a simple local search/Louvain-based algorithm is used (see section 4.1 of [2] ). For modularity clustering, the best known approximation is obtained by solving an SDP [24] . Again, Louvain-based algorithms are used in practice. In the case of k-means, while the best known approximation is a constant obtained via the primal-dual method [10] , in practice Lloyd's heuristic paired with k-means++ seeding is used which has a logarithmic approximation guarantee. For balanced graph partitioning, coarsening combined with a brute-force algorithm is used to obtain the best results in practice. On the other hand, the algorithms for solving the corresponding theoretical formulations, i.e. balanced cut and multiway cut, are entirely different.If D(u, v) ≤ r then Pr h∼H (h(u) = h(v)) ≥ p 1 . 
preliminaries: 
 In this section, we review our conventions and preliminaries. Throughout this paper, we will work in d-dimensional Euclidean space R d . We will use D(•, •) for the Euclidean metric in R d ; that is, for x, y ∈ R d , the distance between x and y is D(x, y):= i (x[i] -y[i]) 2where x[i] and y[i]are the ith coordinates of x and y respectively. Likewise, we let D(x, Y ) := min y∈Y D(x, y) for Y ⊆ R d .While we have described centroid-linkage in a cluster-centric way (see Appendix A for a formal cluster-centric definition), it will be more useful to use an equivalent centroid-centric definition. Specifically, we can equivalently define c-approximate Centroid-Linkage HAC as repeatedly "merging centroids". Initialize the set of all centroids C = P and a weight w u = 1 for each u ∈ C; these weights will encode the cluster sizes of each centroid. Then, until |C| = 1 we do the following: let x, ŷ ∈ C be a pair of centroids satisfyingD(x, ŷ) ≤ c • min x,y∈C D(x, y)where the min is taken over distinct pairs; merge x and ŷ into their weighted midpoints by removing them from C, adding z = (w x x + w ŷ ŷ)/(w x + w ŷ ) to C and setting w z = w x + w ŷ . Also, note that (exact) Centroid-Linkage HAC is just the above with c = 1.We define a dynamic ANNS data structure. Typically, dynamic ANNS requires a lower bound on distances; we cannot guarantee this for centroid HAC. Thus, we make use of β additive error below. u, ū ∈ R d return v ∈ S \ {ū} s.t. D(u, v) ≤ α • D(u, S \ {ū}) + β.If N is a dynamic (α, 0)-approximate NNS data structure then we will simply call it α-approximate. For a set of points S, we will use the notation NonAdaptiveAnn(S) to denote the result of starting with an empty dynamic ANNS data structure and inserting each point in S (in an arbitrary order).We say N succeeds for a set of points S ⊆ R d and a query point u if after starting with an empty set and any sequence of insertions and deletions of points in S we have that the query N .query(u, ū) is correct for any ū (i.e., the returned point satisfies the stated condition in 3). If the data structure is randomized then we say that it succeeds for adaptive updates if with high probability (over the randomness of the data structure) it succeeds for all possible subsets of points and queries. Queries have true distances at most ∆ if the (true) distance from S to any query is always at most ∆.The starting point for our dynamic ANNS data structure is the following data structure requiring non-adaptive updates which alone does not suffice for centroid HAC as earlier described. See Appendix B for a proof. Theorem 1 (Dynamic ANNS for Oblivious Updates, [3, 19] ). Suppose we are given γ > 1 and c, β, ∆ and n where log(∆/β), γ ≤ poly(n) and a dynamically updated set S with at most n insertions. Then, if all queries have true distance at most ∆, we can compute a randomized (O(c), β)-approximate NNS data structure with update, deletion and query times of n 1/c 2 +o (1) • log(∆/β) • d • γ, which for a fixed set of points and query point succeeds except with probability at most exp(-γ).If D(u, v) ≥ cr then Pr h∼H (h(u) = h(v)) ≤ p 2 .We will make use of known locality sensitive hashing schemes for Euclidean distance in R d . For this result, see also [4] . Theorem 6 (Lemma 3.2.3 of [3] ). Given r and c, there is a (r, cr, p 1 , p 2 )-sensitive LSH family withlog(1/p 1 ) log(1/p 2 ) = 1/c 2 + o(1) and log(1/p 2 ) ≤ o(log n).Furthermore, one can sample h ∼ H and given u ∈ R d compute h(u) in time dn o (1) . Theorem 1 (Dynamic ANNS for Oblivious Updates, [3, 19] ). Suppose we are given γ > 1 and c, β, ∆ and n where log(∆/β), γ ≤ poly(n) and a dynamically updated set S with at most n insertions. Then, if all queries have true distance at most ∆, we can compute a randomized (O(c), β)-approximate NNS data structure with update, deletion and query times of n 1/c 2 +o (1) • log(∆/β) • d • γ, which for a fixed set of points and query point succeeds except with probability at most exp(-γ).Proof. We first construct a data structure which succeeds with constant probability for a given query.To do so we begin by turning our LSH hashes from Theorem 6 into a hash function which is an "L ors of K ands". In particular, we letL := n 1/c 2 +o(1) • Θ(log n).and letK := 1 log(1/p 2 )• (Θ(log n) + log L + log log(∆/β)) .Notice that by Theorem 6 we have that log(1/p 2 ) ≤ o(log n) and by assumption we have log(∆/β), γ ≤ poly(n) and so we haveEQUATIONFor each i such that β ≤ 2 i ≤ ∆, we define L-many new hash functions, each of which consists of K-many h from Theorem 6. In particular, for i ∈ [log(∆/β)], l ∈ [L] and k ∈ [K], we let h ikl be a uniformly random h sampled according to Theorem 6 using radius r = 2 i and our input c. Likewise, we let g il : R d → R k be a new hash function defined on u ∈ R d asg il (u) = (h il1 (u), h il2 (u), . . . , h ilk (u)).We let G i := {g il : l ≤ L} be all such hash functions for i and let G = i G i be all such hash functions across scales.By Theorem 6, sampling each constituent hash functions for each g ∈ G can be done in time dn o (1) and so by this and Equation ( 1) initializing all hash functions in G takes timeEQUATIONThis initialization time will be dominated by the time it takes to insert a single point. For each such hash function we will maintain a collection of buckets which partitions points in S.To insert a point into S, we simply compute its hash for each g ∈ G and store this hash in each corresponding bucket. To delete a point from S we simply delete it from all buckets to which it hashes. Lastly, to answer a query on (u, ū), we hash u according to each g ∈ G and then return the v ∈ S \ ū such that for minimum i we have a g ∈ G i such that g(u) = g(v). Since by Theorem 6 evaluating each hash function takes time dn o (1) , we get that all of these operations also take time at most that given by Equation ( 2).Lastly, we argue the (constant probability) correctness of a query for a fixed set of points. Fix query points (u, ū) and a set of points S. Let w ∈ S \ {u} be a fixed point of S \ {u} where we let i w be such that 2 iw ≤ D(u, w) ≤ 2 iw+1 . Say that w succeeds if 1. One Correct Collision: there is a g ∈ G iw such that g(w) = g(u) and 2. No Incorrect Collisions: for all i < i w and g ∈ G i we have that g(u) ̸ = g(w).Observe that if both of the cases happen for every w ∈ S \ {u} then we have that our query is correct. The only non-trivial part of this observation is the fact that if u is within distance β of S \ {ū} then i w = 1 and the returned point is within an additive β of u so long as one correct collision holds.Continuing, we next lower bound the probability of one correct collision. For a fixed g iwl ∈ G iw , we have that g(w) = g(u) iff h iwlk (w) = h iwlk (u) for all k ≤ K. The probability of this for one k ≤ K is, by definition, p 1 . By Theorem 6 we know that log(1/p1) log(1/p2) = 1/c 2 + o(1) and so for a fixed l, we have that this occurs with probability at leastp K 1 = exp(-K • log(1/p 1 )) = exp - log(1/p 1 ) log(1/p 2 ) • (O(log n) + log L + log log(∆/β)) = exp(-(1/c 2 + o(1)) • (O(log n) + log L + log log(∆/β))) ≥ exp -(1/c 2 + o(1)) • O(log n)It follows that the probability that this does not hold for some l is at most (1) .1 -p K 1 L ≤ exp(-L • p K 1 ) ≤ exp -L • exp -(1/c 2 + o(1)) • O(log n) = exp -Ω(log n) -exp[((1/c 2 ) + o(1)) • Ω(log n)] • exp -(1/c 2 + o(1)) • O(log n) = n -ΩNext, we lower bound the probability of no incorrect collisions. Fix an i < i w . Observe that for a given h ilk we have, by definition, that h ilk (u) = h ilk (w) with probability at most p 2 . Thus, we have g il (u) = g il (w) with probability at most p K 2 = exp(-log(1/p 2 )K) = exp(-Ω(log n) -log L -log log(∆/β)). Union bounding over our at most log(∆/β)-many such i and L-many such l, we get that that we have no incorrect collisions except with probability at most 1) . Finally, union bounding over our n-many possible points in S and using the fact that γ ≥ 1 we have that we succeed except with probability at least a fixed constant bounded away from 0 for S and this query. Furthermore, by taking γ-many independent repetitions of this data structure and taking the returned point closest to our query, we increase our runtime by a multiplicative γ and reduce our failure probability to exp(-γ).p K 2 • log(∆/β) • L = exp(-Ω(log n)) = n -Ω( 
dynamic anns with adaptive updates: 
 The main theorem we show in this section is how to construct a dynamic ANNS data structure that succeeds for adaptive updates using one which succeeds for oblivious updates. Theorem 2 (Reduction of Dynamic ANNS from Oblivious to Adaptive). Suppose we are given c, β, ∆, n, s ∈ R d and dynamically updated set S ⊂ R d with at most n insertions, such that all inserted points and query points lie in B s (∆). Moreover, assume that we can compute a dynamic (c, β)-approximate NNS data structure with query, deletion and insertion times T Q , T D and T I which succeeds for S and a fixed query except with probability at most n -O(d log d log(∆/β)) . Then, we can compute a randomized dynamic (c, cβ)-approximate NNS data structure that succeeds for adaptive updates with query, deletion and amortized insertion timesO(log n) • T Q , O(T D ) and O(log n) • T I .We note that our result is slightly stronger than the above: the insertions we make into the oblivious ANNS that we use only occur upon their instantiation, not dynamically.As a corollary of the above reduction and known constructions for dynamic ANNS that work for oblivious updates-namely, Theorem 1 with γ = Θ(d log d • log(∆/β) • log n) and using parameter β ′ = β/c where β ′ is the additive distortion parameter for Theorem 1-we obtain the following. Theorem 3 (Dynamic ANNS for Adaptive Updates). Suppose we are given c, β, ∆, n, s ∈ R d where log(∆/β), d, c ≤ poly(n) and dynamically updated set S ⊂ R d with at most n insertions, such that all inserted and query points lie in B s (∆).Then, we can compute a randomized dynamic (O(c), β)-approximate NNS data structure that succeeds for adaptive updates with query, deletion and amortized insertion timen 1/c 2 +o(1) • log(∆/β) • d 2 .A previous work [19] also provided algorithms that work against an "adaptive updates" but only a set of adaptive updates made against a fixed set of query points. Also, note that if we have a lower bound of δ on the true distance of any query then lowering c by a constant and setting β = Θ(δ) for a suitably small hidden constant gives an O(c)-approximate NNS with similar guarantees.S 0 S 1 S 2 𝒩 0 𝒩 1 𝒩 2 (a) Initial state. S 0 S 1 S 2 𝒩 0 𝒩 1 𝒩 2 (b) Insert point. S 0 S 1 S 2 𝒩 0 𝒩 1 𝒩 2 (c) S1 ← S0 ∪ S1. S 0 S 1 S 2 𝒩 0 𝒩 1 𝒩 2 (d) S2 ← S1 ∪ S2.Figure 3 : Our merge-and-reduce strategy when a point (in green) is inserted.Our algorithm for dynamic ANNS for adaptive updates uses two ingredients. First, we make use of the following "covering nets" to fix our queries to a small set against which we can union bound.Lemma 1. Given s ∈ R d and β, ∆ > 0, there exists a covering netQ ⊆ R d such that 1. Small Size: |Q| ≤ (∆/β) O(d log d) 2. Queries: given u ∈ B s (∆), one can compute u ′ ∈ Q in time O(d) such that D(u, u ′ ) ≤ β.Second, we make use of a "merge-and-reduce" approach. A similar approach was taken by [12] for exact k-nearest-neighbor queries in the plane. Namely, for each i ∈ [O(log n)] we maintain a set S i of size at most 2 i (where all S i partition all inserted points) and a dynamic ANNS data structure N i for S i which only works for oblivious updates. Other than the size constraint, the partition is arbitrary. Informally, we perform deletions, insertions and queries as follows.Deletion: To delete a point we simply delete it from its corresponding N i .Insertion: To insert a point, we insert it into S 0 and for each i we move all points from S i to S i+1 if S i contains more than 2 i points; we update N i accordingly each time we move points; namely, we recompute N i and N i+1 from scratch on S i and S i+1 respectively. See Figure 3 .Query: Lastly, to query a point u we first map this point to a point in our covering net u ′ (as specified by Lemma 1), query u ′ in each of our N i and then return the best output (i.e., the point output by anN i that is closest to u ′ ).Our algorithm is more formally described in pseudo-code in Algorithm 1. Let i be such that u ∈ S i N i .delete(u) and remove u from S i function query(u, ū)Input: β, ∆, n, s ∈ R d , Q ← CoveringNet(s, ∆, β) (computedLet u ′ ∈ Q be such that D(u, u ′ ) = O(β) ▷ computed using Lemma 1 Let v i = N i .query(u ′ , ū) return v = arg min vi D(u ′ , v i )(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. • The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
centroid-linkage hac algorithm: 
 In this section, we give our algorithm for Centroid-Linkage HAC. Specifically, we show how to (with Algorithm 3) reduce approximate Centroid-Linkage HAC to roughly linearly-many dynamic ANNS function calls, provided the ANNS works for adaptive updates.Theorem 4 (Reduction of Centroid HAC to Dynamic ANN for Adaptive Updates). Suppose we are given a set of n points P ⊂ R d , c > 1, ϵ > 0, and lower and upper bounds on the minimum and maximum pairwise distance in P of δ and ∆ respectively. Suppose we are also given a data structure that is a dynamic c-approximate NNS data structure for all queries in Algorithm 3 that works for adaptive updates with insertion, deletion and query times of T I , T D , and T Q .Then there exists an algorithm for c(1 + ϵ)-approximate Centroid-Linkage HAC on n points that runs in time Õ(n• (T I + T D + T Q )) assuming ∆ δ ≤ poly(n) and 1 ϵ ≤ polylog(n).Our final algorithm for centroid HAC is not quite immediate from the above reduction and our dynamic ANNS algorithm for adaptive updates. Specifically, the above reduction requires a c-approximate NNS data structure but in the preceding section we have only provided a (c, β)-approximate NNS data structure for β > 0. However, by leveraging the structure of centroid HAC-in particular, the fact that there is only ever at most one "very close" pair-we are able to show that this suffices. In particular, for a given c, there exists a c 0 and β 0 such that a (c 0 , β 0 )-approximate NNS data structure functions as a c-approximate NNS data structure for all queries in Algorithm 3. This is formalized by Lemma 6. Combining Lemma 6, Theorem 4 and Theorem 3, we get the following.Theorem 5 (Centroid HAC Algorithm). For n points in R d , there exists a randomized algorithm that given any c > 1.01 and lower and upper bounds on the minimum and maximum pairwise distance in P of δ and ∆ respectively, computes a c-approximate Centroid-Linkage HAC with high probability.It has runtime O(n 1+O(1/c 2 ) d 2 ) assuming c, d, ∆ δ ≤ poly(n).In what remains of this section, we describe the algorithm for Theorem 4. Consider a set of points P ∈ R d . Suppose we have a dynamic c-approximate NNS data structure N that works with adaptive updates with query time T Q , insertion time T I , and deletion time T D . Let Q be a priority queue storing elements of the form (l, x, y) where x and y are "nearby" centroids and l = D(x, y). The priority queue supports queuing elements and dequeueing the element with shortest distance, l. At any given time while running the algorithm, say a centroid x is queued if there is an element of the form (l, x, y) in Q. We will maintain a set, C, of active centroids allowing us to check if a centroid is active in constant time. For each active centroid, C will also store the weight of the centroid so we can preform merges and an identifier to distinguish between distinct clusters with the same centroid.Note that any time we store a centroid, including in Q, we will implicitly store this identifier.First, we describe how the algorithm handles merges. To merge two active centroids, x and y, we remove them from C and delete them in N . Let z = w x x + w y y be the centroid formed by merging x and y. We use N to find an approximate nearest neighbor, y * , of z and add to Q the tuple (D(z, y * ), z, y * ). Lastly, we add z to N and C. Pseudo-code for this algorithm is given by Algorithm 2. Crucially, we do not try to update any nearest neighbors of any other centroid at this stage. We will do this work later and only if necessary. Since we are using an approximate NNS, it is possible that the centroid z will be the same point in R d as the centroid of another cluster. We can detect this in constant time using C and we will immediately merge the identical centroids.We now describe the full algorithm using the above merge algorithm; we also give pseudo-code in Algorithm 3. Let ϵ > 0 be a parameter that tells the algorithm how aggressively to merge centroids.To begin, we construct N by inserting all points of P in any order. Then for each p ∈ P , we use N to find an approximate nearest neighbor y ∈ P \ {p} and queue (D(p, y), p, y) to Q. We also initialize C = {p : p ∈ P }. The rest of the algorithm is a while loop that runs until Q is empty. Each iteration of the while loop begins by dequeuing from Q the tuple (l, x, y) with minimum l.(1) If x and y are both active centroids then we merge them.(2) Else if x is not active then we do nothing and move on to the next iteration of the while loop.(3) Else if x is active but y is not, we use N to compute a new approximate nearest neighbor, y * of x.Let l * = D(x, y * ). If l * ≤ (1 + ϵ)l, then we merge x and y * . Otherwise we add (l * , x, y * ) to Q. y ← An approximate nearest neighbor of p returned by N .query(p, p)6: Queue (D(p, y), p, y) to Q 7: while Q is not empty do 8: (l, x, y) ← dequeue shortest distance from Q 9: if x, y ∈ C then 10: Handle-Merge(C, N , Q, x, y) 11: else if x ∈ C then 12:y * ← An approximate nearest neighbor of x returned by N .query(x, x) Queue (l * , x, y * ) to Q13: l * ← D(x, y * ) 14: if l * ≤ (1 + ϵ)l then 15: Handle-Merge(C, N , Q, x, y * ) 16: 
empirical evaluation: 
 We empirically evaluate the performance and effectiveness of our approximate Centroid HAC algorithm through a comprehensive set of experiments. Our analysis on various publicly-available benchmark clustering datasets demonstrates that our approximate Centroid HAC algorithm:(1) Consistently produces clusters with quality comparable to that of exact Centroid HAC, (2) Achieves an average speed-up of 22× with a max speed-up of 175× compared to exact Centroid HAC when using 96 cores in our scaling study. On a single core, it achieves an average speed-up of 5× with a max speed-up of 36×.As demonstrated in previous works, although LSH-based algorithms have strong theoretical guarantees, they tend to perform poorly in practice when compared to graph-based ANNS data structures [18] . For instance, in ParlayANN [18] , the authors find that the recall achievable by LSH-based methods on a standard ANNS benchmark dataset are strictly dominated by state-of-the-art graph-based ANNS data structures.Therefore, for our experiments, we use DiskANN, a widely-used state-of-the-art graph-based ANNS data structure [38] . In particular, we consider the in-memory version of this algorithm from ParlayANN, called Vamana. In a nutshell, the algorithm builds a bounded degree routing graph that can be searched using beam search. Note that this graph is not the k-NN graph of the pointset. The graph is constructed using an incremental construction that adds bidirectional edges between a newly inserted point, and points traversed during a beam search for this point; if a point's degree exceeds the degree bound, the point is pruned to ensure a diverse set of neighbors. See [38] for details.For Centroid HAC, we require the ANNS implementation to support dynamic updates. The implementation provided by ParlayANN currently only supports fast static index building and queries. Recently, FreshDiskANN [37] described a way to handle insertions and deletions via a lazy update approach. However, their approach requires a periodic consolidation step that scans through the graph and deletes inactive nodes and rebuilds the neighborhood of affected nodes, which is expensive.Instead, in our implementation of Centroid HAC, we adopt the following new approach that is simple and practical, and adheres to the updates required by our theoretical algorithm: when clusters u and v merge, choose one of them to represent the centroid, and update its neighborhood N (u) (without loss of generality) to the set obtained by pruning the set N (u) ∪ N (v). The intuition here is that N (u) ∪ N (v) is a good representative for the neighborhood of the centroid of u and v in the routing graph. Further, during search, points will redirect to their current representative centroid (via union-find [14] ), thus allowing us to avoid updating the in-neighbors of a point in the index. We believe application-driven update algorithms, such as the one described here, can help speed-up algorithms that uses dynamic ANNS as a subroutine.We evaluate the clustering quality of our approximate centroid HAC algorithm against ground truth clusterings using standard metrics such as the Adjusted Rand Index (ARI), Normalized Mutual Information (NMI), Dendrogram Purity, and the unsupervised Dasgupta cost [11] . Our primary objective is to assess the performance of our algorithm across various values of ϵ on a diverse set of benchmarks. We primarily compare our results with those obtained using exact Centroid HAC.For these experiments, we consider the standard benchmark clustering datasets iris, digits, cancer, wine, and faces from the UCI repository (obtained from the sklearn.datasets package). We also consider the MNIST dataset which contains images of grayscale digits between 0 and 9, and birds, a dataset containing images of 525 species of birds; see Appendix E for more details.Results. The experimental results are presented in Table 1 , with a more detailed quality evaluation in Appendix E.1. We summarize our results here.We observe that the quality of the clustering produced by approximate Centroid HAC is generally comparable to that of the exact Centroid HAC algorithm. On average, the quality is slightly better in some cases, but we attribute this to noise. In particular, we observe that for the value of ϵ = 0.1, we consistently get comparable quality to that of exact centroid: the ARI and NMI scores are on average within a factor of 7% and 2% to that of exact Centroid HAC, respectively. We also obtained good dendrogram purity score and Dasgupta cost as well, with values within 0.3% and 0.03%, respectively.Next, we evaluate the scalability of approximate Centroid HAC against exact Centroid HAC on large real-world pointsets. We consider the optimized exact Centroid HAC implementation from the fastcluster package [29] . We also implement an efficient version of exact Centroid HAC based on our framework (i.e. setting ϵ = 0 and using exact NNS queries) which has the benefit of using only linear space and supports parallelism in NNS queries. We also implement a bucket-based version of approximate Centroid HAC as a baseline, based on the approach of [1] along with an observation to handle the non-monotonicity of Centroid HAC; details in Appendix E. Results. We now summarize the results of our scalability study; see Appendix E.2 for a more details and plots. Figures 4a and 4b shows the running times on varying slices of the SIFT-1M dataset using one thread and 192 parallel threads, respectively. For the approximate methods, we considered ϵ = 0.1. Figure 4c further compares the heap and bucket based approaches as a function of ϵ.We observe that the bucket based approach is very slow for small values of ϵ due to many redundant nearest-neighbor computations. Yet, at ϵ = 0.1, both algorithms, with 192 threads, obtain speed-ups of up to 175× with an average speed-up of 22× compared to the exact counterparts. However, on a single thread, the bucket based approach fails to scale, while the heap based approach still achieves significant speed-ups of upto 36× with an average speed-up of 5×. Overall, our heap-based approximate Centroid HAC implementation with ϵ = 0.1 demonstrates good quality and scalability, making it a good choice in practice. 
conclusion: 
 In this work we gave an approximate algorithm for Centroid-Linkage HAC which runs in subquadratic time. Our algorithm is obtained by way of a new ANNS data structure which works correctly under adaptive updates which may be of independent interest. On the empirical side we have demonstrated up to 36× speedup compared to the existing baselines. An interesting open question is whether approximate Centroid-Linkage HAC admits a theoretically and practically efficient parallel algorithm. 
experiment statistical significance: 
 Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments? Answer: [Yes] Justification: The algorithms in this paper are deterministic, and hence the quality evaluations are the same across all runs. For the scalability study, we only report the mean running times since the variance in running times is low. Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)? Answer: [NA] Justification: Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected? Answer: [Yes] Justification: Datasets are specified and the licenses are stated in Section 5 and the datasets section of Appendix E. Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. 
 
 Appendix: We formally define Centroid-Linkage HAC in a cluster-centric way. For a set of points X ⊂ R d , we define cent(X) to be the centroid of X as below. Definition 2 (Centroid). Given a set of points X ⊂ R d , the centroid of X is the point x ∈ R d whose ith coordinate isx[i] = y∈X y[i] |X| .Centroid HAC uses distance between centroids as a linkage function. That is, for clusters X and Y , L(X, Y ) = D(cent(X), cent(Y )).c-approximate Centroid-Linkage HAC is defined as follows. We are given input points P ⊂ R d . First, initialize a set of active clusters C = {{p}} p∈P . Then, until |C| = 1 we do the following: let X, Ŷ ∈ C be a pair of active clusters satisfyingL( X, Ŷ ) ≤ c • min X,Y ∈C L(X, Y )where the min is taken over distinct pairs; merge X and Ŷ by removing them from C and adding X ∪ Ŷ to C. (Non-approximate) Centroid-Linkage HAC is just the above with c = 1.In the interest of completion, we give a proof of Theorem 1 (dynamic ANNS for oblivious updates).We will make use of the well-known algorithmic primitive of locality sensitive hashing. Definition 3 (Locality Sensitive Hashing). A family of hash functions H is called (r, cr, p 1 , p 2 )sensitive if for any p, q ∈ R d we haveLemma 1. Given s ∈ R d and β, ∆ > 0, there exists a covering netQ ⊆ R d such that 1. Small Size: |Q| ≤ (∆/β) O(d log d) 2. Queries: given u ∈ B s (∆), one can compute u ′ ∈ Q in time O(d) such that D(u, u ′ ) ≤ β.Proof. The basic idea is simply to take an evenly spaced hypergrid centered at s of radius ∆. More formally, we let The number of points in Q is trivially the number of offsets y which is, in turn,Q :=    s + β √ d • y : y ∈ Z d ∩ -∆ √ d β , ∆ √ d β2 ∆ • √ d β d = ∆ β O(d log d)as desired.Next we analyze how to query. Given a point u, we letu ′ i := s i + β √ d • √ d β (u i -s i )be u i rounded down to the nearest multiple of β √ d after offsetting by s i . Observe that u ′ = (u ′ 1 , u ′ 2 , . . . u ′ d ) ∈ Q by our assumption that u ∈ B s (∆) and that, furthermore, u ′ can be computed from u in time O(d). Lastly, observe that|u i -u ′ i | ≤ β √ dand so we have that the distance between u and u ′ isD(u, u ′ ) = i (u i -u ′ i ) 2 ≤ i β √ d 2 = β.Lemma 2. Algorithm 1 is a (c, cβ)-approximate nearest-neighbor search data structure with high probability for adaptive updates.Proof. Consider a query on (u, ū) with corresponding u ′ ∈ Q in our covering net (i.e., the u ′ returned by Lemma 1 when input u). We refer to the state of our data structure at the beginning of its for loop (i.e., upon receiving a new update or query) as at the beginning of an iteration.First, observe that a simple argument by induction on our insertions and deletions shows that the point set for which N i is a data structure at the beginning of an iteration is S i and that, furthermore, if S is the set of points which have been inserted but not deleted at a given point in time then {S i } i partitions S. Thus, we havemin i D(u ′ , S i \ ū) = D(u ′ , S \ ū).Furthermore, by the guarantees of Lemma 1 and the fact that all points of S and all query points are always contained in B s (∆), we know thatEQUATIONIt follows thatD(u ′ , S \ ū) ≤ D(u, S \ ū) + O(β)and soEQUATIONFurthermore, notice that since the randomness of N i is chosen after we fix S i , we have by assumption that for a fixed u ′ ∈ Q that N i succeeds for a query on (u ′ , ū) with points S i except with probability at most n -Ω(d log d•log(∆/β)) .On the other hand, by Lemma 1 we know that|Q| ≤ ∆ β O(d log d)so taking a union bound over all points in Q we have that with high probability N i succeeds for S i and every query point in Q except with probability at most 1/poly(n). Union bounding over our O(n)-many N i that we ever instantiate, we get that every N i succeeds for its corresponding S i and every query point of Q except with probability at most 1/poly(n). In other words, with high probability we always haveEQUATIONThus, combining the triangle inequality, Equations 3, 4 and 5, we have that with high probability whenever we query point u, the returned point v satisfiesD(u, v) = min i D(u, v i ) ≤ D(u, u ′ ) + min i D(u ′ , v i ) ≤ O(β) + min i D(u ′ , v i ) ≤ O(β) + c • min i •D(u ′ , S i \ ū) ≤ c • D(u, S \ ū) + O(c • β)Choosing our hidden constant appropriately, we get that with high probability our data structure is indeed (c, cβ)-approximate.• T I where T I is the insertion time of the input oblivious dynamic ANNS data structure.Proof. Let n insert ≤ n be the total number of insertions up to some point in time. Each time we instantiate an N i+1 on an S i+1 , it consists of at most 2 i+2 -many points and so by assumption takes us timeT I • 2 (i+2)time to instantiate. On the other hand, the total number of times we can instantiate N i+1 is at most n insert /2 i (since each time we instantiate it we add a new set of at least 2 i points to S i+1 and points only move from S j s to S j+1 s). It follows that the total time to instantiate N i+1 is at mostT I • 2 (i+2) • n insert 2 i = O(T I • n insert ).Applying, the fact that i ≤ log n, we get that the total time for all insertions isO(log n) • T I • n insert .Dividing by our n insert -many insertions, we get an amortized insertion time of O(log n) • T I , as desired.Theorem We now prove that Algorithm 3 gives a (c(1 + ϵ))-approximate HAC where c is the approximation value inherited from the approximate nearest neighbor data structure. Note that a smaller ϵ will give a more accurate HAC but will increase the runtime of the algorithm. Lemma 4. Algorithm 3 gives a c(1 + ϵ)-approximate Centroid-Linkage HAC.Proof. The idea of the proof is to observe that one of the two endpoints of the closest pair will always store an approximate of the minimum distance in Q. In that spirit, suppose that right before a centroid is dequeued, the shortest distance between active centroids is l OPT . Then there must exist active centroids x 1 and x 2 such that D(x 1 , x 2 ) = l OPT . Assume without loss of generality that x 1 was queued before x 2 and x 2 was queued with distance l 2 . Then, when x 2 was queued, x 1 was active sol 2 ≤ c • l OPT .Thus, when a centroid, x, is dequeued, its distance l x obeys the inequalityl x ≤ l 2 ≤ c • l OPT .Then, if x is merged during the while loop, it is merged with a centroid of distance at most(1 + ϵ)l x ≤ c(1 + ϵ)l OPTas desired. If a centriod gets merged without being dequeued, that means it merged with another identical centroid which is consistent with 1-approximate HAC.We now prove the runtime of the algorithm. Say that the maximum pairwise distance between points in P is ∆ and the minimum is δ. We make the following observations about the geometry of the centriods throughout the runtime of Algorithm 3. Observation 1. Consider any point p ∈ P . Throughout the runtime of Algorithm 3, the distance from any centroid to p is at most ∆.Proof. Consider a cluster Q ⊆ P . Since for each q ∈ Q we have D(p, q) ≤ ∆, it follows that D(p, cent(Q)) ≤ ∆.Observation 2. When running Algorithm 3, at any point in time there can only ever be a single centroid queued with distance less than δ.Proof. Initially, every centroid is queued with distance at least δ. Thus, if a centroid is ever queued with a smaller distance it will immediately be dequeued and merged.Lemma 5. Algorithm 3 runs in time O(n• T I + n • T D + n • T Q • log 1+ϵ ( 2∆ δ ) + n log(n) log 1+ϵ ( 2∆ δ )).Proof. The idea of the proof is that each time we dequeue a centroid we either merge it or we increase the distance associated with it by a multiplicative factor of 1 + ϵ which can only happen log 1+ϵ (2∆/δ) times for each of our O(n)-many centroids.Each centroid is inserted into N at most once, either during the initial construction or in Algorithm 2, and deleted at most once. Since there are 2n -1 centroids throughout the runtime, the total time spent on updates to N is O(n• T I + n • T D ).At any given point during the algorithm, a centroid can be in Q at most once. Since only 2n -1 centroids are ever created, there are at most 2n -1 elements in Q at any time. Thus, a single queue or dequeue operation takes O(log(n)) time. In the initial for loop we queue n elements to Q and make n queries to N taking time O(n log(n) + n • T Q ). Each time algorithm 2 gets called it queues one element to Q and makes up to 4 updates to C. Since Algorithm 2 is called at most n -1 times, the total time spent on merges, not including updates to N , is O(n• T Q + n log(n)).Lastly, we consider how many times a centroid x can be dequeued on Line 8. Suppose x is queued with distance l i the ith time it is queued. If l 1 < δ then by Observation 2, x will be immediately dequeued and merged so x is only queued once. Thus, we assume l 1 ≥ δ. If x is dequeued with distance l i then either x is merged and will not be queued again or it is queued again with distancel i+1 > (1 + ϵ)l i > (1 + ϵ) i l 1 .Since l 1 ≥ δ by assumption and l i ≤ 2∆ for all i ≥ 1 by Observation 1, it follows that x can only be queued, and therefore dequeued, log 1+ϵ ( 2∆ δ ) times. Each time, not including merges, there can be one dequeue from Q, one queue to Q, and one query to N . Since there are 2n -1 centroids created by the algorithm, the while loop (not including merges) runs in timeO(n log(n) log 1+ϵ ( 2∆ δ ) + n • T Q • log 1+ϵ ( 2∆ δ )).Combining the previous two lemmas together we get the following. Theorem 4 (Reduction of Centroid HAC to Dynamic ANN for Adaptive Updates). Suppose we are given a set of n points P ⊂ R d , c > 1, ϵ > 0, and lower and upper bounds on the minimum and maximum pairwise distance in P of δ and ∆ respectively. Suppose we are also given a data structure that is a dynamic c-approximate NNS data structure for all queries in Algorithm 3 that works for adaptive updates with insertion, deletion and query times of T I , T D , and T Q .Then there exists an algorithm for c(1 + ϵ)-approximate Centroid-Linkage HAC on n points that runs in time Õ(n• (T I + T D + T Q )) assuming ∆ δ ≤ poly(n) and 1 ϵ ≤ polylog(n).Proof. By Lemma 4 and Lemma 5 we get that Algorithm 3 gives c(1+ ϵ)-approximate HAC in time O(n • T I + n • T D + n • T Q • log 1+ϵ ( 2∆ δ ) + n log(n) log 1+ϵ ( 2∆ δ )). It is only left to show that log 1+ϵ ( 2∆ δ ) ≤ polylog(n).Changing the base of the logarithm we getlog 1+ϵ 2∆ δ = log( 2∆ δ ) log(1 + ϵ) .We have that log( 2∆ δ ) ≤ polylog(n) by the assumption that ∆ δ ≤ poly(n). By the inequality log(1 + x) ≥ x 1+x for all x > -1 we get that if ϵ ≤ 1 1 log(1 + ϵ) ≤ 1 + ϵ ϵ ≤ 2 ϵ = O 1 ϵ ≤ polylog(n)where the final inequality comes from the assumption that 1 ϵ ≤ polylog(n). If ϵ > 1 then1 log(1 + ϵ) ≤ 1 log(2)≤ polylog(n).Combining the above inequalities and the above-stated runtime gives our final runtime.We use Theorem 4 with our results from Section 3 to get our main result. First we prove the following lemma which says that our ANNS that includes an additive term functions as a purely multiplicative ANNS for our algorithm.λ ∈ (1, c), if N is a c λ , δ(λ-1)(1+c)λ -approximate NNS data structure that works for adaptive updates then N is a c-approximate NNS data structure for all queries in Algorithm 3.Proof. Let C be the set of active centroids at the time any given query takes place. We will say a point x is close to the set C if D(x, C) < δ c . We show that if x is not close to C, then N query (x, x) = y will give us the multiplicative result we desire:EQUATION)We now prove the lemma by induction. The base case is the set of neighbors calculated for the initial set of points P . Each p ∈ P is not close to C \ p so by Eq. ( 6),D(p, N query (p, p)) < c • D(p, C \ p)as desired. Now for the induction step, consider some query for a centroid x that has just been dequeued. For x, y ∈ C, say y is a c-approximate nearest neighbor ofx if D(x, y) ≤ c • D(x, C \ x).Assume by induction that all queries up to this point have returned c-approximate nearest neighbors. We know by Observation 2 and the inductive hypothesis that either no points are close to C or x is the lone point that is close to C. In the former case, N returns a c-approximate nearest neighbor by Eq. ( 6). Now consider the latter case. Let y be a centroid such that D(x, y) > c • D(x, C \ x). We will show that N query (x, x) ̸ = y.Define L = D(x, y) + D(x, C \ x). By the triangle inequality and because x is the only point that is close to C, we have L ≥ δ c . Define p y = D(x, y)/L and p S = D(x, C \ x)/L to be the percentage of the distance L from D(x, y) and D(x, C \ x) respectively. It is enough for us to show thatc λ (p S • L) + δ(λ -1) (1 + c)λ < p y • Lfor all L ≥ δ c . Since this equation is linear in L, if it is true for any such L then it is true for L = δ c . Since p y > c • p S , we have p S < 1 1+c and p y > c 1+c . Plugging in L = δ c we getc λ p S • L + δ(λ -1) (1 + c)λ < c λ 1 1 + c • δ c + δ(λ -1) (1 + c)λ = δ + δ(λ -1) (1 + c)λ = δ (1 + c) < p y • δ c = p y • L as desired.Theorem 5 (Centroid HAC Algorithm). For n points in R d , there exists a randomized algorithm that given any c > 1.01 and lower and upper bounds on the minimum and maximum pairwise distance in P of δ and ∆ respectively, computes a c-approximate Centroid-Linkage HAC with high probability.It has runtime O(n 1+O(1/c 2 ) d 2 ) assuming c, d, ∆ δ ≤ poly(n).Proof. Let ϵ = 0.001 and define ĉ = c 1+ϵ > 1.005. For any λ ∈ (1, ĉ), Theorem 3 says we can construct a ĉ λ , δ(λ-1) (1+ĉ)λ -approximate nearest neighbor data structure, N , with query, deletion, andamortized insertion time n O(λ 2 /ĉ 2 ) log ∆(1+ĉ)λ δ(λ-1) d 2 .Note that we moved the big O in Theorem 3 from the approximation guarantee to the runtime by scaling by some constant. By Lemma 6, N is a ĉ-approximate ANNS data structure for all queries in Algorithm 3. By Theorem 4 and our definition of ĉ, we get c-approximate centroid HAC by running Algorithm 3 with N . The runtime is Õ(n 1+O(λ 2 /ĉ 2 ) log ∆(1+ĉ)λ δ(λ-1) d 2 ). Letting λ = 1.005 and by the assumptions that ∆ δ , c ≤ poly(n), all log terms are dominated by the big O in the exponent and we get the runtimeO(n 1+O(1/ĉ 2 ) d 2 ) = O(n 1+O(1/c 2 ) d 2 ).In Algorithm 3, in the worst case we may have to queue a centroid log 1+ϵ ( 2∆ δ ) times. However, in practice a centroid may only be queued a small number of times. Say a tuple (l, x, y) in our priority queue becomes stale if the centroid y is merged. Let Γ be the average number of tuples in our priority queue that become stale with each merge. Our experiments summarized in Table 4 suggest that Γ is a small constant in practice with Γ < 1 for all experiments. The following lemma shows that in this case we lose the log 1+ϵ ( 2∆ δ ) in the runtime. Lemma 7. If the average number of queued element made stale for each merge in Algorithm 3 is Γ, then Algorithm 3 runs in time O(n• T I + n • T D + n • T Q • Γ + n log(n) • Γ).Proof. There are a total of n -1 merges during the runtime of Algorithm 3 so the total number of tuples that are made stale is (n -1)Γ. There are an additional n -1 tuples that get merged. Thus, the total number of tuples that get dequeued throughout the runtime of the algorithm is O(Γ). The rest of the proof is the same as that of Lemma 5 except we pay the cost of de-queuing on Line 8O(n • Γ) times instead of n • log 1+ϵ ( 2∆ δ ).Experimental Setup: We run our experiments on a 96-core Dell PowerEdge R940 (with two-way hyperthreading) machine with 4×2.4GHz Intel 24-core 8160 Xeon processors (with 33MB L3 cache) and 1.5TB of main memory. Our programs are implemented in C++ and compiled with the g++ compiler (version 11.4) with -O3 flag. 2 . mnist [13] is a standard machine learning dataset that consists of 28 × 28 dimensional grayscale images of digits between 0 and 9. Here, each digit corresponds to a cluster. The birds dataset contains 224 × 224 × 3 dimensional images of 525 species of birds. As done in previous works [42] , we pass each image through ConvNet [27] to obtain an embedding. Algorithms Implemented We implement the approximate Centroid HAC algorithm described in Algorithm 3, using the dynamic ANNS implementation based on DiskANN [38, 18] . We denote this algorithm as (1 + ϵ)-Centroid HAC or Centroid ϵ . We emphasize that although we refer to this as (1 + ϵ)-Centroid HAC, it is not necessarily a true (1 + ϵ)-approximate algorithm since we are using a heuristic ANNS data structure with no theoretical guarantees on the approximation factor. The use of (1 + ϵ) here captures only the loss from merging a "near-optimal" pair of clusters, as detailed in Algorithm 3.As our main baseline, we consider the optimized implementation of exact Centroid HAC from the fastcluster4 package. This implementation requires quadratic space since it maintains the distance matrix, and is a completely sequential algorithm. We also implement an efficient version of exact Centroid HAC based on our framework (i.e. setting ϵ = 0 and using exact NNS queries) which has the benefit of using only linear space and supports parallelism in NNS queries.We also implement a bucket-based version of approximate Centroid HAC as a baseline, based on the approach of [1] . The algorithm of [1] runs in rounds processing pairs of clusters whose distance is within the threshold defined by that round; the threshold value scales by a constant factor, resulting in logarithmic (in aspect ratio) number of rounds. In each round, they consider a pair of clusters that are within the threshold distance away, and merge them. For the linkage criteria considered in [1] , the distances between two clusters can only go up as a result of two clusters merging. However, this is not true for Centroid HAC, as discussed in Section 1.Here, we observe that when two clusters merge in a round, the only distances affected are distances of other clusters to the new merged cluster. Thus, we can repeatedly check if the nearest neighbor of this cluster is still within the threshold and merge with it in that case. Thus, in a nutshell, the algorithm filters clusters whose nearest neighbors are within the current threshold in each round and merges clusters by the approach mentioned above. This algorithm requires many redundant NNS queries, as we will see during the experimental evaluation.The output of a hierarchical clustering algorithm is typically a tree of clusters, called a dendrogram, that summarizes all the merges performed by the algorithm. The leaves correspond to the points in the dataset, and each internal node represents the cluster formed by taking all the leaves in its subtree.There is a cost associated to each internal node denoting the cost incurred when merging two clusters to form the cluster associated to that node. Given a threshold value, a clustering is obtained by cutting the dendrogram at certain internal nodes whose associated cost is greater than the threshold. In our experiments, the ARI and NMI scores are calculated by taking the best score obtained over the clusterings formed by considering all possible thresholded-cuts to the output dendrogram. For larger datasets, we consider cuts only at threshold values scaled at a constant factor (i.e. logarithmic-many cuts).Table 3 contains the results obtained for our quality evaluations. We also plot the scores obtained by approximate Centroid HAC as a function of ϵ; see Figure 6 , Figure 7 , Figure 8 , Figure 9 , and Figure 10 .We observe that the quality of the clustering produced by approximate Centroid HAC is generally comparable to that of the exact Centroid HAC algorithm. On average, the quality is slightly better in some cases, but we attribute this to noise. In particular, we observe that for the value of ϵ = 0.1, we consistently get comparable quality to that of exact centroid: the ARI and NMI scores are on average within a factor of 7% and 2% to that of exact Centroid HAC, respectively. We also obtained good dendrogram purity score and Dasgupta cost as well, with values within 0.3% and 0.03%, respectively.Unlike other HAC algorithms, the costs of merging clusters for centroid HAC is not monotone (i.e., non-decreasing). The monotonicity property of HAC algorithms is useful when producing the clusterings for applications from the output dendrgoram: typically the dendrogram is cut at certain thresholds, and the disconnected subtrees obtained are flattened to compute the clusters. Here, cutting the dendrogram for a given threshold returns a well-defined and unique clustering. However, this is not the case with Centroid HAC. Indeed, for a given threshold value, we can obtain different clusterings depending on how we define to cut the dendrogram. For e.g., if two nodes U and V exists in the output dendrogram with associated merge costs (i.e., cost incurred at the time of merging the two clusters to form that node) of U being greater than V , and V is an ancestor of U (which is a possibility for Centroid HAC). Then, given a threshold value cost(V ) < τ < cost(V ), we could either cut at V (or above), or stop somewhere below U .We would like to characterize and study this non-monotonicity of Centroid HAC in real-world datasets.For this, we define the following notion: for nodes U, V in the output dendrogram, the pair (U, V ) is called a δ-inversion if V is an ancestor of U and cost(U ) ≥ (1 + δ)cost(V ), where cost(U ) is the merge cost associated to node U .We calculate the number of such inversions incurred by exact Centroid HAC, and compare with the inversions incurred by our approximate algorithm. Figure 5 shows the number of δ-inversions incurred on the benchmark datasets. We observe that for small values of ϵ, approximate Centroid HAC has similar number of δ inversions as that of exact Centroid HAC.Figures 4a and 4b shows the running times on varying slices of the SIFT-1M dataset using one thread and 192 parallel threads, respectively. For the approximate methods, we considered ϵ = 0.1. Figure 4c further compares the heap and bucket based approaches as a function of ϵ.We observe that the bucket based approach is very slow for small values of ϵ due to many redundant nearest-neighbor computations. Yet, at ϵ = 0.1, both algorithms, with 192 threads, obtain speed-ups of up to 175× with an average speed-up of 22× compared to the exact counterparts. However, on a single thread, the bucket based approach fails to scale, while the heap based approach still achieves significant speed-ups of upto 36× with an average speed-up of 5×. Overall, our heap-based approximate Centroid HAC implementation with ϵ = 0.1 demonstrates good quality and scalability, making it a good choice in practice.We compute the number of stale merges (see Section 1) incurred by our algorithm on various benchmark datasets. See Table 4 . We observe that, on average, the number of such stale queries are incurred about 60% of the time compared to the number of actual merges (i.e. number of points). We also compare the number of nearest-neighbor queries made by the heap and bucket based algorithms. As expected, the number of NNS queries made by the bucket-based algorithms is an order of magnitude higher than that of the heap based approach. This is due to many redundant queries performed at each round by the bucket-based algorithm. However, these queries can be performed in parallel, thus resulting in good speed-ups for this algorithm when run on many cores. This assumption is inherited from the current state of the art in the ANNS area, where the fastest known methods (including the one that we evaluate) do not come with strong theoretical guarantees. Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.38th Conference on Neural Information Processing Systems (NeurIPS 2024).1We note that the notion of ANNS that we use is, in fact, slightly more general than this. 2Other papers sometimes refer to this property as working against an adaptive adversary.4https://pypi.org/project/fastcluster/ 


Decision: This paper gives faster algorithms for centroid linkage clustering.  There has been a general interest in speeding up notoriously unscalable hierarchical clustering methods. The paper gives and algorithm with strong theoretical and empirical evidence. The reviewers agreed these results will be of interest to the clustering community at NeurIPS.

Reviews: review: paper_summary: paper_summary: Paper studied Hierarchical clustering algorithm. HAC is popular method where n points are initially singleton clusters (leaves). At each step, two clusters are merged, and the process is continued until one giant cluster of all n points (root) remains. THe merging process gives a hierarchy or nested clusterings. Depending on which clusters are merged ,the algorithm differs.

Common methods are single linkage where the distance between clusters is in terms of min x in C_1 y in C_2 d(x,y).
Another method which is focus of current undertaking is centroid linkage, which is d(C_1, C_2) = dist(mu_1, mu_2) where mu_i is centroid vector of cluster i. 
Technical point of view: Interestingly, after merging two clusters,the best inter-cluster distance can sometimes decrease, so it is not monotone increasing over time as in single linkage.

Since many lower bounds exist showing that we need theta(n^2) time to do HAC, paper studies approximate HAC, where we merge any apprxoximately closest clusters instead of exact closest clusters.

To do this, authors capture the HAC problem as that of approximate nearest neighbors over dynamic dataset. Each cluster is a point in the ANN index, and when 2 clusters merge, we capture by 2 deletions and 1 insertion. Then using a priority queue where each cluster maintains its closest distance to other clusters, the authors show how to maintain this PQ online, and use this to quickly determine which clusters to merge. 

To make it all work they devise provable new algorithms for dynamic ANN algorithms using LSH which can handle adaptive inserts/deletes (which is needed for this problem). Apparently, prior work only handles non-adaptive.
strength_weakness: strength_weakness: Datasets chosen for evaluation can be more diverse.

The theoretical result can be presented with more rigor -- For example, I could not follow what objective function they are presenting a c-apprxoximation for the HAC problem in Theorem 5.

The technical reason/idea behind why they are able to convert a non-adaptive ANN structure to adaptive using the power-of-two size indices is not explained well, seems like the crux of the paper but missing in detail. Is there precedent for using such ideas to make non-adaptive inputs into adaptive inputs?
questions: questions: In Line 140, don't you want the weighted average  and not the weighted sum?

What is the objective function you use to compare the quality of HAC soluton and prove approximation factor gaurantee?

Why not use more modern-day neural-network inspired clustering datasets for benchmarking, as it will be an important use-case going forward (like the OpenAi embeddings or more BERT-based embeddings).

What is SOTA known for other linkage functions in HAC like single linkage and average linkage? Do your ideas extend there also? Which linkage is most used in practice?

You cite [2, 14] for dynamic non-adaptive ANN, However I cant find the corresponding results in either. Please be more specific.
limitations: limitations: Would be good to mention a more detailed limitations section apart from those mentioned in conclusion / open problems.
paper_summary: paper_summary: This paper considers the design of approximate versions of the centroid linkage method for hierarchical clustering.  The exact version of centroid linkage requires $\Theta(n^2)$ time, but it could be possible to do better if some relaxation is allowed.  This paper shows how to get a $c$-approximate centroid linkage clustering in sub-quadratic time by using a dynamic approximate near-neighbor search data structure that has guarantees against adaptive updates.  This is non-trivial since the distance between centroids is non-monotone: after merging a pair of clusters, the resulting distance between centroids of other clusters could go down.  To complete the result, they construct such a data structure that supports insertion, deletion, and queries against adaptive updates with sub-linear running time.

An extensive experimental evaluation is performed to test a practical variant of the proposed algorithm on standard benchmark datasets against the exact centroid linkage method provided by the fastcluster library, a standard library for hierarchical clustering algorithms).  The proposed method either shows good performance or minimal loss against the exact method for a variety of clustering metrics and significant speedup over the exact method.
strength_weakness: strength_weakness: - The experimental evaluation isn't for the algorithm that is proposed theoretically, but instead for a variant with several practical modifications (e.g., using graph-based ANNS which may not have the same theoretical guarantees as LSH-based ones).

- For the experiments, the authors also consider running the algorithm with 192 cores.  For this case, it seems that fastcluster cannot make use of the additional cores, which seems to be an unfair comparison.  However, they do use an exact version of the new algorithm which does make use of the additional cores, providing a more fair comparison.
questions: questions: Minor comment - the usage of $Q$ for both covering net and priority queue is somewhat confusing.  Please consider modifying this.
limitations: limitations: Limitations are made clear through a clear description of assumptions.  I don't see potential negative societal impact from this work since it develops more efficient versions of already widely-used algorithms.
paper_summary: paper_summary: The paper deals with centroid-linkage agglomerative clustering (centroid HAC), where the distance between two clusters is the distance between their centers. It presents a subquadratic time algorithm (approximate centroid HAC) that, instead of requiring the two closest clusters to be merged at each step, it allows for any two clusters to be merged if their distance is within a factor of c of that of the two closest clusters.
strength_weakness: strength_weakness: W1. It is confusing that the ANNS algorithm used in the experiments is not the same as the one proposed and analyzed in Section 3 of the paper (Algorithm 1). Since the focus of the paper is on clustering, this issue creates confusion regarding the actual contribution of this work. Perhaps the ANNS algorithm (contribution 2) is better to be presented in a forum related to data structures.

W2. It is not clear if the ‘approximate centroid HAC’ idea is novel or not (independently of the ANNS algorithm used in the implementation).
questions: questions: Q1, Q2. See comments W1, W2 above.

Q3. How does the method compare to other approximate HAC methods in terms of NMI? In Figure 3  running time comparison is provided. In Figure 5 some results are provided on the Iris dataset, but I cannot see any general conclusions.

Q4. It is meaningless to compare using small datasets, since the optimal solution can be quickly obtained. More experiments using big would add value to the paper.

Q5. Is it possible to have a (rough) estimate of the increase in execution time as the value of \epsilon decreases?
limitations: limitations: The authors discuss the limitations. I do not see any potential negative societal impact of this work.
paper_summary: paper_summary: The authors give an algorithm that approximates a centroid linkage clustering. Their algorithm is fast both in theory and in practice. The algorithm is based on a novel fully dynamic data structure for nearest neighbors which is of independent interest.
strength_weakness: strength_weakness: The techniques are very much geared towards the centroid linkage function. It could be interesting to mention how much of it could potentially generalize to other linkage functions and where the bottlenecks for full generalization lie.
questions: questions: -Do you have any intuition as to why "worse" approximation factors often seem to give better clustering results (as your Table 1 indicates)?
limitations: limitations: The limitations are addressed in the checklist.

Comments: comment: We thank the reviewer for increasing their score. We will make an effort to add a detailed discussion on this gap between theory and practice in the revised version, which we hope will address their concern.
Thank you again to the reviewer for their suggestions. We will further investigate this direction by comparing with single-linkage clustering and include our findings in the revised version.
I thank the authors for the rebuttal and their efforts to respond to my comments. I am still concerned about the incoherence issue: the ANNS algorithm proposed and analyzed is different from the one used in the experiments. For a reader who wishes to apply the method, the proposed ANNS algorithm is actually irrelevant. I will increase my score from 4 to 6.
I would like to thank the authors for their rebuttal. Their answer to my second question seems plausible, and it would perhaps be interesting to compare the quality of the clusterings to the results of other methods, such as single-linkage to investigate the proposed reason.
Thank you for your response and addressing my questions.  My overall evaluation remains the same.