Fuzzy C Means 算法及其 Python 实现

1. K-Means 算法向 FCM 算法的扩展

K-Means 算法中,如果要将数据集合 X=\{ X_1,X_2,X_3,\dots,X_n \} 划分为 k\,(1 \le k \le n) 个类,使得任意数据对象 X_i 必须属于并且仅属于一个类,同时每一个类至少包含一个数据对象,那么可以用一个 k\times n 的矩阵 U 来表示,矩阵中的任意一个元素 u_{ij} 可以表示为:

    \[{u_{ij}} = \left{ \begin{cases} 1& {X_i \in G_j} \\ 0& {X_i \notin G_j} \end{cases} \right.\]

其中 {G_j}\left( {j = 1,2, \ldots ,k} \right) 表示第 j 个类。并且 U 需要满足如下条件 (1 \le i \le k,\,1 \le j \le n)

    \[\left{ \begin{cases} {u_{ij}} \in \{0,1\} \\ \sum\limits_{i = 1}^k {{u_{ij}}} = 1 \\ \sum\limits_{j = 1}^n {{u_{ij}}} > 0 \end{cases} \right.\]

如果上述矩阵 U 中的元素 u_{ij} 的取值范围不仅仅是 0 或者 1,那么就可以推广到模糊集合上的划分,U 就变成了模糊判定矩阵。此时 {u_{ij}} 需满足:

(1)   \begin{equation*}  \left{ \begin{cases} {u_{ij}} \in [ 0,1 ]} \\ \sum\limits_{i = 1}^k {{u_{ij}}} = 1 \\ \sum\limits_{j = 1}^n {{u_{ij}}} > 0 \end{cases} \right. \end{equation*}

2. 目标函数与聚类中心

K-Means 算法在度量数据对象的非相似性(或者说距离)时一般使用欧几里得距离,要求每个类的聚类中心与数据对象的距离平方之和最小,目标函数可以表示为:

    \[J = \sum\limits_{i = 1}^k {\sum\limits_{j = 1}^n {s_{ij}^2} }\]

    \[{s_{ij}} = Eculid({C_i},{X_j})\]

其中 C_i 表示任意聚类中心,而聚类中心一般取类内所有对象在各属性上的平均值,因此可以表示为:

    \[{C_i} = \frac{{\sum\limits_{j,{X_j} \in {G_i}} {{X_j}} }}{{\sum\limits_{j = 1}^n {{u_{ij}}} }}\]

{G_i}{\kern 1pt} \left( {1 \le i \le k} \right) 表示任意一个类。

将算法推广到模糊集后,Dunn 对样本与类中心之间的距离采用隶属度的平方来加权,Bezdek 则进一步引入了隶属度的加权指数 m 从而得到了新的目标函数:

(2)   \begin{equation*}  J = \sum\limits_{i = 1}^k {\sum\limits_{j = 1}^n {{{\left( {{u_{ij}}} \right)}^m}s_{ij}^2} } \end{equation*}

要使得 (2) 式达到最小值则要求聚类中心 C_i 和隶属度 u_{ij} 满足如下条件:

(3)   \begin{equation*}  {C_i} = \frac{{\sum\limits_{j = 1}^n {u_{ij}^m{X_j}} }}{{\sum\limits_{j = 1}^n {u_{ij}^m} }} \end{equation*}

(4)   \begin{equation*}  {u_{ij}} = \frac{1}{{\sum\limits_{l = 1}^k {{{\left( {\frac{{{u_{ij}}}}{{{u_{lj}}}}} \right)}^{2/\left( {m - 1} \right)}}} }} \end{equation*}

3. FCM 算法计算过程

由于聚类中心和隶属度两个必要条件的存在,FCM 算法的可以表示为一个简单的迭代过程。将 n 个拥有 m 个属性的对象构建成一个 n\times m 的矩阵记为 DFCM 的算法可以描述为:

输入:隶属度加权指数 m,数据矩阵 D,聚类的个数 k,最大迭代次数 iterMax,阈值 threshold
输出:隶属度矩阵 U
算法过程:
Step1: 用区间 [0,1] 内的满足 (1) 式的随机数初始化隶属度矩阵,记为 D^{(0)}。初始化迭代计数器 iterCounter
Step2: 根据 (3) 式和 U^{(iterCounter)} 计算 k 个聚类中心,记做 C^{(iterCounter)}
Step3: 根据 (2) 式计算目标函数的值。如果本次迭代的目标函数值与上一次计算的目标函数值之差小于 threshold 或者迭代计数器 iterCounter 大于等于 iterMax 那么迭代结束,输出结果,算法停止。
Step4: 依据 C^{(iterCounter)} 和 (4) 式计算 U^{(iterCounter+1)}iterCounter++ ,返回Step2。

4. Python 实现代码

import numpy as np

def loadDataFromTxt(fileName):
    dataSet=np.mat(np.loadtxt(fileName,delimiter='\t'))
    return dataSet

def normalization(dataSet):
    colNum=np.shape(dataSet)[1]
    for index in range(colNum):
        colMax=np.max(dataSet[:,index])
        dataSet[:,index]=dataSet[:,index]/colMax
    return dataSet

def initWithFuzzyMat(n,k):
    fuzzyMat=np.mat(np.zeros((k,n)))
    for colIndex in range(n):
        memDegreeSum=0
        randoms=np.random.rand(k-1,1)
        for rowIndex in range(k-1):
            fuzzyMat[rowIndex,colIndex]=randoms[rowIndex,0]*(1-memDegreeSum)
            memDegreeSum+=fuzzyMat[rowIndex,colIndex]
        fuzzyMat[-1,colIndex]=1-memDegreeSum
    return fuzzyMat

def eculidDistance(vectA,vectB):
    return np.sqrt(np.sum(np.power(vectA-vectB,2)))

def calCentWithFuzzyMat(dataSet,fuzzyMat,p):
    n,m=dataSet.shape
    k=fuzzyMat.shape[0]
    centroids=np.mat(np.zeros((k,m)))
    for rowIndex in range(k):
        degExpArray=np.power(fuzzyMat[rowIndex,:],p)
        denominator=np.sum(degExpArray)
        numerator=np.array(np.zeros((1,m)))
        for colIndex in range(n):
            numerator+=dataSet[colIndex]*degExpArray[0,colIndex]
        centroids[rowIndex,:]=numerator/denominator
    return centroids

def calFuzzyMatWithCent(dataSet,centroids,p):
    n,m=dataSet.shape
    c=centroids.shape[0]
    fuzzyMat=np.mat(np.zeros((c,n)))
    for rowIndex in range(c):
        for colIndex in range(n):
            d_ij=eculidDistance(centroids[rowIndex,:],dataSet[colIndex,:])
            fuzzyMat[rowIndex,colIndex]=1/np.sum([np.power(d_ij/eculidDistance(centroid,dataSet[colIndex,:]),2/(p-1)) for centroid in centroids])
    return fuzzyMat

def calTargetFunc(dataSet,fuzzyMat,centroids,k,p):
    n,m=dataSet.shape
    c=fuzzyMat.shape[0]
    targetFunc=0
    for rowIndex in range(c):
        for colIndex in range(n):
            targetFunc+=eculidDistance(centroids[rowIndex,:],dataSet[colIndex,:])**2*np.power(fuzzyMat[rowIndex,colIndex],p)
    return targetFunc

def fuzzyCMean(dataSet,k,p,initMethod=initWithFuzzyMat):
    n,m=dataSet.shape
    fuzzyMat=initWithFuzzyMat(n,k)
    centroids=calCentWithFuzzyMat(dataSet,fuzzyMat,p)
    lastTargetFunc=calTargetFunc(dataSet,fuzzyMat,centroids,k,p)

    fuzzyMat=calFuzzyMatWithCent(dataSet,centroids,p)
    centroids=calCentWithFuzzyMat(dataSet,fuzzyMat,p)
    targetFunc=calTargetFunc(dataSet,fuzzyMat,centroids,k,p)
    while lastTargetFunc*0.99>targetFunc:
        lastTargetFunc=targetFunc
        fuzzyMat=calFuzzyMatWithCent(dataSet,centroids,p)
        centroids=calCentWithFuzzyMat(dataSet,fuzzyMat,p)
        targetFunc=calTargetFunc(dataSet,fuzzyMat,centroids,k,p)
    return fuzzyMat,centroids


if __name__=='__main__':
    dataSet=loadDataFromTxt('data.txt')
    dataSet=normalization(dataSet)
    fuzzyMat,centroids=fuzzyCMean(dataSet,3,2)
    print 'fuzzyMat=\n',fuzzyMat
    print np.sum(fuzzyMat,axis=0)
    print 'centroids=\n',centroids