Matlab笔记-模糊聚类分析原理及实现

23. 模糊聚类分析原理及实现聚类分析，就是用数学方法研究和处理所给定对象，按照事物间的相似性进行区分和分类的过程。

传统的聚类分析是一种硬划分，它把每个待识别的对象严格地划分到某个类中，具有非此即彼的性质，这种分类的类别界限是分明的。

随着模糊理论的建立，人们开始用模糊的方法来处理聚类问题，称为模糊聚类分析。

由于模糊聚类得到了样本数与各个类别的不确定性程度，表达了样本类属的中介性，即建立起了样本对于类别的不确定性的描述，能更客观地反映现实世界。

本篇先介绍传统的两种（适合数据量较小情形，及理解模糊聚类原理）：基于择近原则、模糊等价关系的模糊聚类方法。

（一）预备知识一、模糊等价矩阵定义1设R=(r ij )n ×n 为模糊矩阵，I 为n 阶单位矩阵，若R 满足 i) 自反性：I ≤R （等价于r ii =1）； ii) 对称性：R T =R;则称R 为模糊相似矩阵，若再满足iii) 传递性：R 2≤R （等价于1()nik kj ij k r r r =∨∧≤）则称R 为模糊等价矩阵。

定理1设R 为n 阶模糊相似矩阵，则存在一个最小的自然数k（k <n ）, 使得R k 为模糊等价矩阵，且对一切大于k 的自然数l ，恒有R l =R k . R k 称为R 的传递闭包矩阵，记为t(R). 二、模糊矩阵的λ-截矩阵定义2设A =(a ij )n ×m 为模糊矩阵，对任意的λ∈[0,1], 作矩阵()()ij n mA a λλ⨯=其中，()1, 0, ij ijij a aa λλλ≥⎧=⎨<⎩称为模糊矩阵A 的λ-截矩阵。

显然，A λ为布尔矩阵，且其等价性与与A 一致。

意义：将模糊等价矩阵转化为等价的布尔矩阵，可以得到有限论域上的普通等价关系，而等价关系是可以分类的。

因此，当λ在[0,1]上变动时，由A λ得到不同的分类。

若λ1＜λ2, 则A λ1≥A λ2, 从而由A λ2确定的分类是由A λ1确定的分类的加细。

当λ从1递减变化到0时，A λ的分类由细变粗，逐渐归并，形成一个分级聚类树。

例1设U={u 1, u 2, u 3, u 4, u 5}, 对给定的U 上的模糊等价关系让λ从1到0变化，观察分类过程。

(1) 当λ=1时，110000 01000 00100 00010 00001R⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦分类结果为5类：（每行代表一类，1代表对应元素在该类）{u1}, {u2}, {u3}, {u4}, {u5}(2) 当λ=0.8时，0.810100 01000 10100 00010 00001R⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦分类结果为4类：{u1, u3}, {u2}, {u4}, {u5}(3) 当λ=0.6时，0.610100 01000 10100 00011 00011R⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦分类结果为3类：{u1, u3}, {u2}, {u4, u5}(4) 当λ=0.5时，0.510111 01000 10111 10111 10111R⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦分类结果为2类：{u1, u3, u4, u5}, {u2}(4) 当λ=0.4（R 中的最小值）时，0.41111111111111111111111111R ⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦分类结果为1类：{u 1, u 2, u 3, u 4, u 5}整个动态分类过程如下：（二）基于择近原则的模糊聚类择近原则就是利用贴近度来实现分类操作，贴近度用来衡量两个模糊集A 和B 的接近程度，用N (A ,B )表示。

贴近度越大，表明二者越接近。

设论域有限或者在一定区间，即U={u 1, u 2, …, u n }或U=[a,b], 常用的贴近度有以下三种： (1) 海明贴近度11(,)1|()()|ni i i N A B A u B u n ==--∑1(,)1|()()|d bi i a N A B A u B u u b a=---⎰ (2) 欧氏贴近度1221(,)1[()()]ni iiN A B A u B u=⎫=--⎪⎭∑)122(,)1[()()]dbi iaN A B A u B u u=--⎰(3) 格贴近度(,)()()c cN A B A B A B=∧其中，()1()()ni iiA B A u B u==∨∧.Matlab实现：格贴近度的实现函数fuz_closing.mfunction y=fuz_closing(A,B,type)%要求A与B列数相同的行向量[m,n]=size(A);switch typecase 1 %海明贴近度y=1-sum(abs(A-B))/n;case 2 %欧氏贴近度y=1-(sum(A-B).^2)^(1/2)/sqrt(n);case 3 %格贴近度y1=max(min(ones(m,n)-A,ones(m,n)-B));%ones(m,n)-A等于A^cy2=max(min(A,B));y=min(y1,y2);end例2设某产品的质量等级分为5级，其中一级有5种评判因素u1, u2, u3, u4, u5. 每一等级的模糊集为B1={0.5 0.5 0.6 0.4 0.3}B2={0.3 0.3 0.4 0.2 0.2}B3={0.2 0.2 0.3 0.1 0.1}B4={0.1 0.1 0.2 0.1 0}B5={0.1 0.1 0.1 0.1 0}假设某产品各评判因素的值为A={0.4 0.3 0.2 0.1 0.2}, 问该产品属于哪个等级？代码：A=[0.4 0.3 0.2 0.1 0.2];B=[0.5 0.5 0.6 0.4 0.3;0.3 0.3 0.4 0.2 0.2;0.2 0.2 0.3 0.1 0.1;0.1 0.1 0.2 0.1 0;0.1 0.1 0.1 0.1 0];for i=1:5haiming(i)=fuz_closing(A,B(i,:),1);oushi(i)=fuz_closing(A,B(i,:),2);ge(i)=fuz_closing(A,B(i,:),3);endhaimingoushige运行结果：haiming = 0.7800 0.9200 0.9000 0.8600 0.8400 oushi = 0.5081 0.9106 0.8658 0.6870 0.6422ge = 0.4000 0.3000 0.2000 0.2000 0.1000可见样本A与各等级的格贴近度分别为0.4, 0.3, 0.2, 0.2, 0.1, 故可认为该产品属于B1等级。

若按令两种贴近度判断，该产品属于B2等级。

（三）基于模糊等价关系的模糊聚类一、算法步骤1. 样本数据归一化设X={x 1, x 2, …, x n }为要分类的n 个样本，每个样本有m 个指标，即x i ={ x i 1, x i 2, …, x im }, i =1,2,..,n得到原始数据矩阵X=(x ij )n ×m .由于不同指标的数据量纲不同，为了使数据能够比较，要先对X 做归一化处理。

2. 建立模糊相似矩阵R先建立样本x i 与x j 相似程度r ij , 进而构造模糊相似矩阵R=(r ij )n ×n建立r ij 常用的方法有：(1) 相似系数法①夹角余弦法：mikjkij xx r ⋅=∑②相关系数法：||||miki jk j ij xx x x r -⋅-=∑(2)距离法一般取r ij =1-c (d (x i ,x j ))α, 其中c 和α为适当选取的参数，使得 0≤r ij ≤1. 常用的距离有：①海明距离：1(,)||mi j ik jk k d x x x x ==-∑②欧氏距离：(,)i j d x x =③切比雪夫距离：1(,)max ||i j ik jk k md x x x x ≤≤=- (3) 贴近度法①最大最小法：11()()mikjk k ij mikjk k x x r xx ==∧=∨∑∑②算术平均最小法：11()1()2mikjk k ij m ik jk k xx r x x ==∧=+∑∑③几何平均最小法：11()mikjk k ij mk xx r ==∧=∑3. 求出R 的传递闭包t(R)即改造相似关系为等价关系：令2R R R =, 再令422R R R =, …, 直到满足2l l l R R R =与R l 相等，即为t(R), 仍记为R.4. 选取合适的λ, 利用λ-截矩阵R λ进行分类（参考例1）。

二、Matlab 实现求模糊相似矩阵R 的函数：fuz_distance.mfunction R=fuz_distance(x,type)%x 为归一化的数据矩阵, type 选择计算相似程度的方法 %返回模糊相似矩阵R[n,m]=size(x);%距离法的选择参数c和a, 需要根据具体情况修改以保证R(i,j)属于[0,1]c=0.1;a=1;for i=1:nfor j=1:nswitch typecase 1 %夹角余弦法R(i,j)=(x(i,:)*x(j,:)')/(norm(x(i,:),2)*norm(x(j,:),2));case 2 %相关系数法Dxi=abs(x(i,:)-mean(x(i,:)));Dxj=abs(x(j,:)-mean(x(j,:)));R(i,j)=(Dxi*Dxj')/(norm(Dxi,2)*norm(Dxj,2));case 3 %海明距离法d=sum(abs(x(i,:)-x(j,:)));R(i,j)=1-c*d^a;case 4 %欧氏距离法d=norm(x(i,:)-x(j,:),2);R(i,j)=1-c*d^a;case 5 %切比雪夫距离法d=max(abs(x(i,:)-x(j,:)));R(i,j)=1-c*d^a;case 6 最大最小(贴近度)法R(i,j)=sum(min([x(i,:);x(j,:)]))/sum(max([x(i,:);x(j,:)]));case 7 算术平均最小(贴近度)法R(i,j)=2*sum(min([x(i,:);x(j,:)]))/sum(x(i,:)+x(j,:));case 8 %几何平均最小(贴近度)法R(i,j)=sum(min([x(i,:);x(j,:)]))/sum(sqrt(x(i,:).*x(j,:)));endendend求R的传递闭包t(R)的函数：tran_R.mfunction [B,k]=tran_R(R)%R为模糊相似矩阵, 循环构造满足传递性的t(R)%k为满足R^2k = R^k的最小的自然数kn=length(R);B=zeros(n,n);flag=0;k=1/2;while flag==0B=fco(R,R); %做模糊合成运算k=2*k;if B==Rflag=1;elseR=B; %循环计算R传递闭包endend上面的函数tran_R.m调用函数矩阵模糊合成算子函数：fco.m function B=fco(Q,R)%实现模糊合成算子的计算, 要求Q的列数等于R的行数[n,m]=size(Q);[m,l]=size(R);B=zeros(n,l);for i=1:nfor k=1:lB(i,k)=max(min([Q(i,:);R(:,k)'])); endend求t(R)的λ-截矩阵的函数：fuz_lamda.mfunction y=fuz_lamda(X,m)%用λ-截矩阵将样本分成m类, m≤总样本数lamda=unique(X)'; %根据R中的值取λ值%unique函数取矩阵不重复元素组成向量并从小到大排好序X(find(X<lamda(m)))=0;X(find(X>=lamda(m)))=1;y=X;例3某地区设有11个雨量站，其分布如图所示：10年来各雨量站测得的年降雨量表如下：现因经费问题，希望撤销几个雨量站，问撤销哪些雨量站而不会太多地减少降雨信息？分析：对11个雨量站进行模糊聚类，同一类的只需保留一个即可。