stats.cc

标准的GP源代码,由Andy Singleton维护

// stats.ccc  version to support queue2.cc, from GPQUICK-2.1 chrome.cc
// W.Langdon cs.ucl.ac.uk $Revision: 1.48 $

//Modifications (reverse order):
// WBL 12 Feb 1997  GatherPStats so includes depths>=MAXTREEDEPTH
//                  DisplayPStats and DisplayLPStats to display it
//                  Doubtless output format of DisplayPStats could be better

// WBL 16 Feb 1997  Add costats

// WBL  7 Feb 1997  Add array bound checks (MAXTREEDEPTH isnt bigenought)

// WBL  1 Jan 1997  BUGFIX complete change of meaning implemented on 4-Apr-96
//                  so non-integer scores now work as they should
//                  Ensure nfit has same meaning as in integer scores case

// WBL  1 Jan 1997  Get DisplayDStats working for tournament sizes other than 4

// WBL 31 Dec 1996  BUGFIX avoid hashchain::~hashchain() calling itself

// WBL 31 Dec 1996  BUGFIX remove delete[] table from hashtable::include()
//                  Improve ERROR display in Pop::DisplayPStats

// WBL 16 Dec 1996  Make GatherPStats gather opcode data by tree depth

// WBL  4 Apr 1996  Make DisplayDStats assume integer scores and use blickle
//                  formula. NB covariance of nrank and nfit changed. 
//                  nrank nlonger supported.

// WBL 16 Mar 1996  Make DisplayDupStats record and display detailed stats

// WBL 15 Mar 1996  Add hashcode, Chrome::same and DisplayDupStats

// WBL 26 Feb 1996  Restore version for stack (from chrome.cxx r1.4) and 
//                  splice in using ifndef STACK_LIB (ie not gp.h)

//WBL 22 Sep 95  Replace list.h by prob.h

// WBL 23 May 1995  Add DisplayCStats

// WBL 12 May 1995  Add gathering problem dependant data

// WBL  8 May 1995  Output from DisplayLPStats >10% and make easier to find

// WBL  7 May 1995  Split GatherPStats from DisplayPStats and add 
//                  DisplayLPStats

// WBL  1 May 1995  Add DisplayPStats

// WBL 23 Apr 1995  Add Ntests_run

// WBL 20 Apr 1995  Add TIMINGS

// WBL 19 Feb 1995  Remove unneccessary output (to reduce disk usage)

// WBL 16 Feb 1995  Allow program scores to be of type scoretype

// WBL 12 Feb 1995  OK see if can get working now

// WBL 23 Jan 1995  Partial conversion to problem independance
//                  add gp.h, Primitives
//                  HACK to make COMPILE

// WBL 22 Nov 1994  Refit lost change.
//                  Change Badadf1 to adf1ver

// WBL 13 Sep 1994  Add displaying adf1

// WBL 13 Sep 1994  Make compilable by solaris gcc

// WBL 25 Aug 1994  Remove BestFitness check

// WBL 22 Aug 1994  Add queue2.h
//                  Take hits value from nfitness, rather than Chrome
//                  Add display of mem_used.

// WBL 15 Aug 1994  Add display of hits and total_hits

// WBL	6 Jul 1994  Fix bug whereby may not set last few values of nfit

// WBL  1 Jul 1994  New file based upon chrome.cxx Revision: 1.9
//                  Add pTraceDynamic.
//		    Display primative frequency, covarience and correllation

// WBL 28 Jun 1994  Make nfit take average of cube rather than linear average
//                  Tidy displays in DisplayDStats

// WBL 27 Jun 1994  Correct display of gencount in DisplayStats
//                  Add class stats and make DisplayStats use it.
//                  Complete impementation of DisplayDStats

// WBL 24 Jun 1994  Direct crossover to a particular tree
//                  Protect sqrt from rounding errors

// WBL 14 Jun 1994  Add support for pPopSeed, pTestSeed and PTrace

// WBL 14 Jun 1994  Add ChromeParams::Save() and ChromeParams::Load()
//                  based upon gpcplus3.0 cpv.cc.
//                  Add Pop::DisplayDStats()
//		    Add Pop::Write()

// WBL 16 May 1994  Add Pop::DisplayStats(). Fix MULTREE CrossOver bug

// WBL 12 May 1994  fix % by zero error in rnd() on SETNODE if varnum==0,
//		    Probably only arises if you have the wrong compiler switch

// WBL  5 May 1994  Make pMaxAge = 0 disable age limit,
//		    Add support for P-M random numbers

// GPQUICK
// C++ prefix stack implementation
// Divides GP system into classes:

#include <math.h>	//for DisplayStats
#include <assert.h>
#include "pch.h"
#include "chrome.h"
#include "prob.h"
#ifndef STACK_LIB
//#include "Primitives.h"
#ifdef QUEUE_LIB
#include "queue2.h"
#else
#include "gp.h"
#endif /*QUEUE_LIB*/
#endif

class costats;

class stats {						//WBL
	int cases;// = 0;
	float min;// = FLT_MAX;
	float max;// = -FLT_MAX;
	float sum;// = 0.0;
	double sqr;// = 0.0;

float	Mean() { return sum / cases; };
double	Variance();

public:
	stats(): cases(0),
	         min(FLT_MAX),
	         max(-FLT_MAX),
	         sum(0.0),
	         sqr(0.0){};

void 	include ( float in );
int	Cases() { return cases; };
float	Sum()   { return sum; };
void	display ( char* name, ostream & fout,
		  stats*,     stats*,         costats* );
friend class costats;
};//end class stats

class costats { //WBL
  stats x;
  stats y;
  double product;
public:
  costats(): product(0) {;}
void include(float xx, float yy){x.include(xx); y.include(yy); product+=xx*yy;}
double Covariance() {return product/x.Cases() - x.Mean()*y.Mean();}
double Correlation() {
  double d = sqrt(x.Variance() * y.Variance());
  if ( d > 0 )			//protect rounding error
    return Covariance() / d;
  else
    return FLT_MAX;
}
};//end class costats

double stats::Variance()					//WBL
{ 	double v = (sqr - double(cases*double(Mean())*double(Mean())) )/cases;
	if ( v >= 0 ) return v; else return 0.0; //protect against rounding
}

void stats::include( float in )  			//WBL
{	cases++;
	sum += in;
    	sqr += double(in) * double(in);
	if ( max < in ) max = in;
	if ( min > in ) min = in;
}//end stats::include


//NB covar stats assume primitive was present or wasnt cant use this code
//for covariance of continous things like length
void stats::display ( char* name, ostream & fout = cout,
                    stats* own_stats = NULL, stats* co_stats = NULL,
                    costats* covar = NULL )
                                                      //WBL
{
if ( cases > 0 )
    {	float mean     = Mean();
	double variance = Variance();
	fout << "Mean " << name << "= " << mean;
//	fout << "\tvar = " << variance;
	fout << "\tSD = "<<sqrt(variance);
	if ( covar != NULL)
	     { 	fout << "\tcovar " << covar->Covariance();
		fout << "\tcorel " << covar->Correlation();
	     };
	if ( own_stats != NULL && co_stats != NULL)
	     { 	fout << "\tcovar ";
//
//covariance = 	sum (x(i)*y(i))
//              ---------------  -  mean x * mean y
//               num cases
//
//In our case sum (x(i)*y(i)) = sum because we count y times. 
//The mean x is number of x (ie cases) divided by the popsize
//
		float covar;

		covar = (float)sum - (float) cases * co_stats->Mean();
		covar = covar/(float) own_stats->Cases();
		fout << covar;

// Correlleation co-efficient =             covariance
//				------------------------------------
//				sqrt( own variance * others variance
//
// Own variance includes cases where chrome contains none of the
// opcode, nb including those not in cases
//
// own variance = sqr / popsize - ( sum / popsize )**2
//
		{
		double d = sqrt( own_stats->Variance() * co_stats->Variance());
		if ( d > 0 )			//protect rounding error
		      { fout << "\tcorel ";
			fout << covar / d;
		      }
		} 
	     }
	fout << "\tmin = " << min;
	fout << "\tmax = " << max <<endl;
    }//endif
}//end stats::display

void Pop::DisplayStats(ostream & fout) const //WBL
{
	time_t now;
//	ChromeParams* params;
//	Problem* problem;
	time (&now);

	fout << ctime(&now);
#ifdef STACK_LIB
	fout << "popsize =\t"<< popsize  <<endl; 
	fout << "individuals =\t" <<gencount+1  <<endl;
	fout << "BestFitness =\t" <<*BestFitness  <<endl;
	fout << "BestMember =\t" <<BestMember  <<endl;
	fout << "initeval =\t" <<initeval  <<endl;
	fout << "nexteval =\t" <<nexteval  <<endl;
#endif
	stats fitness;
	stats length;
#ifdef PARETO
	stats hits [last_op+1];
	stats total_hits;
#endif
	stats memory;
#ifdef MINC_ADF
	stats adf1;
#endif
#ifdef TIMINGS
	stats cpu;
#endif
	stats Ntests;
	for (UINT i=0;i<popsize;i++)
	      { fitness.include (pop[i]->nfitness->fvalue);
		length.include (pop[i]->tree_starts[NUM_TREES-1]
                          + pop[i]->SubLen(pop[i]->tree_starts[NUM_TREES-1]));
#ifdef PARETO
		scoretype total = 0;
#endif
#ifndef STACK_LIB
                                                             //force right type
		ptrmyfitnessvalue temp = ptrmyfitnessvalue (pop[i]->nfitness);
#ifdef PARETO
		for (int j=0; j<=last_op; j++)
			{ total += temp->hits[j];
			  hits[j].include (temp->hits[j]);
		        }
		total_hits.include (total);
#endif
		memory.include (temp->mem_used);
#ifdef MINC_ADF
		adf1.include(temp->adf1);
#endif
#ifdef TIMINGS
		cpu.include(temp->cpu_time);
#endif
#ifndef QUEUE_LIB
		Ntests.include(temp->Ntests_run);
#endif /*QUEUE_LIB*/
#endif
	      }
                                                      
	fitness.display( "fitness", fout );
	length.display( "length ", fout );

#ifdef PARETO
	for (int j=0; j<=last_op; j++) 
		{ problem->WriteTreeName(j, fout);
		  fout << "\t";
		  hits[j].display("", fout);
	        }
	fout << "total\t";
	total_hits.display( "", fout );
#endif
	fout << "memory\t";
	memory.display( "", fout );
#ifdef MINC_ADF
	fout << "adf1ver\t";
	adf1.display( "", fout );
#endif
#ifdef TIMINGS
	fout << "cpu\t";
	cpu.display( "", fout );
#endif
	fout << "tests\t";
	Ntests.display( "", fout );

}//end Pop::DisplayStats

#define QARGNUM(op) (problem->funclist[op]->argnum)
void Pop::DisplayDStats(ostream & fout ) const //WBL
//
// Convert chromosome nfitness into a rank, and then into an averag
// rank and then into an approximate chance of selection. Cf roulette
// wheel selection.
//    Gather and display loads of statistics by tree/opcode
//
{
//float nfit [popsize];    //equivelent roulette wheel
float* nfit = new float [popsize];    //equivelent roulette wheel

const int MaxDepth = (params->params[pMaxDepth]!=0)? 
                         params->params[pMaxDepth] : MAXTREEDEPTH;
assert(MaxDepth<=MAXTREEDEPTH);

#ifdef INT_SCORE
{
int score[max_fitness+1];
float nr4[max_fitness+1];
memset(score,0,sizeof(score));

for (int i = 0; i < popsize; i++ ) score[int(pop[i]->nfitness->fvalue)]++;

assert (params->params[pTournSize] == 4);

int sum = 0;
float r4 = 0;

for (i = 0; i <= max_fitness; i++ ) {
	if(score[i] != 0) {
		float old  = r4;
		sum       += score[i];
		float r    = float (sum) / float (popsize);
		r4         = r*r*r*r;
		nr4[i]     = (r4 - old)*float(popsize)/float(score[i]);
	}
}//end for i
for (UINT j = 0; j<popsize; j++) nfit[j] = nr4[int(pop[j]->nfitness->fvalue)];
}
UINT i;
#else
//Chrome* rank [popsize];
//float nrank [popsize];   //normalised rank
Chrome** rank = new Chrome* [popsize];
//float* nrank  = new float   [popsize];   //normalised rank
UINT* pop_index = new UINT [popsize];

//set rank by sorting population. Sort not efficient but good enough
for (UINT i = 0; i < popsize; i++ ) {rank [i] = pop[i]; pop_index[i] = i;}
UINT j;
for ( i = 0; i < popsize; i++ )
      { BOOL swapflag = FALSE;
	for (j = 1; j < popsize; j++ )
	      { if ( rank[j-1]->nfitness->IsBetter(rank[j]->nfitness) )
		      { Chrome* r  = rank [j-1];
			rank [j-1] = rank [j];
			rank [j]   = r;
		        int index       = pop_index [j-1];
			pop_index [j-1] = pop_index [j];
			pop_index [j]   = index;
			swapflag   = TRUE;
                      }
	      }//end loop j
	if (swapflag == FALSE) break;
      }//end loop i

//Nolonger true WBL 25/8/94
//assert( BestFitness->fvalue == rank[popsize-1]->nfitness->fvalue ); //check

//take large population approximation
#define fitnessofrank(rank) (params->params[pTournSize]* \
               pow(double(rank)/double(popsize),params->params[pTournSize]-1))

UINT k = 0;
double sum_normalised_fitness = fitnessofrank(1);

UINT l;
for ( j = 1; j < popsize; j++ )
      { //float J = j;
	//float K = k;
	if ( rank[j]->nfitness->IsBetter(rank[k]->nfitness))
	      { 
//		nrank [ k ] = (J*(J-1)-(K-1)*K)/(2*(J-K));
		nfit  [ pop_index[k] ] = sum_normalised_fitness/(j-k);
		for ( l = k+1; l<j; l++ )
		      { //nrank[ l ] = nrank [ k ];
			nfit [ pop_index[l] ] = nfit  [ pop_index[k] ];
		      }
		k = j;
		sum_normalised_fitness = 0.0;
	      }//end j is not same rank as k
	sum_normalised_fitness += fitnessofrank(j+1);
      }

for ( l = k; l<popsize; l++ )
	      { //float J = popsize; //NB include j
		//float K = k;
//		nrank [ l ] = (J*(J-1)-(K-1)*K)/(2*(J-K));
		nfit  [ pop_index[l] ] = sum_normalised_fitness/(popsize-k);
	      }
//delete [] nrank;
delete [] pop_index;
delete [] rank;

#endif

stats length [NUM_TREES];
costats length_nfit [NUM_TREES];
costats length_fit [NUM_TREES];
struct x
{ stats frequency; /* stats length; stats children;
  stats fitness;   stats nrank;*/  stats nfit;     }; 
x summary;
//x counts [NUM_TREES][problem->funccount];

#define funccount_max 32

x counts [NUM_TREES][funccount_max][MAXTREEDEPTH+1];
assert (problem->funccount<=funccount_max);

for ( i = 0; i < popsize; i++ )
{
//summary.frequency JUNK
//summary.length    JUNK
//summary.children.	include( pop[i]->num_children );
//summary.fitness.	include( pop[i]->nfitness->fvalue );
//summary.nrank.		include( nrank [i] );
summary.nfit.		include( nfit [i] );

int len;
for ( int t = 0; t<NUM_TREES; t++)
      { len = pop[i]->SubLen(pop[i]->tree_starts[t]);
	length[t].  include ( len );
	length_nfit[t].include (len, nfit[i]);
	length_fit[t].include  (len, pop[i]->nfitness->fvalue);

//	int f [problem->funccount];
	int f [funccount_max][MAXTREEDEPTH+1];
	memset(f,0,sizeof(f));
        int tree_stack [EXPRLEN];
        tree_stack[0] = -1;
	int sp = 0;

	{for ( int n = pop[i]->tree_starts[t];
              n < (pop[i]->tree_starts[t] + len);
	      n++ 				)
	{
	        while (tree_stack[sp] == 0) { //remove finished branches
			tree_stack[--sp] -= 1;
		}
		++sp;
		tree_stack[sp] = QARGNUM(pop[i]->expr[n].op);
		if(sp<MAXTREEDEPTH) f[pop[i]->expr[n].op][sp]++;
		else                f[pop[i]->expr[n].op][MAXTREEDEPTH]++;

//counts[t][pop[i]->expr[n].op].length.  include ( len );
//counts[t][pop[i]->expr[n].op].children.include ( pop[i]->num_children );
//counts[t][pop[i]->expr[n].op].fitness. include ( pop[i]->nfitness->fvalue );
//counts[t][rank[i]->expr[n].op].nrank.   include ( nrank[i] );
counts[t][pop[i]->expr[n].op][sp].nfit.    include ( nfit[i] );

	}};//end for each op in this tree

	{for ( int op = 0; op < problem->funccount; op++ )
	   for ( int d = 1; d <= MaxDepth; d++ )
	      { counts[t][op][d].frequency.include ( f[op][d] ); }}

      };//end for each tree
};//end for whole population