// Classification.cpp: implementation of the Classification class.
// 2005.1.13－16：02 By Superman 
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
//#include "DSPLIT.h"
#include "Classification.h"
#include "matrix.h"

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

Classification::Classification()
{

}

Classification::~Classification()
{

}



/******************************************************************
*   函數(shù)名稱：LeastDistance()
*   函數(shù)類型：number_no，結(jié)構(gòu)體
*   函數(shù)功能：最小距離法 ,返回?cái)?shù)字類別和編號(hào)
******************************************************************/
number_no Classification::LeastDistance()
{
	double min=10000000000;
	number_no number_no;
	for(int n=0;n<10;n++)
	{
		for(int i=0;i<pattern[n].number;i++)
		{
			if(pipei(pattern[n].feature[i],testsample)<min)
			{
				//匹配的最小值
				min=pipei(pattern[n].feature[i],testsample);
				number_no.number=n;//樣品類別
				number_no.no=i;//樣品序號(hào)
			}
		}
	}
	return number_no;//返回手寫數(shù)字的類別和序號(hào)
}

/******************************************************************
*   函數(shù)名稱：Fisher_2Classes(int Class0, int Class1)
*   函數(shù)類型：int 
*   參數(shù)說明：Class0,Class1：0～9中的任意兩個(gè)類別
*   函數(shù)功能：兩類Fisher分類器,返回Class0，Class1中的一個(gè)
******************************************************************/
int Classification::Fisher_2Classes(int Class0, int Class1)
{
	double Xmeans[2][25];//兩類的均值
	double S[2][25][25];//樣品類內(nèi)離散度矩陣
	double Sw[25][25];//總類間離散度矩陣
	double Sw_[25][25];//Sw的逆矩陣
	double W[25];//解向量w*
	double difXmeans[25];//均值差
	double X[25];//未知樣品
	double m0,m1;//類樣品均值
	double y0;//閾值y0
	int i,j,k;

	for(i=0;i<2;i++)
		for(j=0;j<25;j++)
			Xmeans[i][j]=0;
	int num0,num1;		//兩類樣品的個(gè)數(shù)
	//兩類樣品特征
	double mode0[200][25],mode1[200][25];
	//兩類樣品個(gè)數(shù)
	num0=40;//pattern[Class0].number;
	num1=40;//pattern[Class1].number;
	for(i=0;i<num0;i++)
	{
		for(j=0;j<25;j++)
		{
			Xmeans[0][j]+=pattern[Class0].feature[i][j];
			mode0[i][j]=pattern[Class0].feature[i][j];
		}
	}

	for(i=0;i<num1;i++)
	{
		for(j=0;j<25;j++)
		{
			Xmeans[1][j]+=pattern[Class1].feature[i][j];	
			mode1[i][j]=pattern[Class1].feature[i][j];
		}
	}
	//求得兩個(gè)樣品均值向量
	for(i=0;i<25;i++)	
	{
		Xmeans[0][i]/=(double)num0;
		Xmeans[1][i]/=(double)num1;
	}
	//求兩類樣品類內(nèi)離散度矩陣
	for(i=0;i<25;i++)
	for(j=0;j<25;j++)
	{
		double s0=0.0,s1=0.0;
		for(k=0;k<num0;k++)
			s0=s0+(mode0[k][i]-Xmeans[0][i])*(mode0[k][j]-Xmeans[0][j]);
		s0=s0/(double)(num0-1);
		S[0][i][j]=s0;//第一類
		for(k=0;k<num1;k++)
			s1=s1+(mode1[k][i]-Xmeans[1][i])*(mode1[k][j]-Xmeans[1][j]);
		s1=s1/(double)(num1-1);
		S[1][i][j]=s1;//第二類		
	}
	//總類間離散度矩陣
	for(i=0;i<25;i++)
	for(j=0;j<25;j++)
	{
		Sw[i][j]=S[0][i][j]+S[1][i][j];
	}
	//Sw的逆矩陣
	for(i=0;i<25;i++)
		for(j=0;j<25;j++)
			Sw_[i][j]=Sw[i][j];	
	double(*p)[25]=Sw_;	
	brinv(*p,25);		//Sw的逆矩陣Sw_
	//計(jì)算w*  w*＝Sw_×(Xmeans0-Xmeans1)
	for(i=0;i<25;i++)
		difXmeans[i]=Xmeans[0][i]-Xmeans[1][i];
	for(i=0;i<25;i++)
		W[i]=0.0;
	brmul(Sw_,difXmeans,25,W);//計(jì)算出W*
	
	//各類樣品均值
	m0=0.0;
	m1=0.0;
	for(i=0;i<num0;i++)
	{
		m0+=brmul(W,mode0[i],25);
	}
	for(i=0;i<num1;i++)
	{
		m1+=brmul(W,mode1[i],25);
	}
	m0/=(double)num0;
	m1/=(double)num1;
	y0=(num0*m0+num1*m1)/(num0+num1);//閾值y0
	
	//對(duì)于任意的手寫數(shù)字X
	for(i=0;i<25;i++)
		X[i]=testsample[i];
	double y;//X在w*上的投影點(diǎn)
	y=brmul(W,X,25);
	if (y>=y0) 
		return Class0;
	else
		return Class1;
}

/******************************************************************
*   函數(shù)名稱：Fisher()
*   函數(shù)類型：int 
*   函數(shù)功能：Fisher分類器,返回手寫數(shù)字的類別
******************************************************************/
int Classification::Fisher()
{
	int i,j,number,maxval,num[10];
	for(i=0;i<10;i++)
		num[i]=0;
	for(i=0;i<10;i++)
		for(j=0;j<i;j++)
			num[Fisher_2Classes(i,j)]++;
	maxval=num[0];
	number=0;
	for(i=1;i<10;i++)
	{
		if(num[i]>maxval)
		{
			maxval=num[i];
			number=i;
		}
	}
	return number;
}

/******************************************************************
*   函數(shù)名稱：BayesErzhishuju()
*   函數(shù)類型：int 
*   函數(shù)功能：基于二值數(shù)據(jù)的Bayes分類器 ,返回手寫數(shù)字的類別
******************************************************************/
int Classification::BayesErzhishuju()
{
	double Pw[10];//先驗(yàn)概率P(wj)=Nj/N
	double P[10][25];//Pj(wi)wi:wi類，j:第j個(gè)特征
	double PXw[10];//類條件概率P(X|wj)
	double PwX[10];//后驗(yàn)概率P(wj|X)

	int i,j;

	//求先驗(yàn)概率
	int n[10];//各類樣品數(shù)
	int N=0;//樣品總數(shù)
	for(i=0;i<10;i++)
	{
		//各數(shù)字類別樣品數(shù)
		n[i]=pattern[i].number;
		N+=n[i];//樣品總數(shù)
	}
	for(i=0;i<10;i++)
		Pw[i]=(double)n[i]/(double)N;//先驗(yàn)概率
	
	//求類條件概率
	for(i=0;i<10;i++)
	{
		for(j=0;j<25;j++)
		{
			int numof1=0;//二值數(shù)據(jù)中1的個(gè)數(shù)
			for(int k=0;k<pattern[i].number;k++)
				numof1+=pattern[i].feature[k][j]>0.1?1:0;
			P[i][j]=(double)(numof1+1)/(double)(n[i]+2);
		}
	}

	for(i=0;i<10;i++)
	{
		double p=1;
		for(int j=0;j<25;j++)
		{
			p*=(testsample[j]>0.1)?P[i][j]:1-P[i][j];
		}
		PXw[i]=p;
	}

	//求后驗(yàn)概率
	double PX=0.0,maxP=0.0;
	int number;
	for(i=0;i<10;i++)
	{
		PX+=Pw[i]*PXw[i];
	}

	for(i=0;i<10;i++)
	{
		PwX[i]=Pw[i]*PXw[i]/PX;
		if(PwX[i]>maxP)
		{
			maxP=PwX[i];
			number=i;
		}
	}
	return number;
}

/******************************************************************
*   函數(shù)名稱：BayesLeasterror()
*   函數(shù)類型：int 
*   函數(shù)功能：最小錯(cuò)誤概率的Bayes分類器 ,返回手寫數(shù)字的類別
******************************************************************/
int Classification::BayesLeasterror()
{
	double X[25];//待測樣品
	double Xmeans[25];//樣品的均值
	double S[25][25];//協(xié)方差矩陣
	double S_[25][25];//S的逆矩陣
	double Pw;//先驗(yàn)概率
	double hx[10];//判別函數(shù)
	
	int i,j,k,n;
	
	for(n=0;n<10;n++)//循環(huán)類別0～9
	{
		int num=pattern[n].number;//樣品個(gè)數(shù)
		//求樣品平均值
		for(i=0;i<25;i++)
			Xmeans[i]=0.0;
		for(k=0;k<num;k++)
		{
			for(i=0;i<25;i++)
				Xmeans[i]+=pattern[n].feature[k][i]>0.10?1.0:0.0;
		}
		for(i=0;i<25;i++)
			Xmeans[i]/=(double)num;
		//求協(xié)方差矩陣
		double mode[200][25];
		for(i=0;i<num;i++)
			for(j=0;j<25;j++)
				mode[i][j]=pattern[n].feature[i][j]>0.10?1.0:0.0;
		for(i=0;i<25;i++)
		for(j=0;j<25;j++)
		{
			double s=0.0;
			for(k=0;k<num;k++)
				s=s+(mode[k][i]-Xmeans[i])*(mode[k][j]-Xmeans[j]);
			s=s/(double)(num-1);
			S[i][j]=s;
		}
		//求先驗(yàn)概率
		int total=0;
		for(i=0;i<10;i++)
			total+=pattern[i].number;
		Pw=(double)num/(double)total;
		//求S的逆矩陣
		for(i=0;i<25;i++)
			for(j=0;j<25;j++)
				S_[i][j]=S[i][j];

		double(*p)[25]=S_;

		brinv(*p,25);//S的逆矩陣
		//求S的行列式
		double (*pp)[25]=S;
		double DetS;
		DetS=bsdet(*pp,25);//S的行列式
		//求判別函數(shù)
		for(i=0;i<25;i++)
			X[i]=testsample[i]>0.10?1.0:0.0;
		for(i=0;i<25;i++)
			X[i]-=Xmeans[i];
		double t[25];
		for(i=0;i<25;i++)
			t[i]=0;
		brmul(X,S_,25,t);
		double t1=brmul(t,X,25);
		double t2=log(Pw);
		double t3=log(DetS+1);
		hx[n]=-t1/2+t2-t3/2;
	}

	double maxval=hx[0];
	int number=0;
	//判別函數(shù)的最大值
	for(n=1;n<10;n++)
	{
		if(hx[n]>maxval)
		{
			maxval=hx[n];
			number=n;
		}
	}
	return number;
}

/******************************************************************
*   函數(shù)名稱：BayesLeastRisk(double loss[10][10])
*   函數(shù)類型：double*
*   參數(shù)說明：double loss[10][10]：損失
*   函數(shù)功能：最小風(fēng)險(xiǎn)的Bayes分類器 ，返回各類的風(fēng)險(xiǎn)值
******************************************************************/
double* Classification::BayesLeastRisk(double loss[10][10])
{
	double X[25];//待測樣品
	double Xmeans[25];//樣品的均值
	double S[25][25];//協(xié)方差矩陣S
	double S_[25][25];//S的逆矩陣
	double P[10];//后驗(yàn)概率
	double Pw;//先驗(yàn)概率
	double hx[10];//判別函數(shù)
	
	int i,j,k,n;
	
	for(n=0;n<10;n++)//
	{
		int num=pattern[n].number;//樣品個(gè)數(shù)
		//求樣品均值
		for(i=0;i<25;i++)
			Xmeans[i]=0.0;
		for(k=0;k<num;k++)
		{
			for(i=0;i<25;i++)
				Xmeans[i]+=pattern[n].feature[k][i]>0.2?1.0:0.0;
		}
		for(i=0;i<25;i++)
			Xmeans[i]/=(double)num;
		//求協(xié)方差矩陣
		double mode[200][25];
		for(i=0;i<num;i++)
			for(j=0;j<25;j++)
				mode[i][j]=pattern[n].feature[i][j]>0.2?1.0:0.0;
		for(i=0;i<25;i++)
		for(j=0;j<25;j++)
		{
			double s=0.0;
			for(k=0;k<num;k++)
				s=s+(mode[k][i]-Xmeans[i])*(mode[k][j]-Xmeans[j]);
			s=s/(double)(num-1);
			S[i][j]=s;
		}
		//求先驗(yàn)概率
		int total=0;//樣品總數(shù)
		for(i=0;i<10;i++)
			total+=pattern[i].number;
		Pw=(double)num/(double)total;
		//求S的逆矩陣
		for(i=0;i<25;i++)
			for(j=0;j<25;j++)
				S_[i][j]=S[i][j];

		double(*p)[25]=S_;

		brinv(*p,25);//S的逆矩陣
		//求S的行列式
		double (*pp)[25]=S;
		double DetS;
		DetS=bsdet(*pp,25);//S的行列式
		//求判別函數(shù)
		for(i=0;i<25;i++)
			X[i]=testsample[i]>0.2?1.0:0.0;
		for(i=0;i<25;i++)
			X[i]-=Xmeans[i];
		double t[25];
		for(i=0;i<25;i++)
			t[i]=0;
		brmul(X,S_,25,t);
		double t1=brmul(t,X,25);
		double t2=log(Pw);
		double t3=log(DetS+1);
		P[n]=-t1/2+t2-t3/2;
	}
	
	for(n=0;n<10;n++)
	{
		double t=0.0;
		for(i=0;i<10;i++)
			t+=loss[n][i]*P[i];
		hx[n]=t;
	}
	
	return (double*)hx;
}

/******************************************************************
*   函數(shù)名稱：Jiangcheng()
*   函數(shù)類型：int 
*   函數(shù)功能：獎(jiǎng)懲算法 ,返回手寫數(shù)字的類別
******************************************************************/
int Classification::Jiangcheng()
{
	double w[10][26];//權(quán)值
	double d[10];//迭代計(jì)算結(jié)果權(quán)矢量
	double x[26];//增1樣品
	double hx[10];//判別函數(shù)
	bool flag;
	int n,i,j,k;
	int c=0;//迭代次數(shù)