K-means聚类算法

1.原文：http://www.cnblogs.com/jerrylead/archive/2011/04/06/2006910.html

K-means也是聚类算法中最简单的一种了，但是里面包含的思想却是不一般。最早我使用并实现这个算法是在学习韩爷爷那本数据挖掘的书中，那本书比较注重应用。看了Andrew Ng的这个讲义后才有些明白K-means后面包含的EM思想。

聚类属于无监督学习，以往的回归、朴素贝叶斯、SVM等都是有类别标签y的，也就是说样例中已经给出了样例的分类。而聚类的样本中却没有给定y，只有特征x，比如假设宇宙中的星星可以表示成三维空间中的点集。聚类的目的是找到每个样本x潜在的类别y，并将同类别y的样本x放在一起。比如上面的星星，聚类后结果是一个个星团，星团里面的点相互距离比较近，星团间的星星距离就比较远了。

在聚类问题中，给我们的训练样本是，每个，没有了y。

K-means算法是将样本聚类成k个簇（cluster），具体算法描述如下：

1、随机选取k个聚类质心点（cluster centroids）为。

2、重复下面过程直到收敛 {

对于每一个样例i，计算其应该属于的类

对于每一个类j，重新计算该类的质心

}

K是我们事先给定的聚类数，代表样例i与k个类中距离最近的那个类，的值是1到k中的一个。质心代表我们对属于同一个类的样本中心点的猜测，拿星团模型来解释就是要将所有的星星聚成k个星团，首先随机选取k个宇宙中的点（或者k个星星）作为k个星团的质心，然后第一步对于每一个星星计算其到k个质心中每一个的距离，然后选取距离最近的那个星团作为，这样经过第一步每一个星星都有了所属的星团；第二步对于每一个星团，重新计算它的质心（对里面所有的星星坐标求平均）。重复迭代第一步和第二步直到质心不变或者变化很小。

下图展示了对n个样本点进行K-means聚类的效果，这里k取2。

K-means面对的第一个问题是如何保证收敛，前面的算法中强调结束条件就是收敛，可以证明的是K-means完全可以保证收敛性。下面我们定性的描述一下收敛性，我们定义畸变函数（distortion function）如下：

J函数表示每个样本点到其质心的距离平方和。K-means是要将J调整到最小。假设当前J没有达到最小值，那么首先可以固定每个类的质心，调整每个样例的所属的类别来让J函数减少，同样，固定，调整每个类的质心也可以使J减小。这两个过程就是内循环中使J单调递减的过程。当J递减到最小时，和c也同时收敛。（在理论上，可以有多组不同的和c值能够使得J取得最小值，但这种现象实际上很少见）。

由于畸变函数J是非凸函数，意味着我们不能保证取得的最小值是全局最小值，也就是说k-means对质心初始位置的选取比较感冒，但一般情况下k-means达到的局部最优已经满足需求。但如果你怕陷入局部最优，那么可以选取不同的初始值跑多遍k-means，然后取其中最小的J对应的和c输出。

下面累述一下K-means与EM的关系，首先回到初始问题，我们目的是将样本分成k个类，其实说白了就是求每个样例x的隐含类别y，然后利用隐含类别将x归类。由于我们事先不知道类别y，那么我们首先可以对每个样例假定一个y吧，但是怎么知道假定的对不对呢？怎么评价假定的好不好呢？我们使用样本的极大似然估计来度量，这里是就是x和y的联合分布P(x,y)了。如果找到的y能够使P(x,y)最大，那么我们找到的y就是样例x的最佳类别了，x顺手就聚类了。但是我们第一次指定的y不一定会让P(x,y)最大，而且P(x,y)还依赖于其他未知参数，当然在给定y的情况下，我们可以调整其他参数让P(x,y)最大。但是调整完参数后，我们发现有更好的y可以指定，那么我们重新指定y，然后再计算P(x,y)最大时的参数，反复迭代直至没有更好的y可以指定。

这个过程有几个难点，第一怎么假定y？是每个样例硬指派一个y还是不同的y有不同的概率，概率如何度量。第二如何估计P(x,y)，P(x,y)还可能依赖很多其他参数，如何调整里面的参数让P(x,y)最大。这些问题在以后的篇章里回答。

这里只是指出EM的思想，E步就是估计隐含类别y的期望值，M步调整其他参数使得在给定类别y的情况下，极大似然估计P(x,y)能够达到极大值。然后在其他参数确定的情况下，重新估计y，周而复始，直至收敛。

上面的阐述有点费解，对应于K-means来说就是我们一开始不知道每个样例对应隐含变量也就是最佳类别。最开始可以随便指定一个给它，然后为了让P(x,y)最大（这里是要让J最小），我们求出在给定c情况下，J最小时的（前面提到的其他未知参数），然而此时发现，可以有更好的（质心与样例距离最小的类别）指定给样例，那么得到重新调整，上述过程就开始重复了，直到没有更好的指定。这样从K-means里我们可以看出它其实就是EM的体现，E步是确定隐含类别变量，M步更新其他参数来使J最小化。这里的隐含类别变量指定方法比较特殊，属于硬指定，从k个类别中硬选出一个给样例，而不是对每个类别赋予不同的概率。总体思想还是一个迭代优化过程，有目标函数，也有参数变量，只是多了个隐含变量，确定其他参数估计隐含变量，再确定隐含变量估计其他参数，直至目标函数最优。

原文：http://www.cnblogs.com/moondark/archive/2012/03/08/2385770.html

聚类算法——K-means（上）

　　首先要来了解的一个概念就是聚类，简单地说就是把相似的东西分到一组，同 Classification (分类)不同，对于一个 classifier ，通常需要你告诉它“这个东西被分为某某类”这样一些例子，理想情况下，一个 classifier 会从它得到的训练集中进行“学习”，从而具备对未知数据进行分类的能力，这种提供训练数据的过程通常叫做 supervised learning (监督学习)，而在聚类的时候，我们并不关心某一类是什么，我们需要实现的目标只是把相似的东西聚到一起，因此，一个聚类算法通常只需要知道如何计算相似度就可以开始工作了，因此 clustering 通常并不需要使用训练数据进行学习，这在 Machine Learning 中被称作 unsupervised learning (无监督学习)。

　　我们经常接触到的聚类分析，一般都是数值聚类，一种常见的做法是同时提取 N 种特征，将它们放在一起组成一个 N 维向量，从而得到一个从原始数据集合到 N 维向量空间的映射——你总是需要显式地或者隐式地完成这样一个过程，然后基于某种规则进行分类，在该规则下，同组分类具有最大的相似性。

　　假设我们提取到原始数据的集合为(x₁, x₂, …, x_n)，并且每个xi为d维的向量，K-means聚类的目的就是，在给定分类组数k（k ≤ n）值的条件下，将原始数据分成k类
S = {S₁, S₂, …, S_k}，在数值模型上，即对以下表达式求最小值：
$\underset{\mathbf{S}} {\operatorname{arg\,min}} \sum_{i=1}^{k} \sum_{\mathbf x_j \in S_i} \left\| \mathbf x_j - \boldsymbol\mu_i \right\|^2$
这里μ_i 表示分类S_i的平均值。

　　那么在计算机编程中，其又是如何实现的呢？其算法步骤一般如下：

1、从D中随机取k个元素，作为k个簇的各自的中心。

2、分别计算剩下的元素到k个簇中心的相异度，将这些元素分别划归到相异度最低的簇。

3、根据聚类结果，重新计算k个簇各自的中心，计算方法是取簇中所有元素各自维度的算术平均数。

4、将D中全部元素按照新的中心重新聚类。

5、重复第4步，直到聚类结果不再变化。

6、将结果输出。

　　用数学表达式来说，

设我们一共有 N 个数据点需要分为 K 个 cluster ，k-means 要做的就是最小化

$\displaystyle J = \sum_{n=1}^N\sum_{k=1}^K r_{nk} \|x_n-\mu_k\|^2$

这个函数，其中 $r_{nk}$ 在数据点 n 被归类到 cluster k 的时候为 1 ，否则为 0 。直接寻找 $r_{nk}$ 和 $\mu_k$ 来最小化并不容易，不过我们可以采取迭代的办法：先固定 $\mu_k$ ，选择最优的 $r_{nk}$ ，很容易看出，只要将数据点归类到离他最近的那个中心就能保证最小。下一步则固定 $r_{nk}$ ，再求最优的 $\mu_k$ 。将对 $\mu_k$ 求导并令导数等于零，很容易得到最小的时候 $\mu_k$ 应该满足：

$\displaystyle \mu_k=\frac{\sum_n r_{nk}x_n}{\sum_n r_{nk}}$

亦即 $\mu_k$ 的值应当是所有 cluster k 中的数据点的平均值。由于每一次迭代都是取到的最小值，因此只会不断地减小（或者不变），而不会增加，这保证了 k-means 最终会到达一个极小值。虽然 k-means 并不能保证总是能得到全局最优解，但是对于这样的问题，像 k-means 这种复杂度的算法，这样的结果已经是很不错的了。

首先 3 个中心点被随机初始化，所有的数据点都还没有进行聚类，默认全部都标记为红色，如下图所示：

iter_00

然后进入第一次迭代：按照初始的中心点位置为每个数据点着上颜色，重新计算 3 个中心点，结果如下图所示：

iter_01

可以看到，由于初始的中心点是随机选的，这样得出来的结果并不是很好，接下来是下一次迭代的结果：

iter_02

可以看到大致形状已经出来了。再经过两次迭代之后，基本上就收敛了，最终结果如下：

iter_04

不过正如前面所说的那样 k-means 也并不是万能的，虽然许多时候都能收敛到一个比较好的结果，但是也有运气不好的时候会收敛到一个让人不满意的局部最优解，例如选用下面这几个初始中心点：

iter_00_bad

最终会收敛到这样的结果：

iter_03_bad

　　整体来讲，K-means算法的聚类思想比较简单明了，并且聚类效果也还算可以，算是一种简单高效应用广泛的 clustering 方法，接下来，我将讨论其代码实现过程。

聚类算法——K-means（下）

　　K-means的源码实现

　　一般情况下，我们通过C++/Matlab/Python等语言进行实现K-means算法，结合近期我刚刚学的C++，先从C++实现谈起，C++里面我们一般采用的是OpenCV库中写好的K-means函数，即cvKmeans2，首先来看函数原型：
　　从OpenCV manual看到的是：
int cvKMeans2(const CvArr* samples, int nclusters,
　　　　　　　 CvArr* labels, CvTermCriteria termcrit,
　　　　　　　 int attempts=1, CvRNG* rng=0,int flags=0,
　　　　　　　 CvArr* centers=0,double* compactness=0);
由于除去已经确定的参数，我们自己需要输入的为：
void cvKMeans2(
  　　const CvArr* samples, //输入样本的浮点矩阵，每个样本一行。
  　　int cluster_count,  //所给定的聚类数目
   　　* labels,    //输出整数向量：每个样本对应的类别标识
   　　CvTermCriteria termcrit //指定聚类的最大迭代次数和／或精度（两次迭代引起的聚类中心的移动距离）
);
其使用例程为：

 1 #ifdef _CH_
 2 #pragma package <opencv>
 3 #endif
 4 
 5 #define CV_NO_BACKWARD_COMPATIBILITY
 6 
 7 #ifndef _EiC
 8 #include "cv.h"
 9 #include "highgui.h"
10 #include <stdio.h>
11 #endif
12 
13 int main( int argc, char** argv )
14 {
15     #define MAX_CLUSTERS 5    //设置类别的颜色，个数（《=5）
16     CvScalar color_tab[MAX_CLUSTERS];
17     IplImage* img = cvCreateImage( cvSize( 500, 500 ), 8, 3 );
18     CvRNG rng = cvRNG(-1);
19     CvPoint ipt;
20 
21     color_tab[0] = CV_RGB(255,0,0);
22     color_tab[1] = CV_RGB(0,255,0);
23     color_tab[2] = CV_RGB(100,100,255);
24     color_tab[3] = CV_RGB(255,0,255);
25     color_tab[4] = CV_RGB(255,255,0);
26 
27     cvNamedWindow( "clusters", 1 );
28 
29     for(;;)
30     {
31         char key;
32         int k, cluster_count = cvRandInt(&rng)%MAX_CLUSTERS + 1;
33         int i, sample_count = cvRandInt(&rng)%1000 + 1;
34         CvMat* points = cvCreateMat( sample_count, 1, CV_32FC2 );
35         CvMat* clusters = cvCreateMat( sample_count, 1, CV_32SC1 );
36         cluster_count = MIN(cluster_count, sample_count);
37 
38         /** generate random sample from multigaussian distribution */
39         for( k = 0; k < cluster_count; k++ )
40         {
41             CvPoint center;
42             CvMat point_chunk;
43             center.x = cvRandInt(&rng)%img->width;
44             center.y = cvRandInt(&rng)%img->height;
45             cvGetRows( points, &point_chunk, k*sample_count/cluster_count,
46                        k == cluster_count - 1 ? sample_count :
47                        (k+1)*sample_count/cluster_count, 1 );
48 
49             cvRandArr( &rng, &point_chunk, CV_RAND_NORMAL,
50                        cvScalar(center.x,center.y,0,0),
51                        cvScalar(img->width*0.1,img->height*0.1,0,0));
52         }
53 
54         /** shuffle samples */
55         for( i = 0; i < sample_count/2; i++ )
56         {
57             CvPoint2D32f* pt1 = (CvPoint2D32f*)points->data.fl + cvRandInt(&rng)%sample_count;
58             CvPoint2D32f* pt2 = (CvPoint2D32f*)points->data.fl + cvRandInt(&rng)%sample_count;
59             CvPoint2D32f temp;
60             CV_SWAP( *pt1, *pt2, temp );
61         }
62 
63         printf( "iterations=%d\n", cvKMeans2( points, cluster_count, clusters,
64                 cvTermCriteria( CV_TERMCRIT_EPS+CV_TERMCRIT_ITER, 10, 1.0 ),
65                 5, 0, 0, 0, 0 ));
66 
67         cvZero( img );
68 
69         for( i = 0; i < sample_count; i++ )
70         {
71             int cluster_idx = clusters->data.i[i];
72             ipt.x = (int)points->data.fl[i*2];
73             ipt.y = (int)points->data.fl[i*2+1];
74             cvCircle( img, ipt, 2, color_tab[cluster_idx], CV_FILLED, CV_AA, 0 );
75         }
76 
77         cvReleaseMat( &points );
78         cvReleaseMat( &clusters );
79 
80         cvShowImage( "clusters", img );
81 
82         key = (char) cvWaitKey(0);
83         if( key == 27 || key == 'q' || key == 'Q' ) // 'ESC'
84             break;
85     }
86 
87     cvDestroyWindow( "clusters" );
88     return 0;
89 }
90 
91 #ifdef _EiC
92 main(1,"kmeans.c");
93 #endif

　　至于cvKmeans2函数的具体实现细节，可参见OpenCV源码

　　下面是Python的实现代码（网上所找）：

 1  #!/usr/bin/python
 2  
 3 from __future__ import with_statement
 4 import cPickle as pickle
 5 from matplotlib import pyplot
 6 from numpy import zeros, array, tile
 7 from scipy.linalg import norm
 8 import numpy.matlib as ml
 9 import random
10  
11 def kmeans(X, k, observer=None, threshold=1e-15, maxiter=300):
12     N = len(X)
13     labels = zeros(N, dtype=int)
14     centers = array(random.sample(X, k))
15     iter = 0
16  
17     def calc_J():
18         sum = 0
19         for i in xrange(N):
20             sum += norm(X[i]-centers[labels[i]])
21         return sum
22  
23     def distmat(X, Y):
24         n = len(X)
25         m = len(Y)
26         xx = ml.sum(X*X, axis=1)
27         yy = ml.sum(Y*Y, axis=1)
28         xy = ml.dot(X, Y.T)
29  
30         return tile(xx, (m, 1)).T+tile(yy, (n, 1)) - 2*xy
31  
32     Jprev = calc_J()
33     while True:
34         # notify the observer
35         if observer is not None:
36             observer(iter, labels, centers)
37  
38         # calculate distance from x to each center
39         # distance_matrix is only available in scipy newer than 0.7
40         # dist = distance_matrix(X, centers)
41         dist = distmat(X, centers)
42         # assign x to nearst center
43         labels = dist.argmin(axis=1)
44         # re-calculate each center
45         for j in range(k):
46             idx_j = (labels == j).nonzero()
47             centers[j] = X[idx_j].mean(axis=0)
48  
49         J = calc_J()
50         iter += 1
51  
52         if Jprev-J < threshold:
53             break
54         Jprev = J
55         if iter >= maxiter:
56             break
57  
58     # final notification
59     if observer is not None:
60         observer(iter, labels, centers)
61  
62 if __name__ == '__main__':
63     # load previously generated points
64     with open('cluster.pkl') as inf:
65         samples = pickle.load(inf)
66     N = 0
67     for smp in samples:
68         N += len(smp[0])
69     X = zeros((N, 2))
70     idxfrm = 0
71     for i in range(len(samples)):
72         idxto = idxfrm + len(samples[i][0])
73         X[idxfrm:idxto, 0] = samples[i][0]
74         X[idxfrm:idxto, 1] = samples[i][1]
75         idxfrm = idxto
76  
77     def observer(iter, labels, centers):
78         print "iter %d." % iter
79         colors = array([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
80         pyplot.plot(hold=False)  # clear previous plot
81         pyplot.hold(True)
82  
83         # draw points
84         data_colors=[colors[lbl] for lbl in labels]
85         pyplot.scatter(X[:, 0], X[:, 1], c=data_colors, alpha=0.5)
86         # draw centers
87         pyplot.scatter(centers[:, 0], centers[:, 1], s=200, c=colors)
88  
89         pyplot.savefig('kmeans/iter_%02d.png' % iter, format='png')
90  
91     kmeans(X, 3, observer=observer)

　　matlab的kmeans实现代码可直接参照其kmeans(X,k)函数的实现源码。

Matlab函数kmeans：K-均值聚类

转自：http://www.aiseminar.cn/html/93/t-693.html

K-means聚类算法采用的是将N*P的矩阵X划分为K个类，使得类内对象之间的距离最大,而类之间的距离最小。

使用方法：
Idx=Kmeans(X,K)
[Idx,C]=Kmeans(X,K)
[Idx,C,sumD]=Kmeans(X,K)
[Idx,C,sumD,D]=Kmeans(X,K)
[…]=Kmeans(…,’Param1’,Val1,’Param2’,Val2,…)

各输入输出参数介绍：

X N*P的数据矩阵
K 表示将X划分为几类，为整数
Idx N*1的向量，存储的是每个点的聚类标号
C K*P的矩阵，存储的是K个聚类质心位置
sumD 1*K的和向量，存储的是类间所有点与该类质心点距离之和
D N*K的矩阵，存储的是每个点与所有质心的距离

[…]=Kmeans(…,'Param1',Val1,'Param2',Val2,…)
这其中的参数Param1、Param2等，主要可以设置为如下：

1. ‘Distance’(距离测度)
‘sqEuclidean’ 欧式距离（默认时，采用此距离方式）
‘cityblock’ 绝度误差和，又称：L1
‘cosine’ 针对向量
‘correlation’ 针对有时序关系的值
‘Hamming’ 只针对二进制数据

2. ‘Start’（初始质心位置选择方法）
‘sample’ 从X中随机选取K个质心点
‘uniform’ 根据X的分布范围均匀的随机生成K个质心
‘cluster’ 初始聚类阶段随机选择10%的X的子样本（此方法初始使用’sample’方法）
matrix 提供一K*P的矩阵，作为初始质心位置集合

3. ‘Replicates’（聚类重复次数）整数

posted @ 2014-11-08 16:57 止战阅读(16066) 评论(0) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

K-means聚类算法

聚类算法——K-means（上）

聚类算法——K-means（下）

Matlab函数kmeans：K-均值聚类

公告