GeneticPopulation.cxx

Go to the documentation of this file.
00001 // @(#)root/tmva $Id: GeneticPopulation.cxx 31458 2009-11-30 13:58:20Z stelzer $    
00002 // Author: Peter Speckmayer
00003 
00004 /**********************************************************************************
00005  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis       *
00006  * Package: TMVA                                                                  *
00007  * Class  : TMVA::GeneticPopulation                                               *
00008  * Web    : http://tmva.sourceforge.net                                           *
00009  *                                                                                *
00010  * Description:                                                                   *
00011  *      Implementation (see header for description)                               *
00012  *                                                                                *
00013  * Authors (alphabetical):                                                        *
00014  *      Peter Speckmayer <speckmay@mail.cern.ch> - CERN, Switzerland              *
00015  *                                                                                *
00016  * Copyright (c) 2005:                                                            *
00017  *      CERN, Switzerland                                                         *
00018  *      MPI-K Heidelberg, Germany                                                 *
00019  *                                                                                *
00020  * Redistribution and use in source and binary forms, with or without             *
00021  * modification, are permitted according to the terms listed in LICENSE           *
00022  * (http://tmva.sourceforge.net/LICENSE)                                          *
00023  **********************************************************************************/
00024 
00025 #include <iostream>
00026 #include <iomanip>
00027 
00028 #include "Rstrstream.h"
00029 #include "TSystem.h"
00030 #include "TRandom3.h"
00031 #include "TH1.h"
00032 #include <algorithm>
00033 
00034 #include "TMVA/GeneticPopulation.h"
00035 #include "TMVA/GeneticGenes.h"
00036 #include "TMVA/MsgLogger.h"
00037 
00038 ClassImp(TMVA::GeneticPopulation)
00039 
00040 using namespace std;
00041    
00042 //_______________________________________________________________________
00043 //                                                                      
00044 // Population definition for genetic algorithm                          
00045 //_______________________________________________________________________
00046 
00047 //_______________________________________________________________________
00048 TMVA::GeneticPopulation::GeneticPopulation(const std::vector<Interval*>& ranges, Int_t size, UInt_t seed) 
00049    : fGenePool(size),
00050      fRanges(ranges.size()),
00051      fLogger( new MsgLogger("GeneticPopulation") )
00052 {
00053    // Constructor
00054    
00055    // create a randomGenerator for this population and set a seed
00056    // create the genePools
00057    //
00058    fRandomGenerator = new TRandom3( 100 ); //please check
00059    fRandomGenerator->Uniform(0.,1.);
00060    fRandomGenerator->SetSeed( seed );
00061 
00062    for ( unsigned int i = 0; i < ranges.size(); ++i )
00063       fRanges[i] = new TMVA::GeneticRange( fRandomGenerator, ranges[i] );
00064 
00065    vector<Double_t> newEntry( fRanges.size() );
00066    for ( int i = 0; i < size; ++i )
00067       {
00068          for ( unsigned int rIt = 0; rIt < fRanges.size(); ++rIt )
00069             newEntry[rIt] = fRanges[rIt]->Random();
00070          fGenePool[i] = TMVA::GeneticGenes( newEntry);
00071       }
00072 
00073    fPopulationSizeLimit = size;
00074 }
00075 
00076 //_______________________________________________________________________
00077 TMVA::GeneticPopulation::~GeneticPopulation() 
00078 {
00079    // destructor
00080    if (fRandomGenerator != NULL) delete fRandomGenerator;
00081 
00082    std::vector<GeneticRange*>::iterator it = fRanges.begin();
00083    for (;it!=fRanges.end(); it++) delete *it;
00084 
00085    delete fLogger;
00086 }
00087 
00088 
00089 
00090 //_______________________________________________________________________
00091 void TMVA::GeneticPopulation::SetRandomSeed( UInt_t seed ) 
00092 {
00093    // the random seed of the random generator
00094    fRandomGenerator->SetSeed( seed );
00095 }
00096 
00097 //_______________________________________________________________________
00098 void TMVA::GeneticPopulation::MakeCopies( int number )
00099 {
00100    // produces offspring which is are copies of their parents
00101    // Parameters:
00102    //         int number : the number of the last individual to be copied
00103    //
00104    
00105    int i=0; 
00106    for (std::vector<TMVA::GeneticGenes>::iterator it = fGenePool.begin();
00107         it != fGenePool.end() && i < number; 
00108         ++it, ++i ) {
00109       GiveHint( it->GetFactors(), it->GetFitness() );
00110    }
00111 }
00112 
00113 //_______________________________________________________________________
00114 void TMVA::GeneticPopulation::MakeChildren()
00115 {
00116    // does what the name says,... it creates children out of members of the
00117    // current generation
00118    // children have a combination of the coefficients of their parents
00119    //
00120 
00121 #ifdef _GLIBCXX_PARALLEL
00122 #pragma omp parallel
00123 #pragma omp for
00124 #endif
00125    for ( int it = 0; it < (int) (fGenePool.size() / 2); ++it )
00126       {
00127          Int_t pos = (Int_t)fRandomGenerator->Integer( fGenePool.size()/2 );
00128          fGenePool[(fGenePool.size() / 2) + it] = MakeSex( fGenePool[it], fGenePool[pos] );
00129       }
00130 }
00131 
00132 //_______________________________________________________________________
00133 TMVA::GeneticGenes TMVA::GeneticPopulation::MakeSex( TMVA::GeneticGenes male, 
00134                                                      TMVA::GeneticGenes female ) 
00135 {
00136    // this function takes two individuals and produces offspring by mixing (recombining) their
00137    // coefficients
00138    //
00139    vector< Double_t > child(fRanges.size());
00140    for (unsigned int i = 0; i < fRanges.size(); ++i) {
00141       if (fRandomGenerator->Integer( 2 ) == 0) {
00142          child[i] = male.GetFactors()[i];
00143       }else{
00144          child[i] = female.GetFactors()[i];
00145       }
00146    }
00147    return TMVA::GeneticGenes( child );
00148 }
00149 
00150 //_______________________________________________________________________
00151 void TMVA::GeneticPopulation::Mutate( Double_t probability , Int_t startIndex, 
00152                                       Bool_t near, Double_t spread, Bool_t mirror ) 
00153 {
00154    // mutates the individuals in the genePool
00155    // Parameters:
00156    //         double probability : gives the probability (in percent) of a mutation of a coefficient
00157    //         int startIndex : leaves unchanged (without mutation) the individuals which are better ranked
00158    //                     than indicated by "startIndex". This means: if "startIndex==3", the first (and best)
00159    //                     three individuals are not mutaded. This allows to preserve the best result of the 
00160    //                     current Generation for the next generation. 
00161    //         Bool_t near : if true, the mutation will produce a new coefficient which is "near" the old one
00162    //                     (gaussian around the current value)
00163    //         double spread : if near==true, spread gives the sigma of the gaussian
00164    //         Bool_t mirror : if the new value obtained would be outside of the given constraints
00165    //                    the value is mapped between the constraints again. This can be done either
00166    //                    by a kind of periodic boundary conditions or mirrored at the boundary.
00167    //                    (mirror = true seems more "natural")
00168    //
00169 
00170    vector< Double_t>::iterator vec;
00171    vector< TMVA::GeneticRange* >::iterator vecRange;
00172 
00173    //#ifdef _GLIBCXX_PARALLEL
00174    // #pragma omp parallel
00175    // #pragma omp for
00176    //#endif
00177    // The range methods are not thread safe!
00178    for (int it = startIndex; it < (int) fGenePool.size(); ++it) {
00179       vecRange = fRanges.begin();
00180       for (vec = (fGenePool[it].GetFactors()).begin(); vec < (fGenePool[it].GetFactors()).end(); ++vec) {
00181          if (fRandomGenerator->Uniform( 100 ) <= probability) {
00182             (*vec) = (*vecRange)->Random( near, (*vec), spread, mirror );
00183          }
00184          ++vecRange;
00185       }
00186    }
00187 }
00188 
00189 
00190 //_______________________________________________________________________
00191 TMVA::GeneticGenes* TMVA::GeneticPopulation::GetGenes( Int_t index )
00192 {
00193    // gives back the "Genes" of the population with the given index.
00194    //
00195    return &(fGenePool[index]);
00196 }
00197 
00198 //_______________________________________________________________________
00199 void TMVA::GeneticPopulation::Print( Int_t untilIndex )
00200 {
00201    // make a little printout of the individuals up to index "untilIndex"
00202    // this means, .. write out the best "untilIndex" individuals.
00203    //
00204 
00205    for ( unsigned int it = 0; it < fGenePool.size(); ++it )
00206       {
00207          Int_t n=0;
00208          if (untilIndex >= -1 ) {
00209             if (untilIndex == -1 ) return;
00210             untilIndex--;
00211          }
00212          Log() << "fitness: " << fGenePool[it].GetFitness() << "    ";
00213          for (vector< Double_t >::iterator vec = fGenePool[it].GetFactors().begin(); 
00214               vec < fGenePool[it].GetFactors().end(); vec++ ) {
00215             Log() << "f_" << n++ << ": " << (*vec) << "     ";
00216          }
00217          Log() << Endl;
00218       }
00219 }
00220 
00221 //_______________________________________________________________________
00222 void TMVA::GeneticPopulation::Print( ostream & out, Int_t untilIndex )
00223 {
00224    // make a little printout to the stream "out" of the individuals up to index "untilIndex"
00225    // this means, .. write out the best "untilIndex" individuals.
00226    //
00227 
00228    for ( unsigned int it = 0; it < fGenePool.size(); ++it ) {
00229       Int_t n=0;
00230       if (untilIndex >= -1 ) {
00231          if (untilIndex == -1 ) return;
00232          untilIndex--;
00233       }
00234       out << "fitness: " << fGenePool[it].GetFitness() << "    ";
00235       for (vector< Double_t >::iterator vec = fGenePool[it].GetFactors().begin(); 
00236            vec < fGenePool[it].GetFactors().end(); vec++ ) {
00237          out << "f_" << n++ << ": " << (*vec) << "     ";
00238       }
00239       out << endl;
00240    }
00241 }
00242 
00243 //_______________________________________________________________________
00244 TH1F* TMVA::GeneticPopulation::VariableDistribution( Int_t varNumber, Int_t bins, 
00245                                                      Int_t min, Int_t max ) 
00246 {
00247    // give back a histogram with the distribution of the coefficients
00248    // parameters:
00249    //          int bins : number of bins of the histogram
00250    //          int min : histogram minimum 
00251    //          int max : maximum value of the histogram
00252    //
00253 
00254    cout << "FAILED! TMVA::GeneticPopulation::VariableDistribution" << endl;
00255 
00256    std::stringstream histName;
00257    histName.clear();
00258    histName.str("v");
00259    histName << varNumber;
00260    TH1F *hist = new TH1F( histName.str().c_str(),histName.str().c_str(), bins,min,max );
00261 
00262    return hist;
00263 }
00264 
00265 //_______________________________________________________________________
00266 vector<Double_t> TMVA::GeneticPopulation::VariableDistribution( Int_t /*varNumber*/ ) 
00267 {
00268    // gives back all the values of coefficient "varNumber" of the current generation
00269    //
00270 
00271    cout << "FAILED! TMVA::GeneticPopulation::VariableDistribution" << endl;
00272 
00273    vector< Double_t > varDist;
00274 
00275    return varDist;
00276 }
00277 
00278 //_______________________________________________________________________
00279 void TMVA::GeneticPopulation::AddPopulation( GeneticPopulation *strangers )
00280 {
00281    // add another population (strangers) to the one of this GeneticPopulation
00282    for (std::vector<TMVA::GeneticGenes>::iterator it = strangers->fGenePool.begin(); 
00283         it != strangers->fGenePool.end(); it++ ) {
00284       GiveHint( it->GetFactors(), it->GetFitness() );
00285    }
00286 }
00287 
00288 //_______________________________________________________________________
00289 void TMVA::GeneticPopulation::AddPopulation( GeneticPopulation &strangers )
00290 {
00291    // add another population (strangers) to the one of this GeneticPopulation
00292    AddPopulation(&strangers);
00293 }
00294 
00295 //_______________________________________________________________________
00296 void TMVA::GeneticPopulation::TrimPopulation()
00297 {
00298    // trim the population to the predefined size
00299    std::sort(fGenePool.begin(), fGenePool.end());
00300    while ( fGenePool.size() > (unsigned int) fPopulationSizeLimit )
00301       fGenePool.pop_back();
00302 }
00303 
00304 //_______________________________________________________________________
00305 void TMVA::GeneticPopulation::GiveHint( std::vector< Double_t >& hint, Double_t fitness )
00306 {
00307    // add an individual (a set of variables) to the population
00308    // if there is a set of variables which is known to perform good, they can be given as a hint to the population
00309    TMVA::GeneticGenes g(hint);
00310    g.SetFitness(fitness);
00311 
00312    fGenePool.push_back( g );
00313 }
00314 
00315 //_______________________________________________________________________
00316 void TMVA::GeneticPopulation::Sort()
00317 {
00318    // sort the genepool according to the fitness of the individuals
00319    std::sort(fGenePool.begin(), fGenePool.end());
00320 }
00321 

Generated on Tue Jul 5 15:25:00 2011 for ROOT_528-00b_version by  doxygen 1.5.1