聚类算法-Kmeans算法的简单实现

Posted 菜鸟升级

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了 聚类算法-Kmeans算法的简单实现相关的知识，希望对你有一定的参考价值。

1. 聚类与分类的区别：

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

 

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

　　假设我们提取到原始数据的集合为(x1, x2, …, xn)，并且每个xi为d维的向量，K-means聚类的目的就是，在给定分类组数k（k ≤ n）值的条件下，将原始数据分成k类 
S = {S1, S2, …, Sk}，在数值模型上，即对以下表达式求最小值：
\\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 表示分类Si 的平均值。

说明：给定K，聚类完成后，对于每个类，对该类中的所有样本点到该类中心点的距离的平方求和，记为sum，对每个类都这样求和，然后把所有类的

sum加起来，总和记为SUM，聚类Kmeans算法就是要使得SUM最小。

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

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 。即上面的J即为最上面的SUM，直接寻找 r_{nk} 和 \\mu_k 来最小化 J 并不容易，不过我们可以采取迭代的办法：先固定 \\mu_k ，选择最优的 r_{nk} ，很容易看出，只要将数据点归类到离他最近的那个中心就能保证 J 最小。下一步则固定 r_{nk}，再求最优的 \\mu_k。将 J 对 \\mu_k 求导并令导数等于零，很容易得到 J 最小的时候 \\mu_k 应该满足：

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

其中 r_{nk} 在数据点 n 被归类到 cluster k 的时候为 1，即上面的分母就是类cluster k中样本点数，而分子为该分类中所有样本点的和

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

 

2. 实现：

注：该实现用的是opencv中的一个函数：cvKMeans2函数

代码来源：http://www.cnblogs.com/moondark/archive/2012/03/08/2385870.html

 

[cpp] view plain copy
 
 在CODE上查看代码片派生到我的代码片
  1. </pre><p></p><p style="margin:15px auto 10px; padding-top:0px; padding-bottom:0px; line-height:22px; text-align:justify"></p><pre code_snippet_id="500422" snippet_file_name="blog_20141029_2_1547356" name="code" class="cpp" style="font-family: \'Lucida Grande\', \'Lucida Sans Unicode\', Helvetica, Arial, Verdana, sans-serif; font-size: 14px; ">#ifdef _CH_  
  2. #pragma package <opencv>  
  3. #endif  
  4.   
  5. #define CV_NO_BACKWARD_COMPATIBILITY  
  6.   
  7. // opencv函数手册：http://wiki.opencv.org.cn/index.php/Cxcore数组操作#.E9.9A.8F.E6.9C.BA.E6.95.B0.E7.94.9F.E6.88.90  
  8. // 代码来源：http://www.cnblogs.com/moondark/archive/2012/03/08/2385870.html  
  9. #include "opencv/cv.h"  
  10. #include "opencv/highgui.h"  
  11. #include <stdio.h>  
  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 ); //8位，3通道  
  18.     CvRNG rng = cvRNG(-1); //4位的有符号数，若4位全为1，则值为-1，即A=1111,A=-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; //聚类的类别数随机生成，在[1,5]范围内  
  33.         int i, sample_count = cvRandInt(&rng)%1000 + 1; //样本个数  
  34.         CvMat* points = cvCreateMat( sample_count, 1, CV_32FC2 ); //生成sample_count行，1列且数据类型为32位浮点双通道的数组，样本点数组  
  35.         CvMat* clusters = cvCreateMat( sample_count, 1, CV_32SC1 );//生成sample_count行，1列且数据类型为32位有符号数单通道的数组，每个样本对应的类别标识数组  
  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.             /* 
  46.             GetRow, GetRows 
  47.             返回数组的一行或在一定跨度内的行 
  48.             CvMat* cvGetRow( const CvArr* arr, CvMat* submat, int row ); 
  49.             CvMat* cvGetRows( const CvArr* arr, CvMat* submat, int start_row, int end_row, int delta_row=1 ); 
  50.             arr 
  51.             输入数组。 
  52.             submat 
  53.             指向返回的子数组头的指针。 
  54.             row 
  55.             被选定行的索引下标，索引下标从0开始。 
  56.             start_row 
  57.             跨度的开始行（包括此行）索引下标，索引下标从0开始。 
  58.             end_row 
  59.             跨度的结束行（不包括此行）索引下标，索引下标从0开始。 
  60.             delta_row 
  61.             在跨度内的索引下标跨步，从开始行到结束行每隔delta_row行提取一行。 
  62.             */  
  63.             cvGetRows( points, &point_chunk, k*sample_count/cluster_count,  
  64.                 k == cluster_count - 1 ? sample_count :  
  65.                 (k+1)*sample_count/cluster_count, 1 );  
  66.   
  67.             /* 
  68.             RandArr 
  69.             用随机数填充数组并更新 RNG 状态 
  70.             void cvRandArr( CvRNG* rng, CvArr* arr, int dist_type, CvScalar param1, CvScalar param2 ); 
  71.             rng 
  72.             被 cvRNG 初始化的 RNG 状态. 
  73.             arr 
  74.             输出数组 
  75.             dist_type 
  76.             分布类型: 
  77.             CV_RAND_UNI - 均匀分布 
  78.             CV_RAND_NORMAL - 正态分布 或者 高斯分布 
  79.             param1 
  80.             分布的第一个参数。如果是均匀分布它是随机数范围的闭下边界。如果是正态分布它是随机数的平均值。 
  81.             param2 
  82.             分布的第二个参数。如果是均匀分布它是随机数范围的开上边界。如果是正态分布它是随机数的标准差。 
  83.             */  
  84.             cvRandArr( &rng, &point_chunk, CV_RAND_NORMAL,  
  85.                 cvScalar(center.x,center.y,0,0),  
  86.                 cvScalar(img->width*0.1,img->height*0.1,0,0));  
  87.         }  
  88.   
  89.         /** shuffle samples 即随机打乱样本*/   
  90.         for( i = 0; i < sample_count/2; i++ )  
  91.         {  
  92.             CvPoint2D32f* pt1 = (CvPoint2D32f*)points->data.fl + cvRandInt(&rng)%sample_count;  
  93.             CvPoint2D32f* pt2 = (CvPoint2D32f*)points->data.fl + cvRandInt(&rng)%sample_count;  
  94.             CvPoint2D32f temp;  
  95.             CV_SWAP( *pt1, *pt2, temp );  
  96.         }  
  97.   
  98.         // 使用Kmeans聚类  
  99.         /* 
  100.         int cvKMeans2(const CvArr* samples, int nclusters, 
  101.         　　　　　　　 CvArr* labels, CvTermCriteria termcrit, 
  102.                 　　　　　　　 int attempts=1, CvRNG* rng=0,int flags=0,  
  103.                         　　　　　　　 CvArr* centers=0,double* compactness=0); 
  104.                                 由于除去已经确定的参数，我们自己需要输入的为： 
  105.                                 void cvKMeans2(  
  106.                                 　　const CvArr* samples, //输入样本的浮点矩阵，每个样本一行。  
  107.                                   　　int cluster_count,  //所给定的聚类数目  
  108.                                     　　* labels,    //输出整数向量：每个样本对应的类别标识  
  109.                                       　　CvTermCriteria termcrit //指定聚类的最大迭代次数和／或精度（两次迭代引起的聚类中心的移动距离） 
  110.                                         );  
  111.         */  
  112.         // 输出迭代数  
  113.         printf( "iterations=%d\\n", cvKMeans2( points, cluster_count, clusters,  
  114.             cvTermCriteria( CV_TERMCRIT_EPS+CV_TERMCRIT_ITER, 10, 1.0 ),  
  115.             5, 0, 0, 0, 0 ));  
  116.         cvZero( img );  
  117.   
  118.         for( i = 0; i < sample_count; i++ ) // 依照样本所属的类别，以不同的颜色画出样本点  
  119.         {  
  120.             int cluster_idx = clusters->data.i[i];  
  121.             ipt.x = (int)points->data.fl[i*2];  
  122.             ipt.y = (int)points->data.fl[i*2+1];  
  123.             cvCircle( img, ipt, 2, color_tab[cluster_idx], CV_FILLED, CV_AA, 0 );  
  124.         }  
  125.   
  126.         cvReleaseMat( &points );  
  127.         cvReleaseMat( &clusters );  
  128.   
  129.         cvShowImage( "clusters", img );  
  130.   
  131.         key = (char) cvWaitKey(0);  
  132.         if( key == 27 || key == \'q\' || key == \'Q\' ) // \'ESC\'  
  133.             break;  
  134.     }  
  135.   
  136.     cvDestroyWindow( "clusters" );  
  137.     return 0;  
  138. }  

效果如下：

 

 

转：http://blog.csdn.net/damotiansheng/article/details/40582211

以上是关于 聚类算法-Kmeans算法的简单实现的主要内容，如果未能解决你的问题，请参考以下文章

spark.mllib源码阅读-聚类算法1-KMeans

求MATLAB实现canopy-kmeans聚类算法的完整代码

算法笔记：Kmeans聚类算法简介

基于R语言的Kmeans聚类算法

Kmeans 聚类 及其python实现

详解 kmeans 聚类算法