K-means算法的实现与应用举例

K-means 算法的实现与应用举例

1 K-means 方法

K-means 算法如下：

S1：初始化，聚类中心k c ,,c c 21，标号集 k I I I 21； S2： 分类：

end

;i I I ;

c x c x j n

i for j j T

j i j i k

j **1*min arg :1

S3：重新计算聚类中心：

end

;x I c k j for

j

I i i

j

j 1

:1

S4：迭代S2-S3，直至收敛。 其matlab 程序见附录1。

2实验

实验1 随机生成300个 44, 之间的二维数对，用K-means 算法将其分为两类（其matlab 程序见附录2），如fig1，由图1(b)可看出，离群点对聚类中心的位置有明显的影响。

实验2 随机生成600个二维数对，其中300个落于以（0,0）为圆心的单位圆，另外300

(a)

(b)

fig1 实验1

个落入（2,2）为圆心的单位圆，用K-means 算法将其分为两类（其matlab 程序见附录2），如fig2(a)，而fig2(b)则为在以上实验的基础上增加了30个干扰点之后的分类图，可见K-means 算法不能对其很好的分类，离群点对聚类中心的位置有明显的影响。

实验3 随机生成600个二维数对，其中300个落于以（0,0）为圆心的单位元，另外300个落入以（0,0）为圆心的单位圆，长半径为3短半径为2的圆盘，用K-means 算法将其分为2类（其matlab 程序见附录2），结果见fig3，可见K-means 算法同样不能对其很好的分类。

3 K-means 算法修正

修正一：实验2中增加离群点后，K-means 算法失效，是因为采用2范数计算距离，使计算出的重心与实际重心存在较大的误差。为减小误差，可以采用1-范数计算距离，或是采用中值代替均值计算新的聚类中心，即

k ,j ,

I i x medium c j i j 1

(a)

(b)

fig2 实验2

fig3 实验3

通过实验可以知道，采用1-范数计算距离实验效果并没有很好的改进，而采用中值计算聚类中心取得较好的效果（matlab 程序见附录3），采用同实验2增加干扰后相同的实验数据用修正后的K-means 算法进行分类，得到实验结果如fig4(a)，而实验3中结果产生的原因则是由于没有考虑数据点自身的结构特征与其他数据点之间关系引起，并且K-means 算法只考虑类内间距最小性并没有考虑类间间距的最大性，即只考虑了类内数据的相似性的最大性并没有考虑类间数据的差异性的最大性，所以单纯的改变聚类中心的选取方法，而没有对相关性（距离）进行本质的重新的定义，并不能对实验3的实验结果很好的改进，如fig4(b)：

4附录

附录1

function [idx,C,D,sumD]=kmeans_mean(X,k) %kmeans2norm K-means clustering. % X:n*p 的数据矩阵 % k:将X 划分为几类

% startM:k*p 的矩阵，初始类中心

% idx:n*1的向量，存储每个点的聚类编号 % C:k*p 的矩阵，存储k 个聚类的中心位置

% sumD:1*k 的和向量，类间所有点与该类中心的距离之和 % D:n*k 的矩阵，每一点与所有中心的距离

%sumK:1*k 的和向量,记录第k 个类中点的个数 [n,p] = size(X); idx=zeros(n,1); C=zeros(k,p); startM=zeros(k,p); D=zeros(n,k);

fig4 中值修正

(a)

(b)

(b)

(a)

%-----------随机生成初始聚类中心

for i=1:k

bi=ceil(i*n/k*rand);

startM(i,:)=X(bi,:);

end

while sum(abs(startM-C))>0

C=startM;startM(:)=0;

sumD=zeros(1,k);

sumK=zeros(1,k);%记录第k个类中点的个数

% count=zeros(1,k);%计数器

% sortC=zeros(k,n);

for i=1:n

%------计算每一点与所有中心的距离------

for j=1:k

D(i,j)=(X(i,:)-C(j,:))*(X(i,:)-C(j,:))';

% D(i,j)=sum(abs(X(i,:)-C(j,:)));

end

%-----------标号------

mini=inf;

for j=1:k

if D(i,j)

mini=D(i,j);

idx(i)=j;

end

end

sumD(idx(i))=sumD(idx(i))+D(i,idx(i));

sumK(idx(i))=sumK(idx(i))+1;

end

%---------计算新的聚类中心---------

%=======求解质心=======

for i=1:n

startM(idx(i),:)=startM(idx(i),:)+X(i,:)/sumK(idx(i));%求解质心end

end

附录2

clear all

% ------随机生成[-1,1]^2的数

% X = [randn(100,2)+ones(100,2);...

% randn(100,2)-ones(100,2)];

% % %------随机生成圆盘内的数

% deg1=(2*rand(300,1)-1)*2*pi;

% deg2=(2*rand(300,1)-1)*2*pi;

% r1=rand(300,1);

% r2=rand(300,1);

% X1=[r1.*cos(deg1) r1.*sin(deg1)];

% X2=[(r2+2).*cos(deg2) (r2+2).*sin(deg2)];

% X=zeros(600,2);

% for i=1:300

% X(2*i-1,:)=X1(i,:);

% X(2*i,:)=X2(i,:);

% end

% %------随机生成两个圆---------

deg1=(2*rand(300,1)-1)*2*pi;

r1=rand(300,1);

deg2=(2*rand(30,1)-1)*2*pi;

r2=0.1*rand(30,1);

deg3=(2*rand(300,1)-1)*2*pi;

r3=rand(300,1);

X1=[r1.*cos(deg1) r1.*sin(deg1)];

X2=[r2.*cos(deg2)+5 r2.*sin(deg2)+5];

X3=[r3.*cos(deg3)+1 r3.*sin(deg3)+1];

X=[X1;X3];

%X=[X1;X2;X3];

%[idx,ctrs,D,sumD]=kmeans_mean(X,2);

[idx,ctrs,D,sumD]=kmeans_medium(X,2);

figure

plot(X(idx==1,1),X(idx==1,2),'r.','MarkerSize',12)

hold on

plot(X(idx==2,1),X(idx==2,2),'b.','MarkerSize',12)

plot(ctrs(:,1),ctrs(:,2),'kx',...

'MarkerSize',12,'LineWidth',2)

plot(ctrs(:,1),ctrs(:,2),'ko',...

'MarkerSize',12,'LineWidth',2)

legend('Cluster 1','Cluster 2','Centroids',...

'Location','NW')

附录3

function [idx,C,D,sumD]=kmeans_medium(X,k)

%kmeans2norm K-means clustering.

% X:n*p的数据矩阵

% k:将X划分为几类

% startM:k*p的矩阵，初始类中心

% idx:n*1的向量，存储每个点的聚类编号

% C:k*p的矩阵，存储k个聚类的中心位置

% sumD:1*k的和向量，类间所有点与该类中心的距离之和% D:n*k的矩阵，每一点与所有中心的距离

%sumK:1*k的和向量,记录第k个类点的个数

[n,p] = size(X);

idx=zeros(n,1);

C=zeros(k,p);

startM=zeros(k,p);

D=zeros(n,k);

%-----------随机生成初始聚类中心

for i=1:k

bi=ceil(i*n/k*rand);

startM(i,:)=X(bi,:);

end

while sum(abs(startM-C))>0

C=startM;startM(:)=0;

sumD=zeros(1,k);

sumK=zeros(1,k);%记录第k个类中点的个数

% count=zeros(1,k);%计数器

% sortC=zeros(k,n);

for i=1:n

%------计算每一点与所有中心的距离------

for j=1:k

% D(i,j)=sqrt(sum((X(i,:)-C(j,:)).^2));

D(i,j)=(X(i,:)-C(j,:))*(X(i,:)-C(j,:))';

%D(i,j)=sum(abs((X(i,:)-C(j,:))));

end

%-----------标号------

mini=inf;

for j=1:k

if D(i,j)

mini=D(i,j);

idx(i)=j;

end

end

sumD(idx(i))=sumD(idx(i))+D(i,idx(i));

sumK(idx(i))=sumK(idx(i))+1;

end

%---------计算新的聚类中心---------

%=======质心作为聚类中心=======

% for i=1:n

% startM(idx(i),:)=startM(idx(i),:)+X(i,:)/sumK(idx(i));%求解质心% end

%=======中位数作为聚类中心=======

for i=1:k

inds=find(idx==i);

tempX=X(inds,:);

%-----------若标号为空则取离所在聚类中心最远的点为新的类中心------

if numel(tempX)==0

dist=0;

for jj=1:n

if D(jj,idx(jj))>dist

dist=D(jj,idx(jj));

end

end

sumK(idx(jj))=sumK(idx(jj))-1;

idx(jj)=i;

D(jj,i)=0;

sumK(i)=1;

tempX=X(jj,:);

end

for j=1:p

tempXj=tempX(:,j);

tempXj=sort(tempXj);

startM(i,j)=tempXj(round(sumK(i)/2));

end

end

end