00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
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
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
00054
00055
00056
00057
00058 fRandomGenerator = new TRandom3( 100 );
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
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
00094 fRandomGenerator->SetSeed( seed );
00095 }
00096
00097
00098 void TMVA::GeneticPopulation::MakeCopies( int number )
00099 {
00100
00101
00102
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
00117
00118
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
00137
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
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170 vector< Double_t>::iterator vec;
00171 vector< TMVA::GeneticRange* >::iterator vecRange;
00172
00173
00174
00175
00176
00177
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
00194
00195 return &(fGenePool[index]);
00196 }
00197
00198
00199 void TMVA::GeneticPopulation::Print( Int_t untilIndex )
00200 {
00201
00202
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
00225
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
00248
00249
00250
00251
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 )
00267 {
00268
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
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
00292 AddPopulation(&strangers);
00293 }
00294
00295
00296 void TMVA::GeneticPopulation::TrimPopulation()
00297 {
00298
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
00308
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
00319 std::sort(fGenePool.begin(), fGenePool.end());
00320 }
00321