聚类之详解FCM算法原理及应用

【之前】 该文的pdf清晰版已被整理上传，方便保存学习，下载地址：

https://download.csdn.net/download/on2way/

FCM原理介绍

FCM分析1
FCM分析2
FCM分析3

首先介绍一下模糊这个概念，所谓模糊就是不确定，确定性的东西是什么那就是什么，而不确定性的东西就说很像什么。比如说把20岁作为年轻不年轻的标准，那么一个人21岁按照确定性的划分就属于不年轻，而我们印象中的观念是21岁也很年轻，这个时候可以模糊一下，认为21岁有0.9分像年轻，有0.1分像不年轻，这里0.9与0.1不是概率，而是一种相似的程度，把这种一个样本属于结果的这种相似的程度称为样本的隶属度，一般用u表示，表示一个样本相似于不同结果的一个程度指标。

基于此，假定数据集为X，如果把这些数据划分成c类的话，那么对应的就有c个类中心为C，每个样本j属于某一类i的隶属度为$u_{ij}$，那么定义一个FCM目标函数（1）及其约束条件（2）如下所示：

$$J=\sum _{i=1}^c\sum _{j=1}^n{u}_{ij}^m||x_j-c_i|{|}^2\text{\hspace{1em}(1)}$$ $$\sum _{i=1}^c{u}_{ij}=1,j=1,2...,n\text{\hspace{1em}(2)}$$
看一下目标函数（式1）而知，由相应样本的隶属度与该样本到各个类中心的距离相乘组成的，m是一个隶属度的因子，个人理解为属于样本的轻缓程度，就像$x^2$与$x^3$这种一样。式（2）为约束条件，也就是一个样本属于所有类的隶属度之和要为1。观察式（1）可以发现，其中的变量有$u_{ij}、c_i$，并且还有约束条件，那么如何求这个目标函数的极值呢？

分析式(3)，先对第一部分的两级求和的$u_{ij}$求导，对求和形式下如果直接求导不熟悉，可以把求和展开如下：
\begin{bmatrix} u_{11}^m||x_1-c_1||2 &\cdots &u_{1j}^m||x_j-c_1||2 &\cdots & u_{1n}^m||x_n-c_1||2 \\vdots &\vdots& \vdots& \vdots&\vdots\u_{i1}^m||x_1-c_i||2 &\cdots &u_{ij}^m||x_j-c_i||2 &\cdots & u_{in}^m||x_n-c_i||2 \\vdots& \vdots& \vdots&\vdots&\vdots\u_{c1}^m||x_1-c_c||2 &\cdots &u_{cj}^m||x_j-c_c||2 &\cdots & u_{cn}^m||x_n-c_c||2\end{bmatrix}
这个矩阵要对$u_{ij}$求导，可以看到只有$u_{ij}$对应的$u_{ij}^m||x_j-c_i|{|}^2$保留，其他的所有项中因为不含有$u_{ij}$，所以求导都为0。那么$u_{ij}^m||x_j-c_i|{|}^2$对$u_{ij}$求导后就为$m||x_j-c_i|{|}^2{u}_{ij}^{m-1}$。

再来看后面那个对$u_{ij}$求导，同样把求和展开，再去除和$u_{ij}$不相关的（求导为0），那么只剩下这一项：${\lambda }_j(u_{ij}-1)$，它对$u_{ij}$求导就是${\lambda }_j$了。

那么最终J对$u_{ij}$的求导结果并让其等于0就是：
$$\frac{\partial J}{\partial u_{ij}}=m||x_j-c_i|{|}^2{u}_{ij}^{m-1}+{\lambda }_j=0$$
这个式子化简下，将$u_{ij}$解出来就是：
$$u_{ij}^{m-1}=\frac{-{\lambda }_j}{m||x_j-c_i|{|}^2}$$
进一步：
$ u_{ij}=(\frac{-\lambda_j}{m||x_j-c_i||^2}){\frac{1}{m-1}}=(\frac{-\lambda_j}{m})^{{\frac{1}{m-1}}(\frac{1}{||x_j-c_i||}{(\frac{2}{m-1})}})\tag{4}$

要解出$u_{ij}$则需要把${\lambda }_j$去掉才行。这里重新使用公式（2）的约束条件，并把算出来的$u_{ij}$代入式（2）中有：

$1=\sum _{i=1}^c{u}_{ij}=\sum _{i=1}^c(\frac{-{\lambda }_j}{m}{)}^{\frac{1}{m-1}}(\frac{1}{||x_j-c_i|{|}^{(\frac{2}{m-1})}})=(\frac{-{\lambda }_j}{m}{)}^{\frac{1}{m-1}}\sum _{i=1}^c(\frac{1}{||x_j-c_i|{|}^{(\frac{2}{m-1})}})$

好了，式子（5）就是最终的$ u_{ij}$迭代公式。

下面在来求J对$c_i$的导数。由公式（2）可以看到只有$\sum _{i=1}^c\sum _{j=1}^n{u}_{ij}^m||x_j-c_i|{|}^2$这一部分里面含有$c_i$，对其二级求和展开如前面所示的，那么它对$c_i$的导数就是：
$\frac{\partial J}{\partial c_i}=\sum _{j=1}^n(-u_{ij}^m\ast 2\ast (x_j-c_i))=0$
即：
$\sum _{j=1}^n(u_{ij}^m{c}_i)=\sum _{j=1}^n(x_j{u}_{ij}^m)$ KaTeX parse error: \tag works only in display equations

好了，公式（6）就是类中心的迭代公式。

我们发现$ u_{ij}$与$c_i$是相互关联的，彼此包含对方，有一个问题就是在fcm算法开始的时候既没有$ u_{ij}$也没有$c_i$，那要怎么求解呢？很简单，程序开始的时候我们随便赋值给$ u_{ij}$或者$c_i$其中的一个，只要数值满足条件即可。然后就开始迭代，比如一般的都赋值给$ u_{ij}$，那么有了$ u_{ij}$就可以计算$c_i$，然后有了$c_i$又可以计算$ u_{ij}$，反反复复，在这个过程中还有一个目标函数J一直在变化，逐渐趋向稳定值。那么当J不在变化的时候就认为算法收敛到一个比较好的结了。可以看到$ u_{ij}$和$c_i$在目标函数J下似乎构成了一个正反馈一样，这一点很像EM算法，先E在M，有了M在E，在M直至达到最优。

还需要说一点的是，当程序结束后，怎么去判断哪个点属于哪个类呢？在结束后，肯定有最后一次计算的U吧，对于每一个点，它属于各个类都会有一个u，那么找到其中的最大的u就认为这个点就属于这一类。基于此一个基础的程序如下：

clc clear close all %% create samples: for i=1:100 x1(i) = rand()*5; %人为保证差异性 y1(i) = rand()*5; x2(i) = rand()*5 + 3; %人为保证差异性 y2(i) = rand()*5 + 3; end x = [x1,x2]; y = [y1,y2]; data = [x;y]; data = data';%一般数据每一行代表一个样本 %plot(data(:,1),data(:,2),'*'); %画出来 %%--- cluster_n = 2;%类别数 iter = 50;%迭代次数 m = 2;%指数 num_data = size(data,1);%样本个数 num_d = size(data,2);%样本维度 %--初始化隶属度u，条件是每一列和为1 U = rand(cluster_n,num_data); col_sum = sum(U); U = U./col_sum(ones(cluster_n,1),:); %% 循环--规定迭代次数作为结束条件 for i = 1:iter %更新c for j = 1:cluster_n u_ij_m = U(j,:).^m; sum_u_ij = sum(u_ij_m); sum_1d = u_ij_m./sum_u_ij; c(j,:) = u_ij_m*data./sum_u_ij; end %-计算目标函数J temp1 = zeros(cluster_n,num_data); for j = 1:cluster_n for k = 1:num_data temp1(j,k) = U(j,k)^m*(norm(data(k,:)-c(j,:)))^2; end end J(i) = sum(sum(temp1)); %更新U for j = 1:cluster_n for k = 1:num_data sum1 = 0; for j1 = 1:cluster_n temp = (norm(data(k,:)-c(j,:))/norm(data(k,:)-c(j1,:))).^(2/(m-1)); sum1 = sum1 + temp; end U(j,k) = 1./sum1; end end end figure; subplot(1,2,1),plot(J); [~,label] = max(U); %找到所属的类 subplot(1,2,2); gscatter(data(:,1),data(:,2),label) 

这里以matlab下的灰度图像为例。灰度图像一图像的角度看是二维的，但是我们知道，决定图像的无非是里面的灰度值。而灰度值就是一个值，所以当我们把图像变成1维，也就是拉成一行或者一列的时候，其实灰度图像就是一个一维数据（上面那个例子生成的随机点是二维的）。

clc clear close all img = double(imread('lena.jpg')); subplot(1,2,1),imshow(img,[]); data = img(:); %分成4类 [center,U,obj_fcn] = fcm(data,4); [~,label] = max(U); %找到所属的类 %变化到图像的大小 img_new = reshape(label,size(img)); subplot(1,2,2),imshow(img_new,[]); 

需要注意的是label出来的是标签类别（1-4），并非真实的灰度，这里不过是把它显示出来就行了。

这个数据库看介绍好像是关于种子分类的，里面共包含3类种子，每类种子通过什么x射线技术等等采集他们的特征，反正最后每个种子共有7个特征值来表示它（也就是说在数据里面相当于7维），每类种子又有70个样本，那么整个数据集就是210*7的样本集。从上面那个地方下载完样本集存为txt文件，并放到matlab工作目录下就可以使用了（注意看看下下来的数据有没有串位的，有的话要手动调整回去）。因为matlab只能显示低于3维的数据，这里有7维，我们现在在二维下显示其中的两维以及正确的分类结果看看什么情况：

clc clear close all data = importdata('data.txt'); %data中还有第8列，正确的标签列 subplot(2,2,1); gscatter(data(:,1),data(:,6),data(:,8)),title('choose:1,6 列') subplot(2,2,2); gscatter(data(:,2),data(:,4),data(:,8)),title('choose:2,4 列') subplot(2,2,3); gscatter(data(:,3),data(:,5),data(:,8)),title('choose:3,5 列') subplot(2,2,4); gscatter(data(:,4),data(:,7),data(:,8)),title('choose:4,7 列') 

组合有限，随便组合几种看看，发现似乎任意两个特征都可以把他们分开，当然还是有一些分不开的，其中最后一个选择特征4,7似乎很好的分开了。

Ok，看过之后我们来试试fcm算法对其进行分类，并计算一下准确率，我们先把7个特征都用上看看：

clc clear close all data = importdata('data.txt'); %data中还有第8列，正确的标签列 [center,U,obj_fcn] = fcm(data(:,1:7),3); [~,label] = max(U); %找到所属的类 subplot(1,2,1); gscatter(data(:,4),data(:,7),data(:,8)),title('choose:4,7列,理论结果') % cal accuracy a_1 = size(find(label(1:70)==1),2); a_2 = size(find(label(1:70)==2),2); a_3 = size(find(label(1:70)==3),2); a = max([a_1,a_2,a_3]); b_1 = size(find(label(71:140)==1),2); b_2 = size(find(label(71:140)==2),2); b_3 = size(find(label(71:140)==3),2); b = max([b_1,b_2,b_3]); c_1 = size(find(label(141:210)==1),2); c_2 = size(find(label(141:210)==2),2); c_3 = size(find(label(141:210)==3),2); c = max([c_1,c_2,c_3]); accuracy = (a+b+c)/210; % plot answer subplot(1,2,2); gscatter(data(:,4),data(:,7),label),title(['实际结果,accuracy=',num2str(accuracy)]) 

这里选择以第1与6维的数据来可视化这个结果。可以看到准确率为0.89524。

试了很多组特征，都没有超过0.89524的，那就所有特征都用上吧。其实这个准确率是可以提高的，我们看到这7个特征似乎有点重复有没有，如果我们把这7个特征采用pca降维到3,4个特征了再去fcm实验呢？可以去试试，有待实验…