K均值聚类（Kmeans）
Sigma = [1, 0; 0, 1]; mu1 = [1, -1]; x1 = mvnrnd(mu1, Sigma, 200); mu2 = [5.5, -4.5]; x2 = mvnrnd(mu2, Sigma, 200); mu3 = [1, 4]; x3 = mvnrnd(mu3, Sigma, 200); mu4 = [6, 4.5]; x4 = mvnrnd(mu4, Sigma, 200); mu5 = [9, 0.0]; x5 = mvnrnd(mu5, Sigma, 200); % obtain the 1000 data points to be clustered X = [x1; x2; x3; x4; x5]; % Show the data point plot(x1(:,1), x1(:,2), 'r.'); hold on; plot(x2(:,1), x2(:,2), 'b.'); plot(x3(:,1), x3(:,2), 'k.'); plot(x4(:,1), x4(:,2), 'g.'); plot(x5(:,1), x5(:,2), 'm.');
save myX %将X存储到文件中，在其他文件中load就可以了
结果如下：
1 % 初始聚类中心 2 mu0=[X(1,:); X(205,:);X(405,:);X(605,:);X(805,:)]; 3 mu1=zeros(5,c); 4 lable=zeros(r,1); 5 % first cluster 6 dist=zeros(5,1); 7 for k=1:r 8 for n=1:5 9 dist(n)=norm(X(k,:)-mu0(n,:)); 10 end 11 lable(k)=find(dist==min(dist)); 12 end 13 14 % X1=[lable X]; 15 sum=zeros(5,2); 16 count=zeros(5,1); 17 18 % 第一次聚类 19 for n=1:5 20 for k=1:r 21 if lable(k)==n 22 sum(n,:)=sum(n,:)+X(k,:); 23 count(n)=count(n)+1; 24 end 25 end 26 mu1(n,:)=sum(n,:)/count(n); 27 end 28 29 % square error 30 e=zeros(5,1); 31 esum=0; 32 for n=1:5 33 for k=1:r 34 if lable(k)==n 35 e(n,:)=e(n,:)+norm(X(k,:)-mu1(n,:)); 36 end 37 end 38 esum=esum+e(n,:); 39 end 40 41 % ***************** start Iteration ************** 42 43 countNotX=0;% 到达一定值后，停止迭代 44 s=0; 45 count1=zeros(5,1); 46 while true 47 s=s+1; 48 %random choose X to update 49 ran=round(1+999*rand(1));%总共1000个样本 50 Xrand=X(ran,:); 51 rou=zeros(5,1); 52 rouMin=inf; 53 newClassLable=0; 54 if count(lable(ran))~=0 55 for n=1:5 56 if lable(ran)==n 57 rou(n)=norm(Xrand-mu1(n,:))*count(n)/(count(n)-1); 58 else 59 rou(n)=norm(Xrand-mu1(n,:))*count(n)/(count(n)+1); 60 end 61 if rou(n)<=rouMin; 62 rouMin=rou(n); 63 newClassLable=n; 64 end 65 end 66 if rouMin<rou(lable(ran)) 67 countNotX=0; 68 % new mu 这里用课件中的公式，不知道哪里错了，总是不行，就用最蠢的办法了 69 for n=1:5 70 for k=1:r 71 if lable(k)==n 72 sum(n,:)=sum(n,:)+X(k,:); 73 count1(n)=count1(n)+1; 74 end 75 end 76 mu1(n,:)=sum(n,:)/count1(n); 77 end 78 % new count 79 count(lable(ran))=count(lable(ran))-1; 80 count(newClassLable)=count(newClassLable)+1; 81 % new esum 82 e(lable(ran),:)=e(lable(ran),:)-rou(lable(ran)); 83 e(newClassLable,:)=e(newClassLable,:)+rouMin; 84 esum=esum+rouMin-rou(lable(ran)); 85 lable(ran)=newClassLable ; 86 else 87 % disp('No') 88 countNotX= countNotX+1; 89 if countNotX==600 90 %when J is not changed for 600 times,reckon OK 91 break 92 end 93 end 94 end 95 end
