clear
clc
kl
s='y';ii=0;ccc=0;
while s~='n'&ccc<10;
ii=ii+1;ccc=ccc+1
clear si sil si2
data1=xlsread('data1.xls');
data2=xlsread('data2.xls');
data3=xlsread('data3.xls');%从workspace读入3类50组数据以备分析.
num=randperm(150);
data=[data1;data2;data3];%xlsread('Y.xls');%%取所有数据为训练数据.
n=length(data);
X=data(num,:);%随机排列训练数据.
si=(num<=50)+((num<=100)&(num>50))*2+(num>100)*3;%标记待训练数据中各数据的类别状态,x∈datai,标为i.
%--------------准备训练数据并标记-----------------
K=3;%input('Please input the K:'); %期望的聚类数.
Min_N=13;%input('Please input the min_N:');%聚类中的最少样本数.
Min_S=.56;%input('Please input the min_S:');%标准偏差参数.
Min_C=1.2;%input('Please input the min_C:');%两聚类中心的最小距离.
L=1;%input('Please input the L:'); %每次迭代允许合并的最大聚类对数.
I0=100;%input('Please input the I:'); %允许迭代的次数.
c=5;%input('Please input the c:'); %初始聚类数.
a0=.5;%input('Please inout the a:'); %分裂距离因子,0<a<1.
m=X(1:c,:);m0=[m';si(1:c)]'; %初始聚类中心.
%--------------设置初始数据及控制参数(1)-------------
count=0;%迭代次数计数器.
%--------------开始ISODATA迭代算法----------------
N=1;N0=0;m1=zeros(c,4);min_size=c;
tic
while ((N(1:min(length(N),length(N0)))~=N0(1:min(length(N),length(N0))))|sum(sum(m(1:min_size,:)~=m1(1:min_size,:))))&count<I0
N0=N;m1=m;
clear d N av_d dd;
cc=(1:c);
for i=1:n
for j=1:c
d(j)=norm(X(i,:)-m(j,:));%计算欧式距离.
end
f=find(d==min(d));
sil(i)=f(1);%分类后的数据类标识,作用同si.
end
%%----------以最小欧式距离将所有样本数据分类(2)------------
for i=1:c
N(i)=length(find(sil==i));
end
%N
N_th=find(N<Min_N);
if isempty(N_th)~=1
for i=1:length(N_th)
TH=find(sil==N_th(i));
sil(TH)=0; %更新后的分类标识.
end
end
c=c-length(N_th); %更新后的聚类数.
cc(N_th)=[];
clear m;
%%------------取消基数小于Min_N的类(3)---------------
for i=1:c
ff=find(sil==cc(i));
m(i,:)=mean(X(ff,:));%计算新的聚类中心.
for j=1:length(ff)
d1(j)=norm(X(ff(j),:)-m(i,:));
end
av_d(i)=mean(d1);%计算各类中样本至中心的平均距离.
end
av_d;
%%--------计算更新的聚类中心和各类样本的平均距离av_d(4),(5)------
to_av_d=(sum(av_d.*N(cc)))/(sum(N(cc)));
%%--------计算所有样本距其相应中心的平均距离to_av_d,(6)--------
clear S
if (c<=K/2)| ((mod(count+1,2)==1)&(c<=2*K))
for i=1:c
Y=X((find(sil==cc(i))),:);
clear mm;
for j=1:length(Y)
mm(j,:)=m(i,:);
end
S(i,:)=sqrt(mean((Y-mm).^2));%每一类各分量的标准偏差S.
end
S_max=max(S'); %每一类中偏差最大的分量S_max
%%----计算每一类各分量的标准偏差s及每一类中偏差最大的分量S_max,(8),(9)---
f3=find(S_max>Min_S);
nu=c;
if isempty(f3)==0
for i=1:length(f3)
if (nu<=K/2)|(N(f3(i))>2*(Min_N+1))&(av_d(f3(i))>to_av_d)
m(c+1,:)=m(f3(i),:)+a0*S_max(f3(i));
m(f3(i),:)=m(f3(i),:)-a0*S_max(f3(i));
N(c+1)=fix(N(f3(i))/2);
N(f3(i))=fix(N(f3(i))/2)+rem(N(f3(i)),2);
c=c+1; %类分裂后新的聚类个数.
end
end
m; %类分裂后新的聚类中心.
end
%%---------------分裂过程(10)------------------
end
%%---------------判断是否需要考虑分裂(7)-------------
for i=1:c-1
for j=i+1:c
dd(i,j)=norm(m(i,:)-m(j,:)); %两中心间距.
end
end
%%------------对所有中心,求两两间距离dd.(11)-------------
clear f4;
if isempty(find(dd>0&dd<Min_C))==0
[f4(:,1),f4(:,2)]=find(dd>0&dd<Min_C);
d2=find((dd<Min_C)&(dd>0));
f5=sortrows([(dd(d2))';(f4(:,1))';(f4(:,2))']');
l=min(L,length(f4(:,1)));
for i=1:l
m(f5(i,2),:)=(m(f5(i,2),:)*N(f5(i,2))+m(f5(i,3),:)*N(f5(i,3)))/(N(f5(i,2))+N(f5(i,3)));
N(f5(i,2))=N(f5(i,2))+N(f5(i,3));
end
i=1:l;
m(f5(i,3),:)=[];%合并后新的中心向量.
N(f5(i,3))=[] %合并后新的类元素个数.
siz_m=size(m);
c=siz_m(1);%合并后新的类个数.
end
%%--------比较dd和Min_C,合并后求新中心和聚类个数.(12)(13)-------
count=count+1;
m_size=size(m);
m1_size=size(m1);
min_size=min(m_size(1),m1_size(1));
end
time(ii)=toc;
I(ii)=count;
%--------------ISODATA迭代算法结束----------------
sil
clear a d6 b N0 m1
for i=1:c
f10=find(sil==i);
f11=si(f10);
for j=1:c
n1(j)=length(find(f11==j));
end
n2=find(n1==max(n1));
a(i,n2(1))=1;
end
a; %分类标识转换矩阵.
[d6(:,1),d6(:,2)]=find(a==1);
for i=1:length(d6(:,1))
b=find(sil==d6(i,1));
sil(b)=d6(i,2)+10;
N0(d6(i,2))=N(d6(i,1));
m1(d6(i,2),:)=m(d6(i,1),:);
end
a_size=size(a);
if a_size(1)<length(N)
N0(a_size(1)+1:length(N))=N(a_size(1)+1:length(N));
end
m=m1;
zero_N=N0==0;
if zero_N<3
N0(zero_N)=[];
m(zero_N,:)=[];
end
%clear N
if length(N0)<3
Num(ii,:)=[N0,zeros(1,3-length(N0))]; %最终各类中元素个数.
else
Num(ii,:)=N0(1:3);
end
sil=sil-10;
%--------------将最终分类结果sil对照si做标识-----------
error_rate(ii)=1-((sum(si==sil))/n); %计算分类错误率.
%M=[mean(data1);mean(data2);mean(data3)]; %理想的各类中心.
M=[mean(data(1:50,:));mean(data(51:100,:));mean(data(101:150,:))]; %理想的各类中心.
for i=1:n
clear d f
for j=1:3
d(j)=norm(X(i,:)-M(j,:));%以理想中心计算欧式距离.
end
f=find(d==min(d));
si2(i)=f(1);%以理想中心分类后的数据类标识.
end
s_error_rate(ii)=1-(sum(si==si2))/n;%理想中心分类错误率.
%m_size=size(m);
%for i=1:m_size(1)
% clear d;
% F=find(sil==i);
% XX=X(F,:);
% for j=1:length(F)
% d(i,j)=(norm(XX(j,:)-m(i,:)))^2;
% end
%end
%J(ii)=sum(sum(d)); %计算分类误差平方和.
%s_J(ii)=sum(sum(dd));%理想中心分类的误差平方和.
%--------------基于最小欧式距离的聚类分析结束-----------
s=a;%input('Anykey to continue,Space to exit.','s')
end
m
Num
I
time
error_rate
%J
%plot(si,'or');title('实际数据分类分布');xlabel('数据次序'),ylabel('类别');
%grid %描绘实际数据分类分布.
%figure(2)
%plot(sil,'or');title('ISODATA算法分类结果');xlabel('数据次序'),ylabel('类别');
%grid %描绘ISODATA算法分类结果.
%figure(3)
%plot3(data1(:,1),data1(:,3),data1(:,4),'or',data2(:,1),data2(:,3),data2(:,4),'ok',data3(:,1),data3(:,3),data3(:,4),'og');
%title('实际数据的三维分量空间分布');xlabel('第一维'),ylabel('第三维'),zlabel('第四维')
%grid %描绘实际数据的三维分量空间分布.
%figure(4)
%f_N=find(Num(:,3)~=0);
%KK=.02:.02:ccc;
K=1:ccc;%(f_N');
%plot(K,Num);title('各类中数据个数');xlabel('a'),ylabel('个数');
%grid %描绘K取不同值时所分类各类中数据个数.
%figure(5)
%plot(K,I);title('迭代次数');xlabel('a'),ylabel('次数');
%grid %描绘K取不同值时所需的迭代次数.
%figure(6)
s_error=0.0667+zeros(1,length(K));
%s_J=s_error-0.067+43.4850;
%plot(K,error_rate,'b',K,s_error,'r');title('错误率');xlabel('a'),ylabel('error_rate');
grid %描绘K取不同值时的分类错误率.
%figure(7)
%plot(K,J,'b',K,s_J,'r');title('误差平方和');xlabel('a'),ylabel('J');
%grid %描绘K取不同值时的误差平方和.
%figure(8)
plot(K,time);title('消耗时间');xlabel('a'),ylabel('time/second');
%grid %描绘K不同时ISODATA算法分类的时间分布.
%figure(9);title('分类错误和正确的个数');xlabel('a'),ylabel('个数');
%hist(J)
%grid %描绘在100次实验中由于聚类中心初�