基于滤波器的聚类算法

Vol.43 No.10Oct. 2010

基于滤波器的聚类算法

张 强，李 成，吴腾飞

(天津大学精密测试技术及仪器国家重点实验室，天津 300072)

摘 要：针对聚类分析精度和效率低的问题，设计了一种聚类算法FBCLUS．应用卷积定理和傅里叶变换，提出了频率滤波法来消除噪声的干扰；提出了单阈值、多阈值幅度滤波法消除噪声和提取不同密度的感兴趣区间；提出一个数学形态学算子提取聚类簇．实验表明：FBCLUS算法能够发现任意形状的聚类；速度快，计算复杂度为O(N)；能够发现不同密度的聚类簇；抗噪声性能强；对网格大小有一定的适应性．FBCLUS算法有很高的聚类精度和效率． 关键词：频率滤波；幅度滤波；数学形态学；聚类

中图分类号：TP391 文献标志码：A 文章编号：0493-2137(2010)10-0884-06

A Filter-Based Clustering Algorithm

ZHANG Qiang，LI Cheng，WU Teng-fei

(State Key Laboratory of Precision Measuring Technology and Instruments，Tianjin University，Tianjin 300072，China)

a clustering algorithm FBCLUS was pro-Abstract：To overcome the problems of low accuracy and low efficiency，

posed. Firstly，convolution theorem and Fourier transform were used to design frequency filters to reduce noise；single threshold and multi-threshold amplitude filters were introduced to reduce noise and distinguish different densityregions of interest．Thirdly，a mathematical morphology clustering operator was designed to discover clusters. The experiments show that FBCLUS is able to detect arbitrarily shaped clusters，and it is very efficient with a complexity of O(N)；it can distinguish clusters of different density；it is insensitive to large amounts of noise；and it is not sen-sitive to the grid size. FBCLUS has high accuracy and efficiency.

Keywords：frequency filter；amplitude filter；mathematics morphology；clustering

聚类分析是一种重要的无指导学习方法，有着广泛的应用领域．在图像处理方面，聚类被用来实现图像的自动分割、人脸的自动识别．在地球物理方面，聚类被用来处理卫星遥感信息．在生物医学方面，聚类被用来分析基因表达的关系．在文本处理方面，聚类被用来实现因特网中文本的自动分类．

当前的聚类算法主要有4类：①划分方法，如文献[1-2]等，其缺点是确定参数K困难，算法对噪声敏感；②层次方法，如文献[3-5]等，其缺陷是聚类质量较低；③密度方法，如文献[6-7]等，其缺点是计算量过大，不适合处理大量数据，对参数过分敏感；④网格的方法，如CLUGD[8-9]等，其缺点是对网格大小过于敏感，生成的簇的形状不是水平的就是垂直的，结果精

收稿日期：2009-06-02；修回日期：2009-10-20.

度低．

针对以上问题，笔者提出了基于滤波器的聚类算法 (filter-based clustering algorithm，FBCLUS)．其主要贡献是：①利用卷积定理和傅里叶变换，提出了频率滤波法来消除噪声的干扰；②提出了单阈值、多阈值幅度滤波法消除噪声和提取不同密度的感兴趣区域(可能包含聚类簇的区域)；③提出了一个新的数学形态学算子，从感兴趣区域中提取聚类簇．

FBCLUS算法的特点有：①能够发现任意形状的聚类；②计算复杂度与数据量呈线性关系；③能够区分不同密度性质的聚类簇；④抗噪声能力强；⑤对网格大小有一定的适应性；⑥聚类结果与数据输入顺序无关．

基金项目：中国博士后科学基金资助项目(20090450767)；天津市高等学校科技发展基金资助项目((20080810)．作者简介：张 强（1972— ），男，博士，讲师． 通讯作者：张 强，zq8622@gmail.com.

2010年10月 张 强等：基于滤波器的聚类算法 ·885·

1 FBCLUS聚类算法

FBCLUS算法在思想上借鉴了数字信号处理的原理：提出聚类滤波方法降低噪声水平，通过数学形态学方法提取聚类簇．FBCLUS分为3步：①将数据空间网格化；②进行频率滤波和幅度滤波；③以数学形态学方法发现聚类簇．

FBCLUS的第1步是将数据空间划分为等间距的网格，每个网格保存着其包含数据对象的个数．如果将数据集看作多维度信号，那么网格中数据对象的个数就标志着该网格处信号的强弱，通过信号的滤波和数学形态学处理，最终得到聚类簇．数据空间网格化的目的是浓缩信息，提高后续处理的速度．确定合理的网格大小(或者疏密)对后续处理有重要影响，目前提出的网格聚类算法对网格大小十分敏感，网格大小的变化对聚类结果有严重的影响．笔者提出的FBCLUS算法由于采用了模板滤波和数学形态学聚类，所以对网格的大小有一定的适应性． 1.1 频率滤波

实际应用中获得的数据通常含有大量的噪声对象，它们不但影响聚类的精度，同时也加大了计算量，所以有必要将它们滤除．如果将数据集视为多维度信号，那么聚类簇和噪声所在区域表现为不同的特 征[10]：聚类簇所在的区域信号以低频为主，信号的幅度大，信号的能量高；噪声所在的区域信号以高频为主，信号的幅度小，能量低．

文中提出了使用频率和幅度两种滤波器过滤数据集的方法，以提高聚类处理的精度．首先，使用频率低通滤波器去除噪声所对应的高频信号；然后，使用幅度滤波器滤除残存的小幅度噪声信号．经过两次滤波处理后，噪声水平明显降低，为提高聚类簇的质量提供了保证．为了降低频率滤波的计算复杂度，提出了时域和频域的综合滤波法，避免了对整个空间的变换和反变换．实验表明，该方法有很好的精度和效率．

频率滤波可以直接施加于信号，也可以先对信号做变换(如傅里叶变换、小波变换等)，然后对变换后的结果进行滤波．前者称为时域滤波，因为信号所处的空间是时域空间；后者称为频域滤波，因为变换后的空间是频域空间．在时域中进行频率滤波的优点是计算量低，因为它不需要进行空间变换．但是在时域中信号的频率特性没有被明显展示，所以频率滤波器设计困难．在频域中进行频率滤波的优点是滤波

器设计简单直观，因为信号直接表现为频率的函

数．但是频域滤波需要对数据集进行变换，计算量大．

文献[11]在聚类中采用了频域滤波，使用小波变换实现时域到频域的变换，在频域中使用低通滤波器实现过滤．文献[11]只考虑了不同区域的频率特性而没有考虑幅度特性．其算法优点是抗噪声强、精度高；缺点是需要对整个数据空间做小波变换和反变换，计算量大．

笔者提出的频率滤波方法综合了时域和频域滤波的优点，避免了各自的缺点，不但精度高而且减小了计算量．其核心思想是：在频域中设计低通滤波器；根据卷积定理，使用傅里叶反变换得到时域空间中对应的低通滤波器；在时域中完成滤波．其优点是：①在频域中能够方便地设计出高质量的滤波器，通过傅里叶反变换在时域中得到的对应滤波器具有频域中相同的性能，这是由卷积定理保证的；②使用时域滤波器直接滤波避免了空间的傅里叶变换和傅里叶反变换，减少了计算量．

下面详细讨论频域中低通滤波器的设计和时域中对应滤波器的获得以及时域中滤波的处理．

频域中频率滤波的数学模型为 G(u,v)=H(u,v)F(u,v)

式中：F(u,v)为数据集f(x,y)在频域的傅里叶变换，也就是要被滤波的对象；H(u,v)是频域中滤波器函数，完成对信号的滤波；

G(u,v)是对F(u,v)滤波后的结果．问题的关键是如何确定频率滤波器H(u,v)．以原型法设计频率滤波器，确定滤波器的类型，再确定滤波器的参数．采用原型法的最大好处是现有的滤波器理论成熟，可以很容易地设计出符合要求的滤波器．主要考虑理想低通滤波器、巴特沃思低通滤波器和高斯低通滤波器[12]．以这3种滤波器为设计原型是因为它们的频率特性曲线涵盖了从尖锐陡峭(理想低通滤波器)到圆滑平缓(高斯低通滤波器)的整个范围．理想低通和高斯低通滤波器是两个极端，巴特沃思低通滤波器则介于二者之间，为过渡类型．

在频域中完成滤波器的设计之后需要将它变换到时域中，进而完成滤波．卷积定理提供了理论依据，即频域的乘积对应于时域的卷积 F(u,v)H(u,v)⇔f(x,y)∗h(x,y)

由卷积定理和傅里叶反变换就可以获得需要的

时域滤波器h(x,y)，

从而可以直接对时域中的数据滤波．现在时域空间是大小为R×T的网格，

应采用离散形式的卷积运算，即

·886· 天 津 大 学 学 报 第43卷 第10期

f(x,y)∗h(x,y)=

1

R−1T−1

RT

∑∑

f(i,j)h(x−i,y−j)

i=0j=0

式中：

h(x,y)的大小(作用范围)是R×T个网格．为了提高处理速度需要将h(x,y)做得小一些，文献[12] 指出以小尺寸的h(x,y)在时域滤波的效率高于频域滤波，同时有很好的滤波效果．下面采用大小为r×t的

滤波器，并用模板实现与网格的卷积，即

f(x,y)∗h(x,y)=1r−1trt∑∑−1

f(i,j)h(x−i,y−j)=

i=0j=0

w󰀁+wrt1f1+w2f2+rtfrt=∑wi=1

ifi 式中：

wi为掩模的系数；fi是与该系数对应的网格信号强度(即网格中含有数据对象的个数)．

以模板计算卷积的过程是在网格中逐个网格地移动模板，在每一个网格(x，y)处滤波器的响应由模板的掩模系数与相应网格乘积之和给出．对于一个

3×3的模板(见图1)，在网格(x，y)处卷积结果为

f(x,y)∗h(x,y)=w−1,−1f(x−1,y−1)+…+

w0,0f(x,y)+…+w1,1f(x+1,y+1)

计算中要求网格(x，y)与模板中心w0,0对应．对于r×t的模板，

r和t都要求是奇数． w−1,−1 w−1,0w−1,1 w0,−1 w0,0w0,1 w1,−1

w1,0 w1,1 图1 3×3的模板 Fig.1 3×3 mask

经过频率滤波之后高频噪声被大幅度衰减，但还会有些残存，笔者采用幅度滤波法以进一步降低噪声 水平． 1.2 幅度滤波

表1为DBSCAN和FBCLUS的运行时间．从表1可知，噪声区域的信号强度明显小于聚类簇区域，

表1

DBSCAN和FBCLUS的运行时间 Tab.1

Run time of DBSCAN and FBCLUS 对象个数 运行时间/s DBSCAN FBCLUS 100 000 69.2 0.87

200 000 195.3 1.71

300 000 576.8

2.53 500 000

1 277.9

4.37

所以文中提出通过信号幅度的差异来进一步滤除噪

声．幅度滤波的思想是确定一个幅度阈值ξ，所有幅度小于ξ的区域都视为噪声加以滤除，其数学表达为

⎧⎨

cell(x,y)∈ROI⎩cell(x,y)∉ROI ifif Ifrefiltered(x,y)≥ξIfrefiltered

(x,y)＜ξ 式中cell(x,y)为中心坐标为(x，y)的网格，Ifrefiltered(x,y)为网格cell(x,y)经过频率滤波后的信号

强度；ROI(regions of interest)表示后续聚类处理的感兴趣区域．如果网格(x，y)的信号强度高于幅度阈值ξ，

则表明该网格可能属于某个聚类簇，需要算法进一步分析处理；反之，表明该网格信号强度太小，则属于噪声，要加以滤除．

幅度滤波的关键是如何确定合理的幅度阈值ξ，笔者以直方图方式选择幅度阈值．通常情况下，经过频率滤波后的信号幅度是非负的实数，呈现为连续分布．为了构建直方图，需要对其离散化，采用均匀量化，将信号幅度转化为非负整数．然后分别计算各个幅度值对应网格数量的百分数，其数学处理方法是：设非空网格的总数为M，幅度值为i的网格数量是

m分比为Pm

i，则幅度值为i的网格数量百i=i度值I为横轴，以M

．以幅

Pi为纵轴生成幅度直方图．幅度直方图一般呈多峰、谷分布，因为噪声信号的强度很低，

所以选取第一个明显的谷底对应的幅度值作为阈值

ξ即可．

通过以上的幅度滤波，整个空间被分为ROI区域(可能属于某个聚类簇)和噪声区域，所以这种幅

度滤波也可以看作是一种二值划分问题．

为了进一步提高聚类分析的精确性，又提出了多阈值的幅度滤波方法，它对应于多值划分问题，其核心思想是：呈多峰、谷分布的幅度直方图中每一个峰及其附近隆起部分都反映某一种特定的信号强度分布，同时也对应于特定的数据点密度分布．而两个峰之间的谷底则是两种不同密度分布的分割点，选取这些谷底对应的幅度值作为幅度阈值，就可以将数据空间划分为不同密度的区域，每一种密度区域都可能对应不同的物理性质．多幅度阈值的数学表达为

⎧⎪cell(x,y)∈ΡΟΙ0ifIfrefiltered(x,y)＜ξ0

⎪

cell(x,y)∈ΡΟΙ

1ifξ0≤Ifrefiltered(x,y)＜ξ1 ⎪⎪

⎨

󰀁⎪cell(x,y)∈ROIif ξIξ ⎪

i

i−1≤frefiltered(x,y)＜i⎪󰀁

⎪⎩cell(x,y)∈ROIzifξz−1≤Ifrefiltered(x,y)＜ξz

式中0≤ξ0＜ξ1＜󰀁＜ξi＜󰀁＜ξz，是一组选定的阈值．区域ROI1到ROIz是不同密度的感兴趣区域，对应于不同的ROI，需要分别处理来发现不同密度的聚

类簇．区域ROI0是噪声区域，要加以滤除．通过多阈值幅度滤波，FBCLUS算法能够发现不同密度的聚

2010年10月 张 强等：基于滤波器的聚类算法 ·887·

类簇．

1.3 数学形态学聚类

经过频率和幅度滤波之后，得到ROI区域，对于二值幅度滤波，ROI区域只有一个，而对于多值幅度滤波，ROI区域有多个，每一个ROI对应的区域密度不同，需要分别处理．下面任务是从每个ROI区域中提取由相连网格构成的聚类簇．

笔者提出了一个数学形态学算子完成聚类簇的搜索，其伪代码见算法，其中⊕代表膨胀运算，Δ为膨胀运算的结构元素．算子处理过程：第1步，以ROI区域作为种子集合seeds，并从中任意取出一个没有

被处理的网格p作为起始点X0i进行第i轮迭代．第2

步，以起始点为中心进行膨胀，并将膨胀后新得到的

网格作为中心继续膨胀，直到没有发现新的网格为止．为了保证膨胀不会超出ROI区域，需要将新得到的网格与seeds取交．第3步，将上面膨胀得到的结果归属于一个新的聚类簇Ci，并从种子集合中去除已经处理过的部分，如果种子集合不为空，则返回第1步继续第i+1轮迭代．此外，这一步还滤除了体积过小的聚类簇．第4步，当迭代结束后输出得到的聚类簇． 算法 i=1

seeds=ROI

WHILE(seeds≠Φ)

∀p∈seeds，

X0i=p k=1

Xki=(Xki−1⊕Δ)∩seeds

Add =Xki-Xki−1 WHILE(add≠Φ)

k=k+1

Xki=(Xki−1∪(add⊕Δ))∩seeds

add=Xki-Xki−1 END WHILE Ci=Xki

seeds=seeds-(Ci∩seeds) i = i+1

if Ci is too small then ignore it END WHILE OUTPUT {Ci}

END

笔者提出的数学形态学算子本身具有一定的降噪能力，它主要消除那些体积过小相互集聚的网格，如图2所示．图2(a)中小黑点是一些高密度网格的

聚集体，数学形态学算子能够将它们滤除(见图

2(b)．图2(a)中小黑点出现的原因是：为了提高效率，对频率和幅度滤波均采用小尺寸的模板，从而导致抗局部高密度噪声能力的下降．但是，数学形态学算子只能区分体积大小的不同，而不能区分密度的疏密，所以必须与频率和幅度滤波相结合才能取得较好的效果．

(a) 降噪前 (b) 降噪后

图2 新的数学形态学算子的降噪效果

Fig.2 Filtering effect of new mathematical morphological

operator

1.4 算法参数的选择

滤波器模板的大小对算法的计算复杂度有很大的影响，模板越大，计算量越大，文献[12]指出在使用小模板的情况下，卷积的计算量会明显小于傅里叶变换．通常模板大小需要由实验决定，经验表明：在数据点密集区域模板中能够包含10个以上的数据点就能基本满足要求．在某些情况下，小尺寸模板会导致局部高密度噪声的滤除不净(见图2)，这可以由提出的数学形态学算子解决．

膨胀运算结构元素的大小决定了算法对提取聚类簇时的黏合性的强弱．当网格划分过密时，聚类区域相对分散，要采用较大的结构元素；反之，当网格划分稀疏时，聚类簇区域相对集中，要采用较小的结构元素．可以通过试探法决定结构元素的大小：首先结构元素的大小取一个较小的值，然后逐步增大结构元素的大小，此时得到的聚类簇个数会逐步减小，当增加结构元素的大小而聚类簇个数变化不大的时候，就找到了合适的大小．

对于多阈值幅度滤波得到的多个ROI区域，由于不同类聚类簇的密度差异较大，需要分别确定结构元素的大小，然后再提取聚类簇，这样就能够保证类内样本较为稀疏的时候，算法判别正确．

滤波器模板和膨胀运算结构元素的形状在低维度下选择矩形，以保证高的聚类质量；在高维度下选择十字形，以减少计算量．因为矩形的模板或者运算结构中包含网格的数目与维度呈指数增长，而十字形中包含网格的数目与维度呈线性增长，所以十字结构可以降低高维度处理的计算量．

·888· 天 津 大 学 学 报 第43卷 第10期

1.5 算法性能分析

为了讨论方便，规定数据集中对象的个数为N，非空网格的数目为M，滤波模板中网格的数量为W，滤波得到的ROI网格的总数目为R，膨胀结构元素中网格的数量为D，FBCLUS算法第1步空间网格化的计算复杂度为O(N)，第2步空间过滤的计算复杂度为O(MW)，第3步数学形态学的计算复杂度为O(RD)，整体复杂度为O(N+MW+RD)．通常MW2 实验测试

2.1 精度测试

测试算法发现任意形状聚类簇的能力和抗噪声性能．测试数据集和聚类结果如图3所示，其中数据集含有噪声对象．

(a) 数据集 (b) 聚类结果

图3 精度测试的数据集和聚类结果

Fig.3 Data set and clustering result of accuracy test

2.2 网格适应性测试

测试算法对网格尺寸的适应性．测试数据集如图4(a)所示，使用两种网格进行聚类．第1种网格是50×50，第2种网格是1,000×1,000，网格边长相差20倍．目前提出的网格聚类算法一般只能在50×50的网格中正常工作，例如CLUGD[8]算法在提取聚类簇时只考虑相邻网格的联通性，当网格划分过于密集时，高密度网格之间被空网格分割开来，其连通性遭到破坏，从而导致算法不能正常工作．而FBCLUS算法在两种网格下均能得到图4(b)的聚类结果．

(a) 数据集 (b) 聚类结果

图4 网格适应性测试的数据集和聚类结果

Fig.4 Data set and clustering result of grid size test

2.3 含障碍物聚类的测试

以含障碍物聚类为例测试算法功能的可扩展性，测试数据集和聚类结果如图5所示，它包含两个障碍物，聚类算法需要区分出被障碍物分割开的聚类簇．本测试将含有障碍物的网格标记为膨胀禁止区，为了防止禁止区的泄漏采用3×3的十字形膨胀结构元素．聚类结果表明，算法区分出左侧被障碍物完全分割开的聚类簇．

(a) 数据集 (b) 聚类结果

图5 含障碍物的测试数据集和聚类结果

Fig.5 Data set and clustering result of obstacle test

2.4 效率测试

测试算法的效率，数据集由随机数产生器生成：每个数据集含有5个聚类簇和10%的噪声．测试计算机为神舟承运M715D，算法FBCLUS和DBSCAN[13]使用JAVA语言编写，处理时间见表1．测试结果表明FBCLUS的速度明显高于DBSCAN．

3 结 语

笔者设计的基于滤波器的聚类算法FBCLUS提

出以频率滤波器和幅度滤波器来消除噪声并获得感兴趣区域ROI，设计一个新的数学形态学算子从ROI中提取聚类簇．理论分析和测试结果表明，FBCLUS算法能够发现任意形状的聚类；计算复杂度与数据量呈线性关系；能够区分不同密度性质的聚类簇；抗噪声能力强；对网格大小有一定的适应性． 参考文献：

［1］ Har-Peled S，Kushal A. Smaller corsets for k-median and

k-means clustering[J]. Discrete and Computational Ge-ometry，2007，37(1)：3-19.

［2］ Rosa La，Nehorai P S，Eswaran A，et al. Detection of uter-ine MMG contractions using a multiple change point es-timator and the k-means cluster algorithm[J]. IEEE Transactions on Biomedical Engineering，2008，55(2)：453-467.

2010年10月 张 强等：基于滤波器的聚类算法 ·889· ［3］ Loewenstein Y，Portugaly E，Fromer M，et al. Efficient

algorithms for accurate hierarchical clustering of huge datasets：Tackling the entire protein space[J]. Bioinfor-matics，2008，24(13)：i41-i49.

［4］ Lee H. E，Park K H，Bien Z Z. Iterative fuzzy clustering

algorithm with supervision to construct probabilistic fuzzy rule base from numerical data[J]. IEEE Transac-tions on Fuzzy Systems，2008，16(1)：263-277.

［5］ Santos J M，Sa de J M，Alexandre L A，et al. LEGClust：A

clustering algorithm based on layered entropic sub-graphs[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence，2008，30(1)：62-75.

［6］ Duan Lian， Xu Lida，Guo Feng，et al. A local-density

based spatial clustering algorithm with noise[J]. Infor-mation Systems，2007，32：978-986.

［7］ Li Jianhua，

Behjat L，Kennings A. Net cluster：A net-reduction-based clustering preprocessing algorithm for partitioning and placement[J]. IEEE Transactions on Computer-Aided Design of Integrated Circuits and Sys-tems，2007，26(4)：669-679.

［8］ Sun Zhiwei，Zhao Zheng. A fast clustering algorithm

based on grid and density[C]//2005 Canadian Confer-ence on Electrical and Computer Engineering. Saskatoon

Sask，Canada：IEEE Press，2005：2276-2279.

［9］ Rasheed，Tinku Javaid，Usman Meddour，et al. An effi-cient stable clustering algorithm for scalable mobile multi-hop networks[C]// 4th IEEE Consumer Communi-cations and Networking Conference. Las Vegas，USA，2007：89-94

［10］ Han Jiawei，Kamber M. Data Mining Concepts and Tech-niques[M]. 2nd ed. San Francisco，CA：Morgan Kauf-mann Publishers，2006.

［11］ 战立强，刘大昕，张健沛. 基于小波滤波的时间序列子

序列聚类[J]. 计算机工程与应用，

2007，43(10)：4-7. Zhan Liqiang，Liu Daxin，Zhang Jianpei. Time series sub-sequence clustering based on wavelet filters[J]. Com-puter Engineering and Applications，2007，43(10)：4-7(in Chinese).

［12］ Gonzalez R C，Woods R E. Digital Image Process-ing[M]. 3th ed. New Jersey，USA：Pearson Prentice Hall，2008.

［13］ Ester M，Kriegel H P，Sander J，et a1. A density-based

algorithm for discovering clusters in large spatial data-bases [C]//Proceedings of 2nd International Conference on Knowledge Discovery and Data Mining. Poland：ACM Press，1996：226-231.

.