Uniform Concentration Bounds toward a Unified Framework for Robust Clustering

Clustering is a fundamental task in unsupervised learning, which seeks to discover naturally occurring groups within a dataset. Among a plethora of algorithms in the literature, center-based methods remain widely popular, with $k$-means (MacQueen,, 1967; Lloyd,, 1982) as the most prominent example. Given $n$ data points $\mathcal{X}=\{\boldsymbol{X}_{i}:i=1,\dots,n\}\subset\mathop{\mathbb{R}}\nolimits^{p}$, $k$-means seeks to partition the data into $k$ mutually exclusive and exhaustive groups by minimizing the within cluster variance. Representing the cluster centroids $\boldsymbol{\Theta}=\{\boldsymbol{\theta}_{1},\dots,\boldsymbol{\theta}_{k}\}\subset\mathop{\mathbb{R}}\nolimits^{p}$, the $k$-means problem is formulated as the minimization of the objective function

where $d(\cdot,\cdot)$ is a dissimilarity measure. Taking the squared Euclidean distance $d(\boldsymbol{x},\boldsymbol{y})=\|\boldsymbol{x}-\boldsymbol{y}\|_{2}^{2}$ yields the classical $k$-means formulation.

Clustering via $k$-means is well-documented to be sensitive to initialization (Vassilvitskii and Arthur,, 2006), prone to getting stuck in poor local optima (Zhang et al.,, 1999; Xu and Lange,, 2019; Chakraborty et al.,, 2020), and fragile against linearly non-separable clusters (Ng et al.,, 2001) or in the presence of even a single outlying observation (Klochkov et al.,, 2020). Researchers continue to tackle these drawbacks of $k$-means while preserving its simplicity and interpretability (Banerjee et al.,, 2005; Aloise et al.,, 2009; Teboulle,, 2007; Ostrovsky et al.,, 2013; Zhang et al.,, 2020; Telgarsky and Dasgupta,, 2012). Although often successful in practice, few of these approaches are grounded in rigorous theory or provide finite-sample statistical guarantees. The available consistency results are mostly asymptotic in nature (Chakraborty and Das,, 2020; Chakraborty et al.,, 2020; Terada,, 2014, 2015), lacking convergence rates or operating under restrictive assumptions on the relation between $p$ and $n$ (Paul et al., 2021a, ). For example, the recent large sample analysis (Paul et al., 2021a, ) of the hard Bregman $k$-means algorithm (Banerjee et al.,, 2005) obtains an asymptotic error rate of $O(\sqrt{\log n/n})$ under restrictive assumptions on the relation between $n$ and $p$. The approach by Teboulle, (2007) focuses on a unified framework from an optimization perspective, while Balcan et al., (2008); Telgarsky and Dasgupta, (2012) focus on different divergence-based methods to better understand the underlying feature-space.

The presence of outliers only further complicates matters. Outlying data are common in real applications, but many of the aforementioned approaches are fragile to deviations from the assumed data generating mechanism. On the other hand, recent work on robust clustering methods (Deshpande et al.,, 2020; Fischer et al.,, 2020; Klochkov et al.,, 2020) do not integrate the practical advances surveyed above, and similarly lack finite-sample guarantees. To bridge this gap, the Median of Means (MoM) literature provides a promising and attractive framework to robustify center-based clustering methods against outliers. MoM estimators are not only insensitive to outliers, but are also equipped with exponential concentration results under the mild condition of finite variance (Lugosi and Mendelson,, 2019; Lecué and Lerasle,, 2020; Bartlett et al.,, 2002; Lerasle,, 2019; Laforgue et al.,, 2019). Recently, near-optimal results for mean estimation (Minsker,, 2018), classification (Lecué et al.,, 2020), regression (Mathieu and Minsker,, 2021; Lugosi and Mendelson,, 2019), clustering (Klochkov et al.,, 2020; Brunet-Saumard et al.,, 2022), bandits (Bubeck et al.,, 2013) and optimal transport (Staerman et al.,, 2021) have been established from this perspective.

Under the MoM lens, we propose a unified framework for robust center-based clustering. Our treatment considers a family of Bregman loss functions not restricted to the classical squared Euclidean loss. We explore the exact sample error bounds for a general class of algorithms under mild regularity assumptions that apply to popular existing approaches, which we show to be special cases. The proposed framework allows for the data to be divided into two categories: the set of inliers ($\mathcal{I}$) and the set of outliers ($\mathcal{O}$). The inliers are assumed to be independent and identically distributed (i.i.d.), while absolutely no assumption is required on $\mathcal{O}$, allowing outliers to be unboundedly large, dependent on each other, sampled from a heavy-tailed distribution and so on. Our contributions within the MoM framework make use of Rademacher complexities (Bartlett and Mendelson,, 2002; Bartlett et al.,, 2005) and symmetrization arguments (Vapnik,, 2013), powerful tools that often find use in the empirical process literature but, in our view, are underexplored in the context of robust clustering.

The paper is organized as follows. In Section 2, we identify a general centroid-based clustering framework which encompasses $k$-means (MacQueen,, 1967), $k$-harmonic means (Zhang et al.,, 1999), and power $k$-means (Xu and Lange,, 2019) as special cases, to name a few. We show how this framework is made robust via Median of Means estimation, yielding an array of center-based robust clustering methods. Within this framework, we derive uniform deviation bounds and concentration inequalities under standard regularity conditions through bounding Rademacher complexity by metric entropy via Dudley’s chaining argument in Section 3. The analysis newly reveals the convergence rate for popularly used clustering methods such as $k$-harmonic means and power $k$-means, matching the known rate results for $k$-means, and elegantly carries over to their MoM counterparts. We then implement and empirically assess the resulting algorithms through simulated and real data experiments. In particular, we find that power $k$-means (Xu and Lange,, 2019) under the MoM paradigm outperforms the state-of-the-art in the presence of outliers.

We consider the problem of partitioning a set of $n$ data points $\mathcal{X}=\{\boldsymbol{X}_{i}:i=1,\dots,n\}\subset\mathop{\mathbb{R}}\nolimits^{p}$ into $k$ mutually exclusive clusters. In a center-based clustering framework, we represent the $j^{\text{th}}$ cluster by its centroid $\boldsymbol{\theta}_{j}\in\mathop{\mathbb{R}}\nolimits^{p}$ for each $j\in\{1,\dots,k\}$. To quantify the notion of “closeness", we allow the dissimilarity measure to be any Bregman divergence. Recall any differentiable, convex function $\phi:\mathop{\mathbb{R}}\nolimits^{p}\to\mathop{\mathbb{R}}\nolimits$ generates the Bregman divergence $d_{\phi}:\mathop{\mathbb{R}}\nolimits^{p}\times\mathop{\mathbb{R}}\nolimits^{p}\to\mathop{\mathbb{R}}\nolimits_{\geq 0}$ ($\mathop{\mathbb{R}}\nolimits_{\geq 0}$ denoting the set of non-negative reals) defined as

For instance, $\phi(\boldsymbol{u})=\|\boldsymbol{u}\|_{2}^{2}$ generates the Euclidean distance. Without loss of generality, one may assume $\phi(\mathbf{0})=\nabla\phi(\mathbf{0})=0.$ In this paradigm, clustering is achieved by minimizing the objective

Here $\Psi_{\boldsymbol{\alpha}}:\mathop{\mathbb{R}}\nolimits^{k}_{\geq 0}\to\mathop{\mathbb{R}}\nolimits_{\geq 0}$ is a component-wise non-decreasing function (such as a generalized mean) of the dissimilarities $\{d_{\phi}(\boldsymbol{X},\boldsymbol{\Theta}_{j})\}_{j=1}^{k}$ which satisfies $\Psi(\mathbf{0})=0$. The hyperparameter $\boldsymbol{\alpha}\in\mathcal{A}\subseteq\mathop{\mathbb{R}}\nolimits^{q}$ is specified by the user, and we will additionally assume that $\Psi$ is Lipschitz. For intuition, we begin by showing how this setup includes several popular clustering methods.

Here we analyze properties of the proposed framework (1), with complete details of all proofs in the Appendix. Denoting $\mathcal{M}$ the set of all probability measures $P$ with support on $[-M,M]^{p}$, i.e. $P([-M,M]^{p})=1,\,\forall P\in\mathcal{M}$, we first make standard assumptions that the data are i.i.d. with bounded components (Ben-David,, 2007; Chakraborty et al.,, 2020; Paul et al., 2021a, ).

Let $P_{n}$ be the empirical distribution based of the data $\boldsymbol{X}_{1}\dots,\boldsymbol{X}_{n}$. That is, for any Borel set $A$, $P_{n}(A)=\frac{1}{n}\sum_{i=1}^{n}\mathbbm{1}\{\boldsymbol{X}_{i}\in A\}$. For notational simplicity, we write $\mu f:=\int fd\mu$. Appealing to the Strong Law of Large Numbers (SLLN) (Athreya and Lahiri,, 2006), $P_{n}f_{\boldsymbol{\Theta}}\xrightarrow{a.s.}Pf_{\boldsymbol{\Theta}}$. Suppose $\widehat{\boldsymbol{\Theta}}_{n}$ be a (global) minimizer of (1) and $\boldsymbol{\Theta}^{\ast}$ be the global minimizer of $Pf_{\boldsymbol{\Theta}}$. Since the functions $P_{n}f_{\boldsymbol{\Theta}}$ and $Pf_{\boldsymbol{\Theta}}$ are close to each other as $n$ become large, we can expect that their respective minimizers, $\widehat{\boldsymbol{\Theta}}_{n}$ and $\boldsymbol{\Theta}^{\ast}$ are also close to one another. To show that $\widehat{\boldsymbol{\Theta}}_{n}$ converges to $\boldsymbol{\Theta}^{\ast}$ as $n\to\infty$, we consider bounding the uniform deviation $\sup_{\boldsymbol{\Theta}}|P_{n}f_{\boldsymbol{\Theta}}-Pf_{\boldsymbol{\Theta}}|$. Towards establishing such bounds, we will posit two regularity assumptions on $\Psi_{\boldsymbol{\alpha}}(\cdot)$ and $\phi(\cdot)$, beginning with a $\tau_{\boldsymbol{\alpha},k}$-Lipschitz condition on $\Psi_{\boldsymbol{\alpha}}$.

We also assume a weaker form of a standard condition that the gradient $\nabla\phi(\cdot)$ is $H_{p}$-Lipschitz (Telgarsky and Dasgupta,, 2013); unlike their work, note we do not additionally require strong convexity of $\phi$.

These conditions are mild, and can be seen to hold for all of the aforementioned popular special cases. For instance, taking $\phi(\boldsymbol{u})=\|\boldsymbol{u}\|^{2}$ yields the squared Euclidean distances, and we see that $\nabla\phi(\boldsymbol{u})=2\boldsymbol{u}$ so that A3 is satisfied with constant $H_{p}=2$. Now, let $\Psi(\boldsymbol{x})=\min_{1\leq j\leq k}x_{j}$ as in classical $k$-means, and denote $j^{\ast}\in\mathop{\rm argmin}\nolimits_{1\leq j\leq k}y_{j}$. For any vectors $\boldsymbol{x},\boldsymbol{y}\in\mathop{\mathbb{R}}\nolimits_{\geq 0}^{p}$,

Thus, $\Psi$ is clearly non-negative and componentwise non-decreasing, and satisfies A2. Similarly, if we take $\Psi_{s}(\boldsymbol{x})=M_{s}(\boldsymbol{x})$, the conditions are met for power $k$-means: it is again non-negative and non-decreasing in its components, and satisfies A2 with $\tau_{\boldsymbol{\alpha},k}=k^{-1/s}$ due to Beliakov et al., (2010).

To validate and assess our proposed framework, we now turn to an empirical comparison of the proposed and peer clustering methods. We evaluate clustering quality under the Adjusted Rand Index (ARI) (Hubert and Arabie,, 1985), with values ranging between $0$ and $1$ and $1$ denoting a perfect match with the ground truth. Though it is not feasible to exhaustively survey center-based clustering methods, we consider a broad range of competitors, comparing to $k$-means (MacQueen,, 1967), Partition Around Medoids (PAM) (Kaufman and Rousseeuw,, 2009), $k$-medians (Jain,, 2010), Robust kmeans++ (RKMpp) (Deshpande et al.,, 2020), Robust Continuous Clustering (RCC) (Shah and Koltun,, 2017), Bregman Trimmed $k$-means (BTKM) (Fischer et al.,, 2020), MoM $k$-means (MOMKM) (Klochkov et al.,, 2020) and a novel MoM variant of power $k$-means (MOMPKM) implied under the proposed framework. An open-source implementation of the proposed method and code for reproducing the simulation experiments are available at https://github.com/SaptarshiC98/MOMPKM.

We consider two thorough simulated experiments below, with additional simulations and large-scale real data results in the Appendix. While the extended comparisons are omitted for space considerations, they convey the same trends as the studies below. In all settings, we generate data in $p=5$ dimensions, varying the number of clusters and the outlier percentages. True centers are spaced uniformly on a grid with $\theta_{k}=\frac{k-1}{10}$ and observations are drawn from Gaussians around their ground truth centers with variance $0.1$. Because we generate Gaussian data, we focus here on the Euclidean case, and do not consider other Bregman divergences in the present empirical study. In all experiments, we take the number of partitions $L$ to be roughly double the number of outliers, and set the hyperparameters $\eta$ and $\alpha$ to be $1.02$ and $1$ respectively by default. Results are averaged over 20 random restarts, and all competitors are initialized and tuned according to the standard implementations described in their original papers.