二维块状数据的聚类数量与边界确定

Determination of the cluster number and cluster boundary for two-dimensional blob data

摘要

针对经典K-means算法在划分块状（斑点）数据聚类结构时存在不确定性以及易导致聚类结果失真的问题，提出聚类数量K、聚类中心以及聚类边界的确定方法。基于数据特征的样本均值和标准差，利用先验法则估计数据的最大聚类数K_max，在[1, K_max]范围内计算不同K值所对应的各种聚类性能指标以评估聚类质量。通过寻求性能指标的最大值以确定数据所固有的聚类数量K_best；将K_best作为K-means算法的输入参数，利用期望最大化（E-M）算法迭代估计出相应的聚类中心；以同一聚类中最远的数据点到该聚类中心的欧氏距离为半径，给出各个聚类的圆形边界。基于不同样本数量、不同聚类数的二维斑点随机数据的仿真实验，验证了本文方法的有效性和鲁棒性。

Abstract

To solve the problem that the classic K-means algorithm divides blob dataset into cluster structure with uncertainty and may cause distortion of clustering result, this paper proposes a set of methods for determining the number of clusters K, the cluster centers and cluster boundaries. Based on the sample mean and standard deviation of the data features, the maximum number of clusters K_max of the dataset is estimated using the prior rule. For different numbers within range of K∈[1,K_max], the multiple clustering performance indicators are calculated, and the inherent number of clusters K_best is determined by finding the maximum value of each indicator. With K_best as the input parameter of the K-means algorithm, the corresponding cluster centers are iteratively updated using the expectation-maximization (E-M) algorithm. The circular boundaries for each cluster are defined with the Euclidean distance from the farthest data point in the same cluster to the cluster center as the radius. The validity and robustness of the proposed methods are verified based on the simulation experiment of two-dimensional blob random data with different sample sizes and different number of clusters.

Graphical abstract

关键词

K-means算法 / 聚类性能指标 / 聚类数量 / 聚类中心 / 聚类半径

Key words

K-means algorithm / clustering performance indicator / cluster number / cluster center / cluster radius

0 引言

在大数据时代，数据挖掘^[1]通过抽取数据间关系，辨识与提取数据中隐藏的信息，实现知识发现^[2]并解决应用中的业务问题。数据中隐含的基本信息即为数据分组的结构。对数据结构的辨识技术包括监督分类^[3]与非监督分类。聚类分析（clustering）^[4]是探索性数据分析方法，属于非监督学习范畴。聚类分析能够揭示相似数据样本点组成的簇（集群、聚类）结构和模式^[5]，数据簇划分的结果可以为监督分类器提供数据的类标签（class label）。

K-means聚类算法^[6]计算复杂度低且易于执行，已成为经典算法。该算法由Forgy、MacQueen分别于1965年、1967年独立开发^[7]。以找到密集型数据集的最优分组策略^[8]为目标，众多学者对其展开研究。STEINLY等^[9]利用K-means聚类的额外信息研究算法的局部以及全局最优解，并对其稳定性进行了分析。SELIM等^[10]通过将聚类问题转化为非凸数学规划，为K-means这类算法的有限收敛性提供了严格证明。在可微的条件下，K-means算法可收敛至库恩-塔克点（Kuhn-Tucker point）^[10]。

尽管K-means算法具备诸多优良特性，但存在以下局限：一是迭代优化过程不能保证收敛到全局最优^[7]，聚类结果对随机初始化的聚类中心（质心）敏感^[8]；二是算法执行依赖于预设聚类数量^[8-9]，当数据集未知时，主观设定的簇数将导致聚类结果不可靠；三是数据中的离群值（outlier）和噪声（noise）可能会扭曲数据簇形状^[7]。

为解决上述问题，对经典K-means算法缺陷的改进方案^[11]和算法变体^[12]应运而生。对于算法的随机初始质心问题，利用最远距离选择原则^[13]确定最初的2个聚类中心位置，其余初始质心则采用最小-最大距离规则进行指定与划分^[13]，保证数据簇的可分离性。算法的输入参数（即簇数量）与数据结构的适应程度，是决定聚类结果是否有意义以及是否具有可解释性的关键。在综合考虑数据簇内部聚合性以及簇之间分离性的基础上，聚类有效性评价函数^[14]被提出并应用于K-means算法中，实现了对数据潜在簇数量的估计。对于高维数据的聚类分析问题，可采用主分量分析（PCA）^[13,15]或线性判别分析（LDA）^[16]对数据进行预处理，并通过轮廓分析法^[13,15-16]实现数据聚类结果的可视化验证。

对K-means聚类而言，簇中心可根据期望最大化（expectation-maximization，E-M）算法^[17]逐步迭代获得，因此，对算法的聚类效果进行实时评价具有重要意义。本文从数据特征的统计量出发，采用不同的聚类性能指标来评价聚类质量，通过估计聚类（簇）数量、确定聚类质心和聚类边界，有效解决经典K-means算法存在的问题。

1 聚类数量估计

设p维数据的样本数量为n，则观测数据矩阵（数据集）为 X ={x_ik |i=1, 2,…, n; k=1, 2, … , p}，其中， X 的每行代表一个样本，即 x_i =[x_i1, x_i2,…, x_ip ]；每列为一个特征，即 x_k =[x_1k, x_2k, …, x_nk ]^T。

在聚类分析中，任意2个样本点之间的相异性（dissimilarity）通常采用欧氏距离^[9]来度量，即

d\left({\mathit{x}}_i,{\mathit{x}}_j\right)={\left({\displaystyle \sum _{k=\mathrm{1}}^p}|x_{ik}-x_{jk}{|}^{\mathrm{1}/\mathrm{2}}\right)}^{\mathrm{2}}，（1）

式中，x_ik 和x_jk 分别表示数据样本 x_i 和 x_j 的第k个特征（属性）值。

K-means聚类算法的目标函数为平方误差和（sum-of-squared-error，SSE）^[7,9]，即

I_{\mathrm{S}\mathrm{S}\mathrm{E}}={\displaystyle \sum _{i=\mathrm{1}}^n}{\displaystyle \sum _{j=\mathrm{1}}^K}\left(w_{ij}||{\mathit{x}}_i-{\mathit{m}}_j|{|}_{\mathrm{2}}^{\mathrm{2}}\right)，（2）

式中： m_j (j=1,2,…,K)为第j个聚类C_j 的质心；K为聚类（簇）数量；w_ij 为权重；||∙||₂为L₂范数，即欧氏距离。

式（2）中，若 x_i 属于聚类C_j，权重w_ij 取1，否则取0。这说明SSE本质上是数据簇内的SSE（within cluster SSE）。

K-means算法通过寻求目标函数I_SSE的最小值以实现对数据的聚类。因此，需确定式（2）中的2个重要参数，即K和 m_j。

任意数据集 X 的聚类数K为正整数^[13,16]。考虑2种极端情况：①若整个数据集的所有样本都归属同一簇，则最小簇数K_min=1；②若每个样本独立构成一簇，则最大簇数K_max=n。聚类数的范围[K_min,K_max]=[1,n]，实际应用中通常满足K_{\mathrm{m}\mathrm{a}\mathrm{x}}\ll n。

在数据分析中，不同特征的量纲与取值范围往往存在显著差异，为此需计算每个特征列的样本统计量，即

m_k=\frac{\mathrm{1}}{n}{\displaystyle \sum _{i=\mathrm{1}}^n}{x}_{ik}k=\mathrm{1,2},\cdots ,p，（3）

{\sigma }_k=\sqrt[]{\frac{\mathrm{1}}{n-\mathrm{1}}{\displaystyle \sum _{i=\mathrm{1}}^n}{\left(x_{ik}-m_k\right)}^{\mathrm{2}}}\hspace{1em}k=\mathrm{1,2},\cdots ,p，（4）

式中：m_k 、σ_k 分别为第k个特征列所有样本的均值、标准差。

基于数据的正态假说^[18]，以3σ准则^[19]（即先验法则）推断数据集的聚类结构。对于具有p个特征的 X，其最大簇数K_max的计算式为

K_{\mathrm{m}\mathrm{a}\mathrm{x}}=\mathrm{c}\mathrm{e}\mathrm{i}\mathrm{l}\left\{\underset{k}{\mathrm{m}\mathrm{a}\mathrm{x}}\left(m_k+\mathrm{3}{\sigma }_k\right)\right\}，（5）

式中，ceil{∙}为向上取整函数。

通过式（5）得到数据集的簇数范围 [1,K_max]。同理采用K_{\mathrm{m}\mathrm{i}\mathrm{n}}=\mathrm{c}\mathrm{e}\mathrm{i}\mathrm{l}\left\{\underset{k}{\mathrm{m}\mathrm{a}\mathrm{x}}\left(m_k-\mathrm{3}{\sigma }_k\right)\right\}求出数据集 X 的最小簇数，该公式计算的K_min可能为负值，考虑到簇数最小值为1，本文统一设定K_min=1。

由式（2）可知，I_SSE与簇数K具有相关性。在样本容量n给定的情况下，当簇数K增大时，样本被分配至各簇中的间距随之减小，即I_SSE随着K的增大而减少。在识别聚类数的问题中，肘部法（elbow method）^[15-16]是一种主流的直观判据，其核心及逻辑是在K从大到小的变化过程中，找到I_SSE开始快速增长所对应的K（即SSE曲线的肘部），从而估计簇数。然而，肘部法的应用具有较大的局限性，特别是当簇数较大时，对SSE曲线肘部的判断具有明显的不确定性。因此，为了客观准确评估真实聚类数，需引入多维度的统计性能指标体系，综合考量聚类质量。

1.1 聚类质量评价指标

对不同聚类算法所产生的聚类结果，分别从算法的运行时间和聚类结果的有效性进行综合评价。常用的性能指标如下。

（1）兰德指数（Rand index，RI）兰德指数通过样本对分类一致性进行度量^[20]，定义为

I_{\mathrm{R}\mathrm{I}}=\frac{a+b}{C_n^{\mathrm{2}}}，（6）

式中：C_n^{\mathrm{2}}为总样本对数；a、b分别为真实标签与聚类结果中同为一簇和不同为一簇的正确样本对数。

I_RI的取值范围为[0,1]，I_RI越接近1，表示聚类精度越高。但该指标未剔除随机聚类的影响，仅适用于简单的分类场景。

（2）调整兰德指数（adjusted Rand index，ARI）

当数据簇的数量大于2，且与样本总数处于相同数量级时，随机划分策略对应的I_RI接近0，难以区分有效聚类与随机划分。为了抵消这种随机影响，引入调整兰德指数^[21]，其表达式为

I_{\mathrm{A}\mathrm{R}\mathrm{I}}=\frac{I_{\mathrm{R}\mathrm{I}}-E\left(I_{\mathrm{R}\mathrm{I}}\right)}{\mathrm{m}\mathrm{a}\mathrm{x}\left(I_{\mathrm{R}\mathrm{I}}\right)-E\left(I_{\mathrm{R}\mathrm{I}}\right)}，（7）

式中，E\left(I_{\mathrm{R}\mathrm{I}}\right)为随机聚类下兰德指数的数学期望。

I_ARI考虑了机会概率，其取值范围为[-1,1]。I_ARI越接近1，表示聚类结果越准确；I_ARI越接近0，表示聚类结果与随机选择相当；I_ARI越接近-1，表示聚类结果与真实数据簇完全相反。一般而言，I_ARI取值范围限定于[0,1]。

基于互信息的分数（mutual information-based score，MI）对于样本数为n的数据集X，假设存在两种簇划分方案U和V，采用信息熵量化簇划分的不确定性，两类方案的信息熵定义分别为

H\left(U\right)=-P\left(i\right)\mathrm{l}\mathrm{o}\mathrm{g}\left[P\left(i\right)\right]，（8）

H\left(V\right)=-P\mathrm{\text{'}}\left(j\right)\mathrm{l}\mathrm{o}\mathrm{g}\left[P\left(j\right)\right]，（9）

式中：P\left(i\right)=\left|U_i\right|/n为从U中随机选取的样本落入U_i 类的概率；P\mathrm{\text{'}}\left(j\right)=\left|V_j\right|/n为从V中随机选取的样本落入V_j 类的概率。

U和V之间的互信息分数MI^[22]定义为

I_{\mathrm{M}\mathrm{I}}\left(U,V\right)={\displaystyle \sum _{i=\mathrm{1}}^{|U|}}{\displaystyle \sum _{j=\mathrm{1}}^{|V|}}P\left(i,j\right)\mathrm{l}\mathrm{o}\mathrm{g}\left[\frac{P\left(i,j\right)}{P\left(i\right)P\mathrm{\text{'}}\left(j\right)}\right]，（10）

式中，P\left(i,j\right)=\left|U_i\bigcap V_j\right|/n为随机选取的样本同时属于U_i 和V_j 的概率。

I_MI取值范围为[0,1]。I_MI接近0，表示聚类效果较差；I_MI接近1，表示聚类效果较好。

（4）归一化互信息（normalized mutual information，NMI）

在MI的基础上，进一步定义归一化互信息NMI^[22]，即

I_{\mathrm{N}\mathrm{M}\mathrm{I}}\left(U,V\right)=\frac{I_{\mathrm{M}\mathrm{I}}\left(U,V\right)}{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}\left[H\left(U\right),H\left(V\right)\right]}，（11）

式中，mean[∙]是熵的广义均值。

与I_MI相比，I_NMI不受数据真实类别数量的影响，更为稳健。I_NMI为0，表示U和V没有互信息；I_NMI为1，表示U和V完全相关。

（5）调整的互信息（adjusted mutual information，AMI）

MI的期望为

E\left[I_{\mathrm{M}\mathrm{I}}\left(U,V\right)\right]={\displaystyle \sum _{i=\mathrm{1}}^{|U|}}{\displaystyle \sum _{j=\mathrm{1}}^{|V|}}{\displaystyle \sum _{n_{ij}={\left(a_i+b_j-n\right)}^{+}}^{\mathrm{m}\mathrm{i}\mathrm{n}\left(a_i,b_j\right)}}\frac{n_{ij}}{n}log\left(\frac{n{n}_{ij}}{a_i{b}_j}\right)\cdot

\frac{a_i!b_j!\left(n-a_i\right)!\left(n-b_j\right)!}{n!n_{ij}!\left(a_i-n_{ij}\right)!\left(b_j-n_{ij}\right)!\left(n-a_i-b_j+n_{ij}\right)!}，（12）

式中：a_i为U_i 中的元素数量，a_i=\left|U_i\right|；b_j为V_j 中的元素数量，b_j=\left|V_j\right|；n_ij =(a_i +b_j -n)⁺为同时落入U_i 和V_j 中的元素数量，上标 “+”表示取正值。

通过式（12），得到调整的互信息AMI^[22]为

I_{\mathrm{A}\mathrm{M}\mathrm{I}}=\frac{I_{\mathrm{M}\mathrm{I}}-E\left(I_{\mathrm{M}\mathrm{I}}\right)}{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}\left[H\left(U\right),H\left(V\right)\right]-E\left(I_{\mathrm{M}\mathrm{I}}\right)}\mathrm{。}（13）

I_NMI和I_AMI的取值依赖于熵的广义均值，不同均值形式对应不同的计算结果。本文采用算术平均值作为广义均值，此时，I_NMI和I_AMI的取值范围均为[0,1]。

同质性（homogeneity index，HI）、完备性（completeness index，CI）、V测度（V-measure，VM）

在给定数据真实聚类划分的情况下，基于条件熵可定义一些直观的评价指标。参考文献[23]，定义2个理想的评价指标^[23]，即同质性HI和完备性CI。

HI用于衡量每个簇中是否只包含单一类别的成员（数据样本），其定义为

I_{\mathrm{H}\mathrm{I}}=\mathrm{1}-\frac{H\left(C|K\right)}{H\left(C\right)}，（14）

式中：H\left(C|K\right)为在给定簇K时类别C的条件熵；H\left(C\right)为类别C的信息熵。

CI用于衡量给定类别的所有成员是否均被分配至同一簇中，其定义为

I_{\mathrm{C}\mathrm{I}}=\mathrm{1}-\frac{H\left(K|C\right)}{H\left(K\right)}，（15）

式中：H\left(K|C\right)为给定类别C时聚类K的条件熵；H\left(K\right)为聚类K的信息熵。

I_HI和I_CI从不同角度量化聚类质量，取值范围均为[0,1]，值越大表示聚类效果越好。在实际应用中，主要观察I_HI与I_CI的差值，差值越小聚类效果越好。为了综合考虑HI和CI，采用VM ^[23]来衡量聚类算法对数据的划分质量，即

I_{\mathrm{V}\mathrm{M}}=\frac{\left(\mathrm{1}+\beta \right)I_{\mathrm{H}\mathrm{I}}{I}_{\mathrm{C}\mathrm{I}}}{\beta I_{\mathrm{H}\mathrm{I}}+I_{\mathrm{C}\mathrm{I}}}，（16）

式中：\beta为权重参数，用于协调指标I_HI与I_CI的关系：\beta<1，则赋予I_HI更多的权重；\beta>1，则赋予I_CI更多的权重。本文取\beta=1，则

I_{\mathrm{V}\mathrm{M}}=\mathrm{2}\frac{I_{\mathrm{H}\mathrm{I}}{I}_{\mathrm{C}\mathrm{I}}}{I_{\mathrm{H}\mathrm{I}}+I_{\mathrm{C}\mathrm{I}}}\mathrm{。}（17）

实际上，I_VM等同于I_MI。特别地，当\beta=1时，I_VM等同于I_NMI。

（7）福克斯-马洛斯指数（Fowlkes-Mallows index，FMI）

FMI最初用于衡量两组聚类结果的相似性，是一种无监督的比较。当已知数据样本的真实聚类时，可将其作为监督自适应指标，本质为聚类精度和召回率的几何平均值^[24]，即

I_{\mathrm{F}\mathrm{M}\mathrm{I}}=\frac{\mathrm{T}\mathrm{P}}{\sqrt[]{\left(\mathrm{T}\mathrm{P}+\mathrm{F}\mathrm{P}\right)\left(\mathrm{T}\mathrm{P}+\mathrm{F}\mathrm{N}\right)}}，（18）

式中：TP表示真阳性，即数据集中真实聚类标签和预测标签属于相同簇的数据点对数；FP表示假阳性，即真实聚类标签同簇而预测标签不同簇的数据点对数；FN表示假阴性，即真实标签不同簇、预测标签同簇的数据点对数。

显然，I_FMI基于真实标签和聚类结果的交集、联合集及簇内/簇间数据点对数的比值计算，其取值范围为[0,1]，值越大，表示聚类效果越好。

（8）轮廓系数（silhouette coefficient，SC）

若数据集固有的聚类标签未知，需基于数据自身的相似性开展聚类性能评估，轮廓系数SC是该场景下的经典无监督评价指标。SC依托同一簇中数据样本点间的欧氏距离，以及与下一个邻近簇中所有数据点间的平均欧氏距离，其取值范围为[0,1]。第i个数据点 x_i 的轮廓系数^[25]为

I_{\mathrm{S}{\mathrm{C}}_{\mathrm{i}}}=\frac{e_i-f_i}{\mathrm{m}\mathrm{a}\mathrm{x}\left\{f_i,e_i\right\}}，（19）

式中：e_i 为数据样本 x_i 与同簇中其他样本点之间的平均欧氏距离；f_i 为样本 x_i 与最邻近簇所有数据点的平均欧氏距离。

实际应用中，单个簇中的轮廓系数常采用该簇中所有样本轮廓系数的平均值。本文采用整个数据集全部样本轮廓系数的平均值。平均轮廓系数越高，算法的聚类性能越好。

1.2 确定聚类数

K-means算法的应用前提是已知数据的聚类数K。聚类有效性（cluster validity）分析^[14]正是建立在准确辨识数据中簇数量K的基础之上。前述各聚类性能指标为数据集潜在簇结构的识别提供了客观标准。同时，上述指标的取值范围均为[0,1]，为数据集潜在聚类（簇）数量的估计提供了可行思路：在聚类数的可能范围[1,K_max]内，分别计算不同K对应的聚类性能指标，通过求解各指标的全局最优值，确定数据集的固有簇数量K_best。

2 聚类中心与聚类边界

K-means算法的输出结果为指定数量K_best的聚类划分结果以及聚类中心位置。算法所产生的数据簇满足以下条件^[7]：一是数据簇的中心为该簇内所有数据样本的平均值位置；二是在同一簇中，每个样本点到其簇中心的距离，均小于到其他簇中心的距离。

显然，随着数据集中样本数量的增加，将样本分配至对应数据簇的组合方案数量将急剧增长。若数据中存在噪声和离群值，穷举式搜索将消耗大量时间，导致算法收敛效率显著下降。

对于大样本数据，K-means算法不采用遍历所有可能结果的搜索方法，而是采用期望最大化（E-M）的迭代算法优化搜索策略。E-M算法是统计学中功能强大的迭代算法^[17]，K-means聚类可视为E-M算法的最简单应用形式。在获得簇数量K_best的前提下，K-means算法的E-M迭代步骤如下。

期望步（E-step）：将数据点划分到距离其最近的聚类中心所对应的簇中。

最大化步（M-step）：将聚类中心更新为当前簇内所有数据样本的平均值。

E-M算法为优化迭代过程，其目标是不断更新K_bset个聚类中心的位置。一般地，每一轮E-step和M-step的迭代，均会获得相较上一轮更优的聚类效果。当算法收敛时，最终得到j个聚类（簇）C_j （j=1, 2, …, K_best）的质心为

m_j=\frac{\mathrm{1}}{n_j}{\displaystyle \sum _{i=\mathrm{1}}^{n_j}}{x}_{ik}\hspace{1em}j=\mathrm{1,2},\cdots ,K_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}，（20）

式中，n_j 为第j个聚类C_j 中数据样本的数量。

经E-M迭代收敛后得到的聚类质心m_j，其位置可能落在某个数据点上，也可能位于多个数据点之间。

在同一簇C_j 中，聚类质心到该簇内最远样本点（可能是噪声）之间的欧氏距离为

d_j\underset{}{=\mathrm{m}\mathrm{a}\mathrm{x}}\left\{d\left(x_{ik},m_j\right)\right\}i\in C_j,j=\mathrm{1,2},\cdots ,K_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\mathrm{。}（21）

d_j 即为簇C_j 的半径。以m_j 为圆心、d_j 为半径，即可得到每个簇C_j 的圆形边界。

3 初始聚类质心改进

在经典K-means算法中，随机选择K_best个数据点作为初始聚类质心，可能引发以下问题：一是初始质心位置过于接近，导致在算法后续迭代中2个及以上的聚类被合并，形成空聚类，即聚类中无数据点；二是虽未产生空聚类，但随机初始化质心选择导致聚类效果较差。为了解决随机初始化问题，将最小相似度的聚类中心选择原则与最大-最小（max-min）相似性规则相结合，针对经典K-means算法提出改进方案。

对于数据集 X，首先计算2个数据样本（假设为 x₁和 x₂）的最小相似度，即

S_{im}\left({\mathit{x}}_{\mathrm{1}}^{\mathrm{0}},{\mathit{x}}_{\mathrm{2}}^{\mathrm{0}}\right)=\mathrm{m}\mathrm{i}\mathrm{n}\left\{S_{im}\left({\mathit{x}}_i,{\mathit{x}}_j\right),i,j=\mathrm{1,2},\cdots ,n\right\}，（22）

式中，上标0表示初始分配。

将 x₁、 x₂作为簇C₁、C₂的初始质心，即

{\mathit{x}}_{\mathrm{1}}^{\mathrm{0}}\in C_{\mathrm{1}}^{\mathrm{0}},{\mathit{x}}_{\mathrm{2}}^{\mathrm{0}}\in C_{\mathrm{2}}^{\mathrm{0}}\mathrm{。}（23）

在此基础上，利用最大-最小相似性规则找到第3个数据样本，作为第3个簇的初始质心，即

\mathrm{m}\mathrm{a}\mathrm{x}\left\{S_{im}\left({\mathit{x}}_{\mathrm{3}}^{\mathrm{0}},{\mathit{x}}_r^{\mathrm{0}}\right),r=\mathrm{1,2}\right\}=

\mathrm{m}\mathrm{i}\mathrm{n}\left\{\left[S_{im}\left({\mathit{x}}_s^{\mathrm{0}},{\mathit{x}}_r^{\mathrm{0}}\right),r=\mathrm{1,2}\right],s\ne \mathrm{1,2}\right\}，（24）

式中：r为初始质心的序号；s为第3个质心的序号；{\mathit{x}}_r^{\mathrm{0}}为初始质心。

以此类推，则可依次获得初始K_best个质心，满足

{\mathit{x}}_{\mathrm{1}}^{\mathrm{0}}\in C_{\mathrm{1}}^{\mathrm{0}},{\mathit{x}}_{\mathrm{2}}^{\mathrm{0}}\in C_{\mathrm{2}}^{\mathrm{0}},\cdots ,{\mathit{x}}_{k_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}}^{\mathrm{0}}\in C_{k_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}}^{\mathrm{0}}\mathrm{。}（25）

通过上述改进，可保证初始聚类质心具有较大的相异性，从而在K-means算法后续迭代过程中，聚类不会被合并而形成空聚类。在聚类分析中，数据样本间的相异性可通过欧式距离进行度量，也可以利用夹角余弦等来计算^[7]。

4 本文方法的实施步骤

针对给定数据集 X （样本容量n、特征维度p），本文方法可估计出 X 的簇数量K_best和聚类中心位置。具体步骤如下。

步骤1 计算数据集 X 各特征维度的样本均值m_k 和标准差σ_k (k=1,2,…,p)。

步骤2 基于式（5），利用m_k 和σ_k 估计出数据集 X 的最大簇数K_max。

步骤3 在[1,K_max]范围内，计算不同聚类数K对应的各聚类性能指标。

步骤4 定位聚类性能指标在区间[1,K_max]内的最大值, 从而确定 X 的最佳聚类（簇）数K_best。

步骤5 对于确定的K_best，利用最小相似度规则以及最大-最小相似性规则选择初始聚类中心。

步骤6 采用E-M算法迭代更新聚类中心位置。

步骤7 基于式（21）计算每个簇的半径，从而获得其圆形边界。

5 仿真实验

由K-means算法的目标函数和迭代过程可知，算法适用于球形（或超球形）的数据簇划分。对于非球形簇的数据集，K-means的聚类结果与数据本身的簇结构相差较大。因此，在仿真时，采用不同样本数量、不同簇数量的球形随机斑点数据。为了聚类结果的可视化，仅采用二维数据进行仿真，将经典K-means算法、文献[13]中改进的K-means算法及本文方法的聚类效果进行对比分析。

仿真平台配置：Intel Core i7-8700 CPU、16 GB内存PC，Windows10系统，Spyder 5运行环境。

5.1 小样本数据的聚类

利用随机数生成器，生成由400个样本点构成的二维Blobs数据集，数据散点图见图1，其中，特征1和特征2为数据维度，无量纲。由图1可见，400个样本点大致形成4个数据簇。根据式（5），估计出该数据集的最大簇数K_max为12。

在K∈[1,12]范围内，计算不同K对应的各聚类性能指标，提取各指标在该范围内的最大值以估计数据集的固有簇数K_best。由于K_best大于2，因此在仿真中采用ARI。另外，由于HI和CI仅从单一维度反映聚类的性能，因此采用VM来评价聚类质量而不采用HI和CI。

小样本数据集聚类性能指标随K的变化规律见图2。各性能指标在K∈[1,12]范围内均存在最大值，且最大值均发生在K为4处，表明该数据集的最佳簇数K_best为4。同时，VM、NMI曲线重合，验证了VM等效于NMI且非常接近于AMI的结论。故后续仿真中，不再绘制NMI曲线。

在K_best为4的条件下，经典K-means算法的聚类结果见图3。与原始数据的散点图相比，经典K-means算法因随机选择数据点作为初始质心，除右下角蓝色数据簇外，其余3个数据簇均出现明显不合理划分：中部原本属于2个簇的数据（绿色）被算法合并为1个簇；左上角数据簇被算法拆分为2个不同的簇（粉色簇和橙色簇），且橙色簇中混入上部绿色簇的个别样本点。从聚类性能指标来看，I_AMI为0.791 41，I_VM为0.793 26，I_FMI为0.793 50，I_ARI为0.641 10，I_SC仅为0.592 98，与理想值1相差较大，运行时间为0.015 6 s，聚类性能较差。

在K_best为4的条件下，分别采用改进K-means算法^[13]及本文方法对该数据集进行聚类分析，聚类结果分别见图4、图5。

由图4可以看出，文献[13]算法的聚类结果与原始数据集固有划分结构匹配度较好，其I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.994 67、0.991 56、0.991 59、0.996 00、0.751 78，运行时间 0.015 7 s。由图5可以看出，本文方法的聚类结果与原始数据簇的划分完全一致。I_ARI、I_AMI、I_VM、I_FMI均达到理想值1，I_SC为0.885 96，接近理想值1，且算法耗时小于0.000 1 s，聚类质量与迭代效率均实现显著提升。

采用本文方法确定数据集的聚类中心后，基于式（21）计算各簇的聚类半径，以聚类中心为圆心，绘制各簇的分界线，见图6。

由图6可以看出，本文方法绘制的圆形聚类边界线精确度较高，说明本文方法能够精准估计和正确识别数据集的簇数量，同时E-M算法为各个聚类中心提供准确定位。

5.2 大样本数据集的验证

为了验证本文算法在“样本容量较大和簇数较多”场景下的聚类性能，开展不同样本规模和不同簇数量的仿真实验。

（1） 样本数为2 100、簇数为7的数据集

随机生成2 100个样本、7个簇的二维数据集，其散点图见图7。

分别计算该数据集2个特征维度的样本均值和标准差，根据式（5）得出该数据集的最大可能簇数K_max为14。在K∈[1,14]范围内，计算不同簇数K对应的各聚类性能指标，并确定各指标的全局最大值，结果见图8。

由图8可知，各聚类性能指标均在K为7处取得最大值，即该数据集的最佳簇数K_best为7。

在K_best为7的条件下，经典K-means算法的聚类结果见图9。该算法将中下方蓝色的2个簇合并为1个簇，将中部上方的1个簇拆分为2个簇（浅蓝色簇和灰色簇），导致数据集固有聚类结构发生扭曲失真。I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.805 36、0.916 66、0.917 04、0.839 76、0.769 11，耗时0.015 7 s，聚类效果较差。

在相同实验条件下，分别采用文献[13]改进K-means算法和本文方法对该数据集进行聚类分析，结果分别见图10、图11。

由图10、图11可以看出，两种算法的聚类结构与原始数据集固有簇结构一致。文献[13]改进算法的I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.949 36、0.969 05、0.969 33、0.958 26、0.813 55，耗时0.015 6 s。本文方法的I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.995 56、0.993 75、0.993 78、0.996 19、0.895 46，相对更接近理想值1，且耗时小于0.000 1 s，收敛性能明显提高，聚类质量显著提升。

根据上述计算结果，绘制本文方法对该数据集的聚类分界线，见图12。尽管样本量增加导致相邻簇数据点存在交叉，本文方法仍能准确识别不同聚类的边界线。

（2） 样本数为7 000、簇数为10的数据集

随机生成7 000个样本、10个簇的二维数据集，其散点图见图13。

根据式（5），估计该数据集的最大簇数K_max为19。在K∈[1,19]范围内，各类聚类性能指标随K的变化规律见图14，所有指标在K为10时取得最大值，因此确定该数据集的最佳簇数K_best为10。

在K_best为10的条件下，经典K-means算法的聚类结果见图15。算法将右上部2个簇合并为1个簇（橙色簇），将左上部1个簇分裂为2个簇（蓝色簇和玫粉色簇）。I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.867 76、0.950 65、0.950 78、0.884 00、0.822 06，耗时0.046 9 s。

分别采用文献[13]改进算法与本文方法对该数据集进行聚类分析，结果分别见图16、图17。两种算法的数据划分与数据集固有的聚类结构完全一致。其中，文献[13]改进算法的I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.936 27、0.954 32、0.956 76、0.947 63、0.834 66，耗时0.044 4 s；本文方法的I_ARI、I_AMI、I_VM、I_FMI、I_SC分别为0.999 05、0.998 68、0.998 68、0.999 14、0.834 66，耗时0.015 6 s。本文方法的I_SC超过0.83，I_ARI、I_FMI、I_AMI、I_VM均高度接近理想值1，且运行时间仅为经典K-means算法的30%，收敛性能与聚类质量显著提升。

对比图11、图17的聚类性能指标可知，随着数据集样本数量和聚类簇数量的增加，本文方法的聚类性能不仅未下降，反而显著提升，表明该方法对于数据集的样本容量和聚类（簇）数均具有较强的鲁棒性。

采用本文方法绘制10个数据簇的边界线，结果如图18所示。

由图18可知，尽管样本规模较大导致顶部2个簇的数据点存在空间交叉，本文方法仍能较准确地给出各数据簇的边界。

6 结论

针对经典K-means算法的固有缺陷，提出一种改进K-means方法。采用3σ准则估计数据的簇数区间，以区间内多维度聚类性能指标的全局最大值确定数据集潜在的簇数；通过最小相似度原理、最大-最小相似度规则解决经典K-means算法的随机初始化问题。以E-M算法迭代收敛的簇质心为圆心、簇内距离质心最远样本点的欧氏距离为半径，构建各簇的圆形边界曲线，有效抑制噪声与离群值对聚类结果的干扰。仿真结果表明，本文方法可为球形数据集提供高精度的聚类划分，为监督分类器提供未知数据集的类标签。

本研究仅针对球形结构数据开展验证，存在一定局限。后续将引入数据变换技术，将本文方法的适用范围拓展至更复杂结构数据集的聚类任务。

参考文献

[1] BANDARU S, NG A H C, DEB K.Data mining methods for knowledge discovery in multi-objective optimization:part A-Survey[J].Expert Systems with Applications,2017,70:139-159.

[2] BUI A T.Dimension reduction with prior information for knowledge discovery[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2024,46(5):3625-3636.

[3] CZARNOWSKI I, JĘDRZEJOWICZ P.Supervised classification problems-taxonomy of dimensions and notation for problems identification[J].IEEE Access,2021,9:151386-151400.

[4] SINGH J, SINGH D.A comprehensive review of clustering techniques in artificial intelligence for knowledge discovery:taxonomy,challenges, applications and future prospects[J].Advanced Engineering Informatics, 2024,62:102799.

[5] MAHNOOR, SHAFI I, CHAUDHRY M,et al.A review of approaches for rapid data clustering:challenges, opportunities,and future directions [J].IEEE Access,2024,12:138086-138120.

[6] SUBASI A.Machine learning techniques[M]//Practical Machine Learning for Data Analysis Using Python.Amsterdam:Elsevier,2020: 91-202.

[7] XU R, WUNSCH-Ⅱ D C.Clustering[M].New Jersey:IEEE Press & WILEY,2009:67-73.

[8] STEINLEY D. K-means clustering:a half-century synthesis[J].British Journal of Mathematical and Statistical Psychology,2006,59:1-34.

[9] STEINLEY D.Stability analysis in K-means clustering[J].British Journal of Mathematical and Statistical Psychology,2008,61:255-273.

[10] SELIM S Z, ISMAIL M A. K-means-type algorithms:a generalized convergence theorem and characterization of local optimality[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,1984, PAMI-6(1): 81-87.

[11] HE X S, HE F, FAN Y P,et al.An effective clustering scheme for high-dimensional data[J].Multimedia Tools and Applications,2024, 83:45001-45045.

[12] IKOTUN A M, EZUGWU A E, ABUALIGAH L,et al.K-means clustering algorithms:a comprehensive review, variants analysis, and advances in the era of big data[J].Information Sciences,2023,622: 178-210.

[13] 何选森,何帆,于海澜.K-均值算法的初始化改进与聚类质量评估[J]. 西安工程大学学报,2024,38(6):114-123.

[14] HE Xuansen, HE Fan, YU Hailan.Initialization improvement and clustering quality evaluation of K-means algorithm[J].Journal of Xi’an Polytechnic University,2024,38(6):114-123.

[15] 何选森,何帆,徐丽,等.K-Means 算法最优聚类数量的确定[J].电子科技大学学报,2022,51(6):904-912.

[16] HE Xuansen, HE Fan, XU Li,et al.Determination of the optimal number of clusters in K-means algorithm[J].Journal of University of Electronic Science and Technology of China,2022,51(6):904-912.

[17] 何选森,何帆,樊跃平,等.高维数据聚类数量可视化确定模式[J].沈阳航空航天大学学报,2024,41(3):71-84.

[18] HE Xuansen, HE Fan, FAN Yueping,et al.Visualized determination mode for clustering quantity of high-dimensional data[J].Journal of Shenyang Aerospace University,2024,41(3):71-84.

[19] 何帆,何选森,刘润宗,等.数据降维与K-均值聚类的质量评估[J].重庆理工大学学报(自然科学),2024,38(1):131-141.

[20] HE Fan, HE Xuansen, LIU Runzong,et al.Data dimensionality reduction and clustering quality evaluation of K-means clustering[J].Journal of Chongqing University of Technology (Natural Science),2024,38(1):131-141.

[21] IWASAKI Y, SASAKI Y, NAGATA T,et al.Dynamic mode decomposition based on expectation-maximization algorithm for simultaneous system identification and denoising[J].Mechanical Systems and Signal Processing,2025,223:111864.

[22] XIE J R, HONG T, LAING T,et al.On normality assumption in residual simulation for probabilistic load forecasting[J].IEEE Transactions on Smart Grid,2017,8(3):1046-1053.

[23] CREIGHTON J H C.A first course in probability models and statistical inference[M].New York,NY:Springer New York,1994:140-156.

[24] HUBERT L, ARABIE P.Comparing partitions[J].Journal of Classification,1985,2:193-218.

[25] STEINLEY D.Properties of the hubert-arabie adjusted rand index[J]. Psychological Methods,2004,9(3):386-396.

[26] STREHL A, GHOSH J.Cluster ensembles-a knowledge reuse framework for combining multiple partitions[J].Journal of Machine Learning Research,2002,3:583-617.

[27] ROSENBERG A, HIRSCHBERG J.V-measure:a conditional entropy-based external cluster evaluation measure[C]//Proceedings of the 2007 Joint Conference on Empirical Methods in Natural Language Processing and Computational Natural Language Learning (EMNLP-CoNLL).Prague,Czech RepublicAssociation for Computational Linguistics,2007:410-420.

[28] FOWLKES E B, MALLOWS C L.A method for comparing two hierarchical clusterings[J].Journal of the American Statistical Association,1983,78(383):553-569.

[29] BAGIROV A M, ALIGULIYEV R M, SULTANOVA N.Finding compact and well-separated clusters:clustering using silhouette coefficients[J].Pattern Recognition,2023,135:109144.

基金资助

广东省普通高校特色创新项目（自然科学）(2024KTSCX131)