A unified approach for cluster-wise and general noise rejection approaches for k-means clustering

Noise rejection in regression analysis and C-means-type clustering

Many machine learning tasks such as regression analysis are formulated in a framework of least mean squares (LMS) proposed by Legendre or Gauss (Legendre, 1805; Gauss, 1809), which minimizes the sum of the squared residuals to fit models to a dataset. However, since the LMS criterion is strongly affected by noise objects and has the lack of robustness, various robust estimation methods have been proposed to reduce the influence of noise objects. Least absolute values (LAV) (Edgeworth, 1887) is a criterion that minimizes the sum of the absolute values of the residuals to reduce the influence of large residuals. M-estimator (Huber, 1964; Huber, 1981) is one of the most widely used robust estimators, which replaces the square function in LMS by a symmetric function with a unique minimum at zero that reduces the influence of large residuals. Least median of squares (LMedS) (Hampel, 1975; Rousseeuw, 1984) minimizes the median of the squared residuals. Least trimmed squares (LTS) (Rousseeuw & Leroy, 1987) minimizes the sum of the squared residuals up to h-th objects in ascending list of residuals.

Since C-means-type clustering methods are generally formulated based on the minimization of the within-cluster sum-of-squared-error, the above-mentioned robust estimation methods are promising approaches to noise in the cluster structure (Kaufmann & Rousseeuw, 1987; Dubes & Jain, 1988). In C-means-type clustering, the distance between object and its nearest cluster center is identified as the residual. Thus, in GNR, objects distant from any cluster center are rejected as noise. For instance, trimmed C-means (TCM; trimmed k-means, TKM) (Cuesta-Albertos, Gordaliza & Matrán, 1997; Garcia-Escudero & Gordaliza, 1999) introduces LTS criterion to HCM. TCM calculates the new cluster center by using objects up to h-th in ascending list of the distances to their nearest cluster centers. As a result, objects more than a certain distance away are rejected as noise. Noise rejection in C-means-type clustering is also well discussed in the context of fuzzy C-means (FCM) (Dunn, 1973; Bezdek, 1981). In noise fuzzy C-means (NFCM) (Davé, 1991; Davé & Krishnapuram, 1997), a single noise cluster is introduced in addition to the intended regular clusters and objects distant from any cluster center are absorbed in the noise cluster. Another approach to noise is CNR. For instance, possibilistic C-means (PCM) (Krishnapuram & Keller, 1993; Krishnapuram & Keller, 1996) considers cluster-wise noise rejection, in which each cluster is extracted independently while rejecting objects distant from its center. The membership values are interpreted as degrees of possibility of the object belonging to clusters. PCM represents typicality as absolute membership to clusters rather than relative membership by eliminating the sum-to-one constraint. Fuzzy possibilistic C-means (FPCM) (Pal, Pal & Bezdek, 1997) uses both relative typicalities (memberships) and absolute typicalities. Possibilistic fuzzy C-means (PFCM) (Pal et al., 2005) is a hybridization of FCM and PCM using both probabilistic memberships of FCM and possibilistic memberships of PCM.

In this study, we show that GNR, CNR, and their combinations are realized by the linear function threshold-based object-cluster assignment in the proposed LiFTCM. The above-mentioned approaches introduce various mechanisms to realize GNR and CNR. In contrast, the linear function threshold-based assignment essentially includes GNR, CNR, and their combinations in compliance with RCM standards without any extra formulation.

Generalized approaches to hard, fuzzy, noise, possibilistic, and rough clustering

Maji & Pal (2007a) proposed rough-fuzzy C-means (RFCM) as a hybrid algorithm of FCM and RCM. RFCM is formulated so that objects in the lower areas have crisp memberships and objects in the boundary areas have FCM-based fuzzy memberships. Furthermore, Maji & Pal (2007b) proposed rough-fuzzy possibilistic C-means (RFPCM)  based on possibilistic fuzzy C-means (PFCM) (Pal et al., 2005). Masson and Denœux proposed evidential C-means (ECM) (Masson & Denoeux, 2008) as one of the evidential clustering (EVCLUS) (Denoeux & Masson, 2003; Denoeux & Masson, 2004) methods based on the Dempster-Shafer theory of belief functions (evidence theory). Evidential clustering considers the basic belief assignment, which indicates the membership (mass of belief) of each object to each subset of clusters with the probabilistic constraints that derive credal partition. Credal partition can represent hard and fuzzy partitions with a noise cluster considering assignments to a singleton and the empty set. Possibilistic and rough partitions are represented by using the plausibility function and the belief function (Denoeux & Kanjanatarakul, 2016).

Although RFCM and RFPCM provide interesting perspectives on the handling of the uncertainty in the boundary area, the object-cluster assignment is different from that of RCM and transform into different types of approach. Although the credal partition in ECM has high expressiveness including hard, noise, possibilistic, and rough clustering, the object-cluster assignment and cluster center calculation of ECM do not boil down to those of RCM. In contrast to the above-mentioned approaches, the formulation of the proposed LiFTCM is fully compliant with RCM standards. This study reveals that RCM itself inherently includes GNR, CNR, and their combinations as well as rough clustering aspects without any extra formulation.

Hard C-means and noise rejection

Let U = {x₁, …, x_i, …, x_n} be a set of n objects, where each object x_i = (x_i1, …, x_ij, …, x_ip)^⊤ is a p-dimensional real feature vector. In C-means-type methods, C(2 ≤ C < n) represents the number of clusters, and C clusters composed of similar objects are extracted from U. Each cluster has a prototypical point (cluster center), which is a p-dimensional vector b_c = (b_c1, …, b_cj, …, b_cp)^⊤. Let u_ci be the degree of belonging of object i to cluster c. Let d_ci = ||x_i − b_c|| be the distance between the cluster center b_c and the object i.

The optimization problem of HCM (MacQueen, 1967) is given by (1)${\displaystyle \text{min.}{J}_{HCM}=\sum _{c=1}^C\sum _{i=1}^n{u}_{ci}{d}_{ci}^2,}$ (2)${\displaystyle \text{s.t.}{u}_{ci}\in \left\{0,1\right\},\forall c,i,}$ (3)${\displaystyle \sum _{c=1}^C{u}_{ci}=1,\forall i.}$ HCM minimizes the total within-cluster sum-of-squared-error (Eq. (1)) under the Boolean domain constraints (Eq. (2)) and the sum-to-one constraints across clusters (Eq. (3)).

HCM first initializes cluster centers and then alternately updates u_ci and b_c until convergence by using the following update rules: (4)${\displaystyle u_{ci}=\left\{\begin{array}{cc}1\hfill & \left(c=\underset{1\le l\le C}{argmin}{d}_{li}\right),\hfill \\ 0\hfill & \left(\text{otherwise}\right),\hfill \end{array}\right.}$ (5)${\displaystyle {\mathit{b}}_c=\frac{\sum _{i=1}^n{u}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{u}_{ci}}.}$ There are various strategies for initializing cluster centers. A naive strategy is to choose C objects as initial cluster centers from U by simple random sampling without replacement. Alternatively, there are strategies that set the initial cluster centers away from each other to reduce initial value dependencies and improve clustering performance, such as KKZ (Katsavounidis, Kuo & Zhang, 1994) and k-means++ (Arthur & Vassilvitskii, 2007).

Generalized rough C-means

In RCM-type methods, which are rough set clustering schemes, membership in the lower, upper, and boundary areas of each cluster represents positive, possible, and uncertain belonging to the cluster, respectively (Lingras & West, 2004; Peters, 2006; Peters et al., 2013; Ubukata, Notsu & Honda, 2017). GRCM is constructed based on a heuristic scheme, not an objective function.

In every iteration, the membership ${\overline{u}}_{ci}$ of object i to the upper area of cluster c is first calculated as follows: (8)${\displaystyle d_i^{min}=\underset{1\le l\le C}{min}{d}_{li},}$ (9)${\displaystyle {\overline{u}}_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le \alpha d_i^{min}+\beta \right),\hfill \\ 0\left(\text{otherwise}\right),\hfill \end{array}\right.}$ where α (α ≥ 1) and β (β ≥ 0) are user-defined parameters that adjust the volume of the upper areas. GRCM assigns each object to the upper area of not only its nearest cluster but also of other relatively nearby clusters using a linear function of the distance to its nearest cluster as a threshold. Larger α and β imply larger clustering roughness and larger overlap of the upper areas of the clusters. Figure 1 shows the linear function threshold T and the allowable range of d_ci (gray area) in GRCM.

The membership ${\underset{¯}{u}}_{ci}$ and ${\hat{u}}_{ci}$ of object i to the lower and boundary areas, respectively, is calculated using ${\overline{u}}_{ci}$ as follows: (10)${\displaystyle {\underset{¯}{u}}_{ci}=\left\{\begin{array}{c}1\left({\overline{u}}_{ci}=1\wedge \sum _{l=1}^C{\overline{u}}_{li}=1\right),\hfill \\ 0\left(\text{otherwise}\right),\hfill \end{array}\right.}$ (11)${\displaystyle {\hat{u}}_{ci}=\left\{\begin{array}{c}1\left({\overline{u}}_{ci}=1\wedge \sum _{l=1}^C{\overline{u}}_{li}\ne 1\right),\hfill \\ 0\left(\text{otherwise}\right)\hfill \end{array}\right.}$ (12)${\displaystyle ={\overline{u}}_{ci}-{\underset{¯}{u}}_{ci}.}$

GRCM represents each cluster by the three regions. Therefore, the new cluster center is determined by the aggregation of the centers of these regions. The cluster center b_c is calculated by the convex combination of the centers of the lower, upper, and boundary areas of the cluster c: (13)${\displaystyle {\mathit{b}}_c=\left\{\begin{array}{c}\frac{\sum _{i=1}^n{\overline{u}}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{\overline{u}}_{ci}}\left(\sum _{i=1}^n{\underset{¯}{u}}_{ci}=0\vee \sum _{i=1}^n{\hat{u}}_{ci}=0\right),\hfill \\ \underset{¯}{w}\frac{\sum _{i=1}^n{\underset{¯}{u}}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{\underset{¯}{u}}_{ci}}+\overline{w}\frac{\sum _{i=1}^n{\overline{u}}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{\overline{u}}_{ci}}+\hat{w}\frac{\sum _{i=1}^n{\hat{u}}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{\hat{u}}_{ci}}\left(\text{otherwise}\right),\hfill \end{array}\right.}$ (14)${\displaystyle \underset{¯}{w},\overline{w},\hat{w}\ge 0,}$ (15)${\displaystyle \underset{¯}{w}+\overline{w}+\hat{w}=1,}$ where $\underset{¯}{w}$, $\overline{w}$, and $\hat{w}$ are user-defined parameters that represent the impact of the centers of the lower, upper, and boundary areas, respectively. Ubukata, Notsu & Honda (2017) suggest $\hat{w}=0$ because the centers of the boundary areas tend to cause instability in the calculations and poor classification performance.

HCM

In HCM, each object is assigned to the cluster whose center is nearest to the object. This assignment can be interpreted as assigning object i to cluster c if d_ci is equal to (or less than) $d_i^{min}$, that is, (16)${\displaystyle u_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le d_i^{min}\right),\hfill \\ 0\left(\text{otherwise}\right).\hfill \end{array}\right.}$ This is the case α = 1 and β = 0 in the linear function threshold $\alpha d_i^{min}+\beta$ for the assignment of upper area in GRCM (Eq. (9)). Figure 2A shows the linear function threshold T and the allowable range of d_ci in HCM. The allowable range is limited to the case $d_{ci}=d_i^{min}$. We note that if there are multiple nearest cluster centers for an object, HCM requires certain tie-breaking rules for satisfying the sum-to-one constraints, such as exclusive assignment based on cluster priority and uniform assignment by distributing the membership, depending on the implementation. However, in the present linear function threshold-based assignment, an object has membership 1 with respect to all its nearest clusters.

The calculation of u_ci in HCM can be represented by that of ${\overline{u}}_{ci}$ in GRCM. The lower and boundary areas are not used in HCM. Thus, the cluster center calculation of HCM is consistent with that of GRCM only using the upper areas, that is, $\overline{w}=1$ in Eq. (13). Therefore, GRCM($\alpha =1,\beta =0,\overline{w}=1$) represents HCM.

GNR

In GNR, a condition that the distance is less than δ is imposed in addition to the threshold-based HCM assignment (Eq. (16)) to reject noise objects over δ distant from any clusters. For each object i to be assigned to the cluster c, d_ci must be equal to (or less than) $d_i^{min}$, and equal to or less than the noise distance δ, that is, (17)${\displaystyle u_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le d_i^{min}\wedge d_{ci}\le \delta \right),\hfill \\ 0\left(\text{otherwise}\right).\hfill \end{array}\right.}$ This assignment can also be approximated using the linear function threshold by setting $\alpha =\frac{\delta -\epsilon }{\delta }$ and β = ε, where ε → +0, that is, (18)${\displaystyle u_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le \frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon \right),\hfill \\ 0\left(\text{otherwise}\right).\hfill \end{array}\right.}$ Equation (18) implies that u_ci = 1 if $d_{ci}\le d_i^{min}$ and d_ci ≤ δ. Thus, Eq. (18) approaches the update rule Eq. (17). In order to show that Eqs. (17) and (18) are equivalent, we show that the condition $d_{ci}\le d_i^{min}\wedge d_{ci}\le \delta$ and the condition $d_{ci}\le \frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon$ are equivalent, under the condition δ > 0 and ε → +0.

proof.

(1) $d_{ci}\le d_i^{min}\wedge d_{ci}\le \delta$ (Assumption)

(2) $d_{ci}\le d_i^{min}$ (Conjunction elimination: (1))

(3) d_ci ≤ δ (Conjunction elimination: (1))

(4) $d_i^{min}\le d_{ci}$ (Definition: (Eq. (8)))

(5) $d_i^{min}\le \delta$ (Transitivity: (3), (4))

(6) $\frac{\epsilon }{\delta }{d}_i^{min}+\frac{\delta -\epsilon }{\delta }{d}_i^{min}\le \frac{\epsilon }{\delta }\delta +\frac{\delta -\epsilon }{\delta }{d}_i^{min}$ (Multiply by $\frac{\epsilon }{\delta }$ and add $\frac{\delta -\epsilon }{\delta }{d}_i^{min}$ in both sides of (5))

(7) $d_i^{min}\le \frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon$ (Deformation: (6))

(8) $d_{ci}\le \frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon$ (Transitivity: (2), (7)) □

proof.

(1) $d_{ci}\le \frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon$ (Assumption)

(2) $d_{ci}\le d_i^{min}$ (From (1) and ε → +0)

(3) $d_i^{min}\le d_{ci}$ (Definition: (Eq. (8)))

(4) $d_{ci}\le \frac{\delta -\epsilon }{\delta }{d}_{ci}+\epsilon$ (From (1), (3))

(5) δd_ci ≤ δd_ci − εd_ci + δε (Multiply by δ in both sides of (4))

(6) d_ci ≤ δ (Deformation: (5))

(7) $d_{ci}\le d_i^{min}\wedge d_{ci}\le \delta$ (Conjunction introduction: (2), (6)) □

Hence, (Eq. (17)) induces (Eq. (18)), and vice versa.

Figure 2B shows the linear function threshold T and the allowable range of d_ci (gray area) in GNR. Since the intersection of the two lines $y=\frac{\delta -\epsilon }{\delta }{d}_i^{min}+\epsilon$ and $y=d_i^{min}$ is (δ, δ), if d_ci > δ, object i is never assigned to cluster c. If $d_i^{min}\le \delta$, the threshold approaches the HCM-based nearest assignment. These characteristics are consistent with those of GNR.

Similar to HCM, in GNR, the cluster centers are calculated only using the upper areas. Therefore, GRCM($\alpha =\frac{\delta -\epsilon }{\delta },\beta =\epsilon ,\overline{w}=1$) represents GNR.

Smooth transition between GNR and CNR by tuning linear function threshold

In reference to the threshold-based assignment of GNR, i.e., Eq. (18), we construct the following rule using a parameter t ∈ [0, δ_c]: (19)${\displaystyle u_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le \frac{{\delta }_c-t}{{\delta }_c}{d}_i^{min}+t\right),\hfill \\ 0\left(\text{otherwise}\right).\hfill \end{array}\right.}$ If t = 0, then Eq. (19) reduces to Eq. (16) of HCM. If t = ε, where ε → +0, then Eq. (19) changes to Eq. (18) of GNR. If t = δ_c, then Eq. (19) reduces to Eq. (7) of CNR. If t ∈ (0, δ_c), then Eq. (19) represents the combinations of GNR and CNR. Thereby, smooth transition between HCM, GNR, and CNR is realized.

Figure 3 shows the linear function threshold T and the allowable range of d_ci (gray area) in the combinations of GNR and CNR. It can be seen that this linear function can transition between the states shown in Fig. 2 by t.

For practical use, we consider the normalized parameter z ∈ [0, 1]. We let $z=\frac{t}{{\delta }_c}\in \left[0,1\right]$ and replace t in Eq. (19) with zδ_c: (20)${\displaystyle u_{ci}=\left\{\begin{array}{c}1\left(d_{ci}\le \left(1-z\right)d_i^{min}+z{\delta }_c\right),\hfill \\ 0\left(\text{otherwise}\right).\hfill \end{array}\right.}$ Then, z = 0 represents HCM, z → +0 represents GNR, z ∈ (0, 1) represents the combinations of GNR and CNR, and z = 1 represents CNR. By Eq. (20), the threshold value is represented by the convex combination of $d_i^{min}$ and δ_c. That is, HCM, GNR, and CNR can be characterized depending on which of $d_i^{min}$ and δ_c is emphasized as the threshold value.

Schematic diagram

Figure 6 is a schematic diagram of the proposal of this study. Representations of HCM, LRCM, PRCM, GRCM, GNR, CNR, and their combinations by the linear function threshold in LiFTCM with the parameters (α, β), and their relationships are shown. (α, β) = (1, 0) is the default state and represents HCM assignment. Increasing α from 1 and β from 0 increases cluster overlap. Simultaneously increasing α and β increases clustering roughness. This shows combinations of LRCM and PRCM, namely, GRCM. LiFTCM gives an interpretation in 0 ≤ α ≤ 1 in addition to GRCM. As proposed in the smooth transition, when the parameter z is increased from 0 to 1, (α, β) transits from (1, 0) to (0, δ_c), namely, from HCM to CNR via GNR. The parameter z has the effect of changing clustering more possibilistic. Cluster overlap in CNR is attributed to the increase in β in LRCM. The direction in which the destination δ_c is lowered is the direction in which noise objects are more rejected.

Cluster center calculation utilizing probabilistic memberships

RCM-type methods have the problem that even if the number of objects in the boundary area is small, they have unnaturally large impacts on the new cluster center compared to the objects in the lower area, because the cluster center is calculated by the convex combination of these areas. To cope with the problem, Peters proposed πPRCM by introducing the cluster center calculation based on the normalized membership of the membership to the upper area, which satisfies the probabilistic constraint (Peters, 2014; Peters, 2015). “π” is an acronym that stands for “Principle of Indifference,” in which the probability is assigned equally by dividing the number of possible clusters. Ubukata et al. proposed πGRCM (Ubukata et al., 2018) based on GRCM. The proposed LiFTCM has almost the same formulation as GRCM except that the condition α ≥ 1 is relaxed to α ≥ 0. Thus, πLiFTCM can be formulated in a similar manner to πGRCM by introducing the following normalized membership ${\tilde{u}}_{ci}$ of the membership to the upper area and the cluster center calculation based on ${\tilde{u}}_{ci}$: (22)${\displaystyle {\tilde{u}}_{ci}=\frac{{\overline{u}}_{ci}}{\sum _{l=1}^C{\overline{u}}_{li}},}$ (23)${\displaystyle {\mathit{b}}_c=\frac{\sum _{i=1}^n{\tilde{u}}_{ci}{\mathit{x}}_i}{\sum _{i=1}^n{\tilde{u}}_{ci}}.}$ Here, attention should be paid to the following cases. In the case of α < 1, that is, in the case of GNR and CNR, since non-belonging of the object to any cluster is handled and thus the denominator ${\sum }_{l=1}^C{\overline{u}}_{li}$ can become zero, it is necessary to set ${\tilde{u}}_{ci}=0$ for all clusters in such cases.