声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1815|回复: 0

[编程技巧] KMeans和KMedoid 的Matlab实现

[复制链接]
发表于 2016-5-6 10:12 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
KMeans和KMedoid算法是聚类算法中比较普遍的方法,本文讲了其原理和matlab中实现的代码。

1.目标:

       找出一个分割,使得距离平方和最小


2.K-Means算法:

       1. 将数据分为k个非空子集

       2. 计算每个类中心点(k-means中用所有点的平均值,K-medoid用离该平均值最近的一个点)center

       3. 将每个object聚类到最近的center

       4. 返回2,当聚类结果不再变化的时候stop


   复杂度:

       O(kndt)

       -计算两点间距离:d

       -指定类:O(kn)   ,k是类数

       -迭代次数上限:t


3.K-Medoids算法:

       1. 随机选择k个点作为初始medoid

       2.将每个object聚类到最近的medoid

       3. 更新每个类的medoid,计算objective function

       4. 选择最佳参数

       4. 返回2,当各类medoid不再变化的时候stop


    复杂度:

       O((n^2)d)

       -计算各点间两两距离O((n^2)d)

       -指定类:O(kn)   ,k是类数


4.特点:

       -聚类结果与初始点有关(因为是做steepest descent from a random initial starting oint)

       -是局部最优解

       -在实际做的时候,随机选择多组初始点,最后选择拥有最低TSD(Totoal Squared Distance)的那组




Kmeans KMedoid Implementation with matlab:

===================

下面是我用matlab上的实现:
说明:fea为训练样本数据,gnd为样本标号。算法中的思想和上面写的一模一样,在最后的判断accuracy方面,由于聚类和分类不同,只是得到一些 cluster ,而并不知道这些 cluster 应该被打上什么标签,或者说。由于我们的目的是衡量聚类算法的 performance ,因此直接假定这一步能实现最优的对应关系,将每个 cluster 对应到一类上去。一种办法是枚举所有可能的情况并选出最优解,另外,对于这样的问题,我们还可以用 Hungarian algorithm 来求解。具体的Hungarian代码我放在了资源里,调用方法已经写在下面函数中了。下面给出Kmeans&Kmedoid主函数。

Kmeans.m 函数:

  1. <font color="#000000">[cpp] view plaincopyprint?
  2. function [ accuracy,MIhat ] = KMeans( K,mode )  
  3.   
  4. % Artificial Intelligence & Data Mining - KMeans & K-Medoids Clustering  
  5. % Author: Sophia_qing [url=home.php?mod=space&uid=156850]@[/url] ZJU  
  6. % CreateTime: 2012-11-18  
  7. % Function: Clustering  
  8. %  -K: number of clusters  
  9. %  -mode:   
  10. %   1: use kmeans cluster algorithm in matlab  
  11. %   2: k_medroid algorithm: use data points as k centers  
  12. %   3: k_means algorithm: use average as k centers  
  13.   
  14. global N_features;  
  15. global N_samples;  
  16. global fea;  
  17. global gnd;  
  18.   
  19. switch (mode)  
  20.     case 1 &#65533;ll system function KMeans  
  21.         label = kmeans(fea,K);  
  22.         [label,accuracy] = cal_accuracy(gnd,label);  
  23.          
  24.     case 2%use kmedroid method  
  25.         for testcase = 1:10% do 10 times to get rid of the influence from Initial_center  
  26.             K_center = Initial_center(fea,K); %select initial points randomly  
  27.             changed_label = N_samples;  
  28.             label = zeros(1,N_samples);  
  29.             iteration_times = 0;  
  30.             while changed_label~=0  
  31.                 cls_label = cell(1,K);  
  32.                 for i = 1: N_samples  
  33.                     for j = 1 : K  
  34.                         D(j) = dis(fea(i,:),K_center(j,:));  
  35.                     end  
  36.                     [~,label(i)] = min(D);  
  37.                     cls_label{label(i)} = [cls_label{label(i)} i];  
  38.                 end  
  39.                 changed_label = 0;  
  40.                 cls_center = zeros(K,N_features);  
  41.                 for i = 1 : K  
  42.                     cls_center(i,:) = mean(fea(cls_label{i},:));  
  43.                     D1 = [];  
  44.                     for j = 1:size(cls_label{i},2)%number of samples clsutered in i-th class  
  45.                         D1(j) = dis(cls_center(i,:),fea(cls_label{i}(j),:));  
  46.                     end  
  47.                     [~,min_ind] = min(D1);  
  48.                     if ~isequal(K_center(i,:),fea(cls_label{i}(min_ind),:))  
  49.                         K_center(i,:) = fea(cls_label{i}(min_ind),:);  
  50.                         changed_label = changed_label+1;  
  51.                     end  
  52.                 end  
  53.                 iteration_times = iteration_times+1;  
  54.             end  
  55.             [label,acc(testcase)] = cal_accuracy(gnd,label);  
  56.         end  
  57.         accuracy = max(acc);  
  58.          
  59.     case 3%use k-means method  
  60.         for testcase = 1:10% do 10 times to get rid of the influence from Initial_center  
  61.             K_center = Initial_center(fea,K); %select initial points randomly  
  62.             changed_label = N_samples;  
  63.             label = zeros(1,N_samples);  
  64.             label_new = zeros(1,N_samples);  
  65.             while changed_label~=0  
  66.                 cls_label = cell(1,K);  
  67.                 changed_label = 0;  
  68.                 for i = 1: N_samples  
  69.                     for j = 1 : K  
  70.                         D(j) = dis(fea(i,:),K_center(j,:));  
  71.                     end  
  72.                     [~,label_new(i)] = min(D);  
  73.                     if(label_new(i)~=label(i))  
  74.                         changed_label = changed_label+1;  
  75.                     end;  
  76.                     cls_label{label_new(i)} = [cls_label{label_new(i)} i];  
  77.                 end  
  78.                 label = label_new;  
  79.                   
  80.                 for i = 1 : K  %recalculate k centroid  
  81.                     K_center(i,:) = mean(fea(cls_label{i},:));  
  82.                 end  
  83.             end  
  84.              [label,acc(testcase)] = cal_accuracy(gnd,label);  
  85.         end  
  86.         accuracy = max(acc);  
  87. end  
  88.   
  89. MIhat = MutualInfo(gnd,label);  
  90.   
  91.   
  92.     function center = Initial_center(X,K)  
  93.         rnd_Idx = randperm(N_samples,K);  
  94.         center = X(rnd_Idx,:);  
  95.     end  
  96.   
  97.     function res = dis(X1,X2)  
  98.         res = norm(X1-X2);  
  99.     end  
  100.   
  101.     function [res,acc] = cal_accuracy(gnd,estimate_label)  
  102.         res = bestMap(gnd,estimate_label);  
  103.         acc = length(find(gnd == res))/length(gnd);  
  104.     end  
  105. end  

  106. 实验结果分析:
  107. 对上面得到的accuracy进行画图,横坐标为10个数据集,纵坐标为在其上进行聚类的准确率。
  108. 其中,auto为matlab内部kmeans函数。
  109. 画图:

  110. [cpp] view plaincopyprint?
  111. function [  ] = Plot( A,B,C )  
  112. %PLOT Summary of this function goes here  
  113. %   Detailed explanation goes here  
  114. figure;  
  115. k = 1:10;  
  116. plot(k,A,'-r',k,B,'-b',k,C,'-g');  
  117. legend('auto','medoid','means');  
  118.   
  119.   
  120. end  
  121. </font>
复制代码

结果:

5类聚类:

1.jpg

7类聚类:

2.jpg

转自:http://blog.sina.com.cn/s/blog_6c41e2f30101c48o.html






回复
分享到:

使用道具 举报

您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-15 12:11 , Processed in 0.076873 second(s), 21 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表