isodata.cpp

模式识别 ISODATA算法

// ISODATA.cpp : Defines the entry point for the console application.
//
#include <iostream.h>
#include <stdio.h>
#include <math.h>
#define  N 10
#define  eps  0.00001
struct Pointf//模式特征向量
{
	int s;
	float x;
	float y;
};
struct PointZ //聚类中心坐标
{
	float x;
	float y;
};

float CalDistancef(Pointf x1,Pointf x2)//计算模式间的距离
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}

float CalDistanceZ(PointZ x1,PointZ x2)//计算类内距离
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}

float CalDistancefZ(Pointf x1,PointZ x2)//计算类间距离
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}

int main(int argc, char* argv[])
{
	Pointf pts[N]={
		{0,0.0,0.0},{1,3.0,8.0},{2,2.0,2.0},{3,1.0,1.0},{4,5.0,3.0},
		{5,4.0,8.0},{6,6.0,3.0},{7,5.0,4.0},{8,6.0,4.0},{9,7.0,5.0},
	};
	int i,j,m;
	cout<<"样本集为:"<<"\n";
        for(i=0;i<N;i++)
        {
           cout<<"X["<<i<<"]=("<<pts[i].x<<","<<pts[i].y<<"),";
		if((i+1)%5==0)
		{
		   cout<<"\n";
		}
        }
        cout<<"\n";
	int Nc=0;
	cout<<"please input Nc(1-10): ";
	cin>>Nc;
    int Z[N];//选取的初始聚类中心序号
	for(i=0;i<Nc;i++)
	{
		cout<<"输入初始第"<<i<<"聚类中心的序号(0-9):";
		cin>>Z[i];
	}
	
    int    Nj[N]; //记录每个类中元素的个数
	PointZ ZArray[N];  //记录每个聚类中心坐标
	Pointf SAArray[N][N]; //记录每类中模式坐标
	float  DjAv[N];
	float  Deltaj[N][2]; //记录每类的标准差矢量
	float  Deltajmax[N];
	int    DeltajmaxCor[N];
	float  DAv;
	int    Nreal=N;
	int    count=0;   //迭代次数
	float  Dij[N*N/2];
	int    Diji[N];
	int    Dijj[N];
	int    q=0;
	int    p=0;
    float ft;
	int   it;
	int   jt;
	int   flag;
	int ss=0;
	PointZ Ztp;
	PointZ ZArraytp[N];
	int    Nctp;
	char ch;
	int cur=0;
	for(i=0;i<N;i++)
	{
		Nj[i]=0; 
	}

	//选取初始聚类中心
	for(i=0;i<Nc;i++)
	{
		int ihere=Z[i];
		ZArray[i].x=pts[ihere].x;
		ZArray[i].y=pts[ihere].y;
	}

	int K,ThetaN;
	float ThetaS,ThetaC;
	int L,I;

Step1://预置
	cout<<"输入预期聚类中心数目K:";
	cin>>K;
	cout<<"输入每个聚类域中最少的样本数ThetaN:";
	cin>>ThetaN;
	cout<<"输入同一聚类域中样本标准差的最大值ThetaS:";
	cin>>ThetaS;
	cout<<"输入不同聚类域距离最小值ThetaC:";
	cin>>ThetaC;
	cout<<"输入一次可以合并的聚类中心的最多对数L:";
	cin>>L;
	cout<<"输入最大迭代次数I:";
	cin>>I;
Step2:
	for(i=0;i<Nc;i++)
	{
	    Nj[i]=0;
	}
    cout<<"\n";
	cout<<"这是第"<<count+1<<"次归类"<<"\n";	
	for(i=0;i<N;i++)                          //将模式样本归类
	{
	     if(pts[i].s==-1)continue;      //若该点的序号为－1则说明它是被剔除的
	     float ftemp=1.0e+10;
	     int xx=0;
	     float dis;
         for(j=0;j<Nc;j++)                     
	     {
	    	dis=CalDistancefZ(pts[i],ZArray[j]);
		    if(dis<ftemp||fabs(dis-ftemp)<eps)
			{
		      xx=j;
		      ftemp=dis;
			}
	     }
	     SAArray[xx][Nj[xx]].x=pts[i].x;
	     SAArray[xx][Nj[xx]].y=pts[i].y;
	     SAArray[xx][Nj[xx]].s=pts[i].s;	
	     Nj[xx]=Nj[xx]+1;	
	}
	count++;
Step3:
	for(j=0;j<Nc;j++)               //样本数少于ThetaN的类合并
	{
		if(Nj[j]<ThetaN)
		{
			for(i=0;i<Nj[j];i++)
			{
				pts[SAArray[j][i].s].s=-1;
			}
			i=j;
			int tr=j;
			Nreal-=Nj[j];
			while(j<Nc-1)
			{
				
				for(m=0;m<Nj[j+1];m++)
				{
					SAArray[j][m].x=SAArray[j+1][m].x;
					SAArray[j][m].y=SAArray[j+1][m].y;
					SAArray[j][m].s=SAArray[j+1][m].s;
				} 
				j++;			
			}

			while(i<Nc-1)
			{
				Nj[i]=Nj[i+1];
				i++;
			}
			Nc--;
			j=tr;
		}		
	}

Step4:              //修正各聚类中心
	for(j=0;j<Nc;j++)
	{
		float temx=0,temy=0;
		for(i=0;i<Nj[j];i++)
		{
			temx+=SAArray[j][i].x;
			temy+=SAArray[j][i].y;		
		}
		ZArray[j].x=temx/Nj[j];
		ZArray[j].y=temy/Nj[j];
	}

    for(i=0;i<Nc;i++)
	{
		   cout<<"第"<<i<<"个聚类中心是: "<<"("<<ZArray[i].x<<","<<ZArray[i].y<<")"<<endl;
		   cout<<"    包含的元素有: ";
		   for(j=0;j<Nj[i];j++)
		   {
			  cout<<"X["<<SAArray[i][j].s<<"],";
		   }
		   cout<<"\n";
		  
	 }

Step5://计算各聚类域中诸样本与聚类中心的平均距离
	float temp=0.0;
	for(j=0;j<Nc;j++)
	{
		for(i=0;i<Nj[j];i++)
		{
			temp+=CalDistancefZ(SAArray[j][i],ZArray[j]);
		}

		DjAv[j]=temp/Nj[j];
		
		temp=0.0;
	}
Step6://计算全部模式样本对应聚类中心的总平均距离
	DAv=0;
	for(j=0;j<Nc;j++)
	{
		DAv+=Nj[j]*DjAv[j];
	}
	DAv/=Nreal;

Step7:

	if(count>=I) goto Step14;

	if(Nc<=K/2)goto Step8;

	if((count%2==0)||Nc>=2*K) goto Step11;
Step8://计算各聚类中样本距离标准差
	
	for(j=0;j<Nc;j++)
	{
		float temx=0.0,temy=0.0;
		for(i=0;i<Nj[j];i++)
		{
           temx+=(SAArray[j][i].x-ZArray[j].x)*(SAArray[j][i].x-ZArray[j].x);
		   temy+=(SAArray[j][i].y-ZArray[j].y)*(SAArray[j][i].y-ZArray[j].y);
		}
        Deltaj[j][0]=sqrtf(temx/Nj[j]);
		Deltaj[j][1]=sqrtf(temy/Nj[j]);	
		temx=0.0,temy=0.0;	
	}	
Step9://求每个标准差向量中的最大分量 
	for(j=0;j<Nc;j++)
	{
		Deltajmax[j]=Deltaj[j][0]>Deltaj[j][1]?Deltaj[j][0]:Deltaj[j][1];
		DeltajmaxCor[j]=Deltaj[j][0]>Deltaj[j][1]?0:1;
	}
Step10://分裂判断和计算
	for(j=0;j<Nc;j++)
	{
		if(Deltajmax[j]>ThetaS)
		{
			if((DjAv[j]>DAv&&Nj[j]>2*(ThetaN+1))||Nc<=K/2)
			{
				float Garma=0.5;
				PointZ Zj1,Zj2;

				if(DeltajmaxCor[j]==0)
				{
					Zj1.x=ZArray[j].x+Deltajmax[j]*Garma;
					Zj1.y=ZArray[j].y;
					Zj2.x=ZArray[j].x-Deltajmax[j]*Garma;
					Zj2.y=ZArray[j].y;
				}
				else if(DeltajmaxCor[j]==1)
				{
					Zj1.x=ZArray[j].x;
					Zj1.y=ZArray[j].y+Deltajmax[j]*Garma;
					Zj2.x=ZArray[j].x;
					Zj2.y=ZArray[j].y-Deltajmax[j]*Garma;
				}
				
                ZArray[j].x=Zj1.x;
				ZArray[j].y=Zj1.y;

				ZArray[Nc].x=Zj2.x;
				ZArray[Nc].y=Zj2.y;
				Nc++;
				goto Step2;
			}
		}
	}

Step11://计算全部聚类中心的距离
	ss=0;
	for(i=0;i<Nc-1;i++)
	{
		for(j=i+1;j<Nc;j++)
		{
			Dij[ss]=CalDistanceZ(ZArray[i],ZArray[j]);
			Diji[ss]=i;
			Dijj[ss]=j;
			ss++;		
		}
	}
		
Step12: //简单起见，只考虑一次只合并一对聚类中心的情况
	    //找出类间距离最小的
	ft=Dij[0];
	it=Diji[0];
	jt=Dijj[0];
	for(i=1;i<ss;i++)
	{
		if(Dij[i]<ft)
		{
			ft=Dij[i];
			it=Diji[i];
			jt=Dijj[i];
		}
	}
Step13:
	if(ft<ThetaC)
	{
		Ztp.x=(Nj[it]*ZArray[it].x+Nj[jt]*ZArray[jt].x)/(Nj[it]+Nj[jt]);
		Ztp.y=(Nj[it]*ZArray[it].y+Nj[jt]*ZArray[jt].y)/(Nj[it]+Nj[jt]);
		ZArray[it].x=Ztp.x;
		ZArray[jt].y=Ztp.y;

		j=jt;
		while(j<Nc-1)
		{
			ZArray[j].x=ZArray[j+1].x;
			ZArray[j].y=ZArray[j+1].y;
			j++;
		}
		j=jt;
		
		while(j<Nc-1)
		{
			Nj[j]=Nj[j+1];
			j++;
		}
		Nc--;
	}
Step14:
	if(count>=I) 
	{
		count=0;	
		cout<<"共分为"<<Nc<<"类"<<"\n";
		return 0;
		
	}
	else
	{
		goto Step2;	
	}
}