Logo ROOT  
Reference Guide
MethodBoost.cxx
Go to the documentation of this file.
1// @(#)root/tmva $Id$
2// Author: Andreas Hoecker, Joerg Stelzer, Helge Voss, Kai Voss,Or Cohen, Jan Therhaag, Eckhard von Toerne
3
4/**********************************************************************************
5 * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6 * Package: TMVA *
7 * Class : MethodCompositeBase *
8 * Web : http://tmva.sourceforge.net *
9 * *
10 * Description: *
11 * Virtual base class for all MVA method *
12 * *
13 * Authors (alphabetical): *
14 * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15 * Peter Speckmayer <Peter.Speckmazer@cern.ch> - CERN, Switzerland *
16 * Joerg Stelzer <Joerg.Stelzer@cern.ch> - CERN, Switzerland *
17 * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
18 * Jan Therhaag <Jan.Therhaag@cern.ch> - U of Bonn, Germany *
19 * Eckhard v. Toerne <evt@uni-bonn.de> - U of Bonn, Germany *
20 * *
21 * Copyright (c) 2005-2011: *
22 * CERN, Switzerland *
23 * U. of Victoria, Canada *
24 * MPI-K Heidelberg, Germany *
25 * U. of Bonn, Germany *
26 * *
27 * Redistribution and use in source and binary forms, with or without *
28 * modification, are permitted according to the terms listed in LICENSE *
29 * (http://tmva.sourceforge.net/LICENSE) *
30 **********************************************************************************/
31
32/*! \class TMVA::MethodBoost
33\ingroup TMVA
34
35Class for boosting a TMVA method
36
37This class is meant to boost a single classifier. Boosting means
38training the classifier a few times. Every time the weights of the
39events are modified according to how well the classifier performed
40on the test sample.
41
42*/
43
44#include "TMVA/MethodBoost.h"
45
47#include "TMVA/Config.h"
48#include "TMVA/Configurable.h"
49#include "TMVA/DataSet.h"
50#include "TMVA/DataSetInfo.h"
51#include "TMVA/IMethod.h"
52#include "TMVA/MethodBase.h"
53#include "TMVA/MethodCategory.h"
55#include "TMVA/MethodDT.h"
56#include "TMVA/MethodFisher.h"
57#include "TMVA/PDF.h"
58#include "TMVA/Results.h"
59#include "TMVA/Timer.h"
60#include "TMVA/Tools.h"
61#include "TMVA/Types.h"
62
63#include "TMVA/SeparationBase.h"
65#include "TMVA/GiniIndex.h"
66#include "TMVA/CrossEntropy.h"
69
70#include "Riostream.h"
71#include "TRandom3.h"
72#include "TFile.h"
73#include "TMath.h"
74#include "TH1F.h"
75#include "TH2F.h"
76#include "TGraph.h"
77#include "TSpline.h"
78#include "TDirectory.h"
79#include "TTree.h"
80
81#include <algorithm>
82#include <iomanip>
83#include <vector>
84#include <cmath>
85
86
87REGISTER_METHOD(Boost)
88
90
91////////////////////////////////////////////////////////////////////////////////
92
94 const TString& methodTitle,
95 DataSetInfo& theData,
96 const TString& theOption ) :
97 TMVA::MethodCompositeBase( jobName, Types::kBoost, methodTitle, theData, theOption)
98 , fBoostNum(0)
99 , fDetailedMonitoring(kFALSE)
100 , fAdaBoostBeta(0)
101 , fRandomSeed(0)
102 , fBaggedSampleFraction(0)
103 , fBoostedMethodTitle(methodTitle)
104 , fBoostedMethodOptions(theOption)
105 , fMonitorBoostedMethod(kFALSE)
106 , fMonitorTree(0)
107 , fBoostWeight(0)
108 , fMethodError(0)
109 , fROC_training(0.0)
110 , fOverlap_integral(0.0)
111 , fMVAvalues(0)
112{
113 fMVAvalues = new std::vector<Float_t>;
114 fDataSetManager = NULL;
116}
117
118////////////////////////////////////////////////////////////////////////////////
119
121 const TString& theWeightFile)
122 : TMVA::MethodCompositeBase( Types::kBoost, dsi, theWeightFile)
123 , fBoostNum(0)
124 , fDetailedMonitoring(kFALSE)
125 , fAdaBoostBeta(0)
126 , fRandomSeed(0)
127 , fBaggedSampleFraction(0)
128 , fBoostedMethodTitle("")
129 , fBoostedMethodOptions("")
130 , fMonitorBoostedMethod(kFALSE)
131 , fMonitorTree(0)
132 , fBoostWeight(0)
133 , fMethodError(0)
134 , fROC_training(0.0)
135 , fOverlap_integral(0.0)
136 , fMVAvalues(0)
137{
138 fMVAvalues = new std::vector<Float_t>;
139 fDataSetManager = NULL;
141}
142
143////////////////////////////////////////////////////////////////////////////////
144/// destructor
145
147{
148 fMethodWeight.clear();
149
150 // the histogram themselves are deleted when the file is closed
151
152 fTrainSigMVAHist.clear();
153 fTrainBgdMVAHist.clear();
154 fBTrainSigMVAHist.clear();
155 fBTrainBgdMVAHist.clear();
156 fTestSigMVAHist.clear();
157 fTestBgdMVAHist.clear();
158
159 if (fMVAvalues) {
160 delete fMVAvalues;
161 fMVAvalues = 0;
162 }
163}
164
165
166////////////////////////////////////////////////////////////////////////////////
167/// Boost can handle classification with 2 classes and regression with one regression-target
168
170{
171 if (type == Types::kClassification && numberClasses == 2) return kTRUE;
172 // if (type == Types::kRegression && numberTargets == 1) return kTRUE;
173 return kFALSE;
174}
175
176
177////////////////////////////////////////////////////////////////////////////////
178
180{
181 DeclareOptionRef( fBoostNum = 1, "Boost_Num",
182 "Number of times the classifier is boosted" );
183
184 DeclareOptionRef( fMonitorBoostedMethod = kTRUE, "Boost_MonitorMethod",
185 "Write monitoring histograms for each boosted classifier" );
186
187 DeclareOptionRef( fDetailedMonitoring = kFALSE, "Boost_DetailedMonitoring",
188 "Produce histograms for detailed boost monitoring" );
189
190 DeclareOptionRef( fBoostType = "AdaBoost", "Boost_Type", "Boosting type for the classifiers" );
191 AddPreDefVal(TString("RealAdaBoost"));
192 AddPreDefVal(TString("AdaBoost"));
193 AddPreDefVal(TString("Bagging"));
194
195 DeclareOptionRef(fBaggedSampleFraction=.6,"Boost_BaggedSampleFraction","Relative size of bagged event sample to original size of the data sample (used whenever bagging is used)" );
196
197 DeclareOptionRef( fAdaBoostBeta = 1.0, "Boost_AdaBoostBeta",
198 "The ADA boost parameter that sets the effect of every boost step on the events' weights" );
199
200 DeclareOptionRef( fTransformString = "step", "Boost_Transform",
201 "Type of transform applied to every boosted method linear, log, step" );
202 AddPreDefVal(TString("step"));
203 AddPreDefVal(TString("linear"));
204 AddPreDefVal(TString("log"));
205 AddPreDefVal(TString("gauss"));
206
207 DeclareOptionRef( fRandomSeed = 0, "Boost_RandomSeed",
208 "Seed for random number generator used for bagging" );
209
210 TMVA::MethodCompositeBase::fMethods.reserve(fBoostNum);
211}
212
213////////////////////////////////////////////////////////////////////////////////
214/// options that are used ONLY for the READER to ensure backward compatibility
215/// they are hence without any effect (the reader is only reading the training
216/// options that HAD been used at the training of the .xml weight file at hand
217
219{
220
222
223 DeclareOptionRef( fHistoricOption = "ByError", "Boost_MethodWeightType",
224 "How to set the final weight of the boosted classifiers" );
225 AddPreDefVal(TString("ByError"));
226 AddPreDefVal(TString("Average"));
227 AddPreDefVal(TString("ByROC"));
228 AddPreDefVal(TString("ByOverlap"));
229 AddPreDefVal(TString("LastMethod"));
230
231 DeclareOptionRef( fHistoricOption = "step", "Boost_Transform",
232 "Type of transform applied to every boosted method linear, log, step" );
233 AddPreDefVal(TString("step"));
234 AddPreDefVal(TString("linear"));
235 AddPreDefVal(TString("log"));
236 AddPreDefVal(TString("gauss"));
237
238 // this option here
239 //DeclareOptionRef( fBoostType = "AdaBoost", "Boost_Type", "Boosting type for the classifiers" );
240 // still exists, but these two possible values
241 AddPreDefVal(TString("HighEdgeGauss"));
242 AddPreDefVal(TString("HighEdgeCoPara"));
243 // have been deleted .. hope that works :)
244
245 DeclareOptionRef( fHistoricBoolOption, "Boost_RecalculateMVACut",
246 "Recalculate the classifier MVA Signallike cut at every boost iteration" );
247
248}
249
250////////////////////////////////////////////////////////////////////////////////
251/// just registering the string from which the boosted classifier will be created
252
254{
255 fBoostedMethodName = Types::Instance().GetMethodName( theMethod );
256 fBoostedMethodTitle = methodTitle;
257 fBoostedMethodOptions = theOption;
258 TString opts=theOption;
259 opts.ToLower();
260 // if (opts.Contains("vartransform")) Log() << kFATAL << "It is not possible to use boost in conjunction with variable transform. Please remove either Boost_Num or VarTransform from the option string"<< methodTitle<<Endl;
261
262 return kTRUE;
263}
264
265////////////////////////////////////////////////////////////////////////////////
266
268{
269}
270
271////////////////////////////////////////////////////////////////////////////////
272/// initialisation routine
273
275{
276
277 Results* results = Data()->GetResults(GetMethodName(), Types::kTraining, GetAnalysisType());
278
279 results->Store(new TH1F("MethodWeight","Normalized Classifier Weight",fBoostNum,0,fBoostNum),"ClassifierWeight");
280 results->Store(new TH1F("BoostWeight","Boost Weight",fBoostNum,0,fBoostNum),"BoostWeight");
281 results->Store(new TH1F("ErrFraction","Error Fraction (by boosted event weights)",fBoostNum,0,fBoostNum),"ErrorFraction");
282 if (fDetailedMonitoring){
283 results->Store(new TH1F("ROCIntegral_test","ROC integral of single classifier (testing sample)",fBoostNum,0,fBoostNum),"ROCIntegral_test");
284 results->Store(new TH1F("ROCIntegralBoosted_test","ROC integral of boosted method (testing sample)",fBoostNum,0,fBoostNum),"ROCIntegralBoosted_test");
285 results->Store(new TH1F("ROCIntegral_train","ROC integral of single classifier (training sample)",fBoostNum,0,fBoostNum),"ROCIntegral_train");
286 results->Store(new TH1F("ROCIntegralBoosted_train","ROC integral of boosted method (training sample)",fBoostNum,0,fBoostNum),"ROCIntegralBoosted_train");
287 results->Store(new TH1F("OverlapIntegal_train","Overlap integral (training sample)",fBoostNum,0,fBoostNum),"Overlap");
288 }
289
290
291 results->GetHist("ClassifierWeight")->GetXaxis()->SetTitle("Index of boosted classifier");
292 results->GetHist("ClassifierWeight")->GetYaxis()->SetTitle("Classifier Weight");
293 results->GetHist("BoostWeight")->GetXaxis()->SetTitle("Index of boosted classifier");
294 results->GetHist("BoostWeight")->GetYaxis()->SetTitle("Boost Weight");
295 results->GetHist("ErrorFraction")->GetXaxis()->SetTitle("Index of boosted classifier");
296 results->GetHist("ErrorFraction")->GetYaxis()->SetTitle("Error Fraction");
297 if (fDetailedMonitoring){
298 results->GetHist("ROCIntegral_test")->GetXaxis()->SetTitle("Index of boosted classifier");
299 results->GetHist("ROCIntegral_test")->GetYaxis()->SetTitle("ROC integral of single classifier");
300 results->GetHist("ROCIntegralBoosted_test")->GetXaxis()->SetTitle("Number of boosts");
301 results->GetHist("ROCIntegralBoosted_test")->GetYaxis()->SetTitle("ROC integral boosted");
302 results->GetHist("ROCIntegral_train")->GetXaxis()->SetTitle("Index of boosted classifier");
303 results->GetHist("ROCIntegral_train")->GetYaxis()->SetTitle("ROC integral of single classifier");
304 results->GetHist("ROCIntegralBoosted_train")->GetXaxis()->SetTitle("Number of boosts");
305 results->GetHist("ROCIntegralBoosted_train")->GetYaxis()->SetTitle("ROC integral boosted");
306 results->GetHist("Overlap")->GetXaxis()->SetTitle("Index of boosted classifier");
307 results->GetHist("Overlap")->GetYaxis()->SetTitle("Overlap integral");
308 }
309
310 results->Store(new TH1F("SoverBtotal","S/B in reweighted training sample",fBoostNum,0,fBoostNum),"SoverBtotal");
311 results->GetHist("SoverBtotal")->GetYaxis()->SetTitle("S/B (boosted sample)");
312 results->GetHist("SoverBtotal")->GetXaxis()->SetTitle("Index of boosted classifier");
313
314 results->Store(new TH1F("SeparationGain","SeparationGain",fBoostNum,0,fBoostNum),"SeparationGain");
315 results->GetHist("SeparationGain")->GetYaxis()->SetTitle("SeparationGain");
316 results->GetHist("SeparationGain")->GetXaxis()->SetTitle("Index of boosted classifier");
317
318
319
320 fMonitorTree= new TTree("MonitorBoost","Boost variables");
321 fMonitorTree->Branch("iMethod",&fCurrentMethodIdx,"iMethod/I");
322 fMonitorTree->Branch("boostWeight",&fBoostWeight,"boostWeight/D");
323 fMonitorTree->Branch("errorFraction",&fMethodError,"errorFraction/D");
324 fMonitorBoostedMethod = kTRUE;
325
326}
327
328
329////////////////////////////////////////////////////////////////////////////////
330
332{
333 Log() << kDEBUG << "CheckSetup: fBoostType="<<fBoostType << Endl;
334 Log() << kDEBUG << "CheckSetup: fAdaBoostBeta="<<fAdaBoostBeta<<Endl;
335 Log() << kDEBUG << "CheckSetup: fBoostWeight="<<fBoostWeight<<Endl;
336 Log() << kDEBUG << "CheckSetup: fMethodError="<<fMethodError<<Endl;
337 Log() << kDEBUG << "CheckSetup: fBoostNum="<<fBoostNum << Endl;
338 Log() << kDEBUG << "CheckSetup: fRandomSeed=" << fRandomSeed<< Endl;
339 Log() << kDEBUG << "CheckSetup: fTrainSigMVAHist.size()="<<fTrainSigMVAHist.size()<<Endl;
340 Log() << kDEBUG << "CheckSetup: fTestSigMVAHist.size()="<<fTestSigMVAHist.size()<<Endl;
341 Log() << kDEBUG << "CheckSetup: fMonitorBoostedMethod=" << (fMonitorBoostedMethod? "true" : "false") << Endl;
342 Log() << kDEBUG << "CheckSetup: MName=" << fBoostedMethodName << " Title="<< fBoostedMethodTitle<< Endl;
343 Log() << kDEBUG << "CheckSetup: MOptions="<< fBoostedMethodOptions << Endl;
344 Log() << kDEBUG << "CheckSetup: fMonitorTree=" << fMonitorTree <<Endl;
345 Log() << kDEBUG << "CheckSetup: fCurrentMethodIdx=" <<fCurrentMethodIdx << Endl;
346 if (fMethods.size()>0) Log() << kDEBUG << "CheckSetup: fMethods[0]" <<fMethods[0]<<Endl;
347 Log() << kDEBUG << "CheckSetup: fMethodWeight.size()" << fMethodWeight.size() << Endl;
348 if (fMethodWeight.size()>0) Log() << kDEBUG << "CheckSetup: fMethodWeight[0]="<<fMethodWeight[0]<<Endl;
349 Log() << kDEBUG << "CheckSetup: trying to repair things" << Endl;
350
351}
352////////////////////////////////////////////////////////////////////////////////
353
355{
356 TDirectory* methodDir( 0 );
357 TString dirName,dirTitle;
358 Int_t StopCounter=0;
359 Results* results = Data()->GetResults(GetMethodName(), Types::kTraining, GetAnalysisType());
360
361
362 InitHistos();
363
364 if (Data()->GetNTrainingEvents()==0) Log() << kFATAL << "<Train> Data() has zero events" << Endl;
365 Data()->SetCurrentType(Types::kTraining);
366
367 if (fMethods.size() > 0) fMethods.clear();
368 fMVAvalues->resize(Data()->GetNTrainingEvents(), 0.0);
369
370 Log() << kINFO << "Training "<< fBoostNum << " " << fBoostedMethodName << " with title " << fBoostedMethodTitle << " Classifiers ... patience please" << Endl;
371 Timer timer( fBoostNum, GetName() );
372
373 ResetBoostWeights();
374
375 // clean boosted method options
376 CleanBoostOptions();
377
378
379 // remove transformations for individual boosting steps
380 // the transformation of the main method will be rerouted to each of the boost steps
381 Ssiz_t varTrafoStart=fBoostedMethodOptions.Index("~VarTransform=");
382 if (varTrafoStart >0) {
383 Ssiz_t varTrafoEnd =fBoostedMethodOptions.Index(":",varTrafoStart);
384 if (varTrafoEnd<varTrafoStart)
385 varTrafoEnd=fBoostedMethodOptions.Length();
386 fBoostedMethodOptions.Remove(varTrafoStart,varTrafoEnd-varTrafoStart);
387 }
388
389 //
390 // training and boosting the classifiers
391 for (fCurrentMethodIdx=0;fCurrentMethodIdx<fBoostNum;fCurrentMethodIdx++) {
392 // the first classifier shows the option string output, the rest not
393 if (fCurrentMethodIdx>0) TMVA::MsgLogger::InhibitOutput();
394
396 fBoostedMethodName.Data(), GetJobName(), Form("%s_B%04i", fBoostedMethodTitle.Data(), fCurrentMethodIdx),
397 DataInfo(), fBoostedMethodOptions);
399
400 // suppressing the rest of the classifier output the right way
401 fCurrentMethod = (dynamic_cast<MethodBase*>(method));
402
403 if (fCurrentMethod==0) {
404 Log() << kFATAL << "uups.. guess the booking of the " << fCurrentMethodIdx << "-th classifier somehow failed" << Endl;
405 return; // hope that makes coverity happy (as if fears I might use the pointer later on, not knowing that FATAL exits
406 }
407
408 // set fDataSetManager if MethodCategory (to enable Category to create datasetinfo objects) // DSMTEST
409 if (fCurrentMethod->GetMethodType() == Types::kCategory) { // DSMTEST
410 MethodCategory *methCat = (dynamic_cast<MethodCategory*>(fCurrentMethod)); // DSMTEST
411 if (!methCat) // DSMTEST
412 Log() << kFATAL << "Method with type kCategory cannot be casted to MethodCategory. /MethodBoost" << Endl; // DSMTEST
413 methCat->fDataSetManager = fDataSetManager; // DSMTEST
414 } // DSMTEST
415
416 fCurrentMethod->SetMsgType(kWARNING);
417 fCurrentMethod->SetupMethod();
418 fCurrentMethod->ParseOptions();
419 // put SetAnalysisType here for the needs of MLP
420 fCurrentMethod->SetAnalysisType( GetAnalysisType() );
421 fCurrentMethod->ProcessSetup();
422 fCurrentMethod->CheckSetup();
423
424
425 // reroute transformationhandler
426 fCurrentMethod->RerouteTransformationHandler (&(this->GetTransformationHandler()));
427
428
429 // creating the directory of the classifier
430 if(!IsSilentFile())
431 {
432 if (fMonitorBoostedMethod) {
433 methodDir=GetFile()->GetDirectory(dirName=Form("%s_B%04i",fBoostedMethodName.Data(),fCurrentMethodIdx));
434 if (methodDir==0) {
435 methodDir=BaseDir()->mkdir(dirName,dirTitle=Form("Directory Boosted %s #%04i", fBoostedMethodName.Data(),fCurrentMethodIdx));
436 }
437 fCurrentMethod->SetMethodDir(methodDir);
438 fCurrentMethod->BaseDir()->cd();
439 }
440 }
441
442 // training
443 TMVA::MethodCompositeBase::fMethods.push_back(method);
444 timer.DrawProgressBar( fCurrentMethodIdx );
445 if (fCurrentMethodIdx==0) MonitorBoost(Types::kBoostProcBegin,fCurrentMethodIdx);
446 MonitorBoost(Types::kBeforeTraining,fCurrentMethodIdx);
447 TMVA::MsgLogger::InhibitOutput(); //suppressing Logger outside the method
448 if (fBoostType=="Bagging") Bagging(); // you want also to train the first classifier on a bagged sample
449 SingleTrain();
451 if(!IsSilentFile())fCurrentMethod->WriteMonitoringHistosToFile();
452
453 // calculate MVA values of current method for all events in training sample
454 // (used later on to get 'misclassified events' etc for the boosting
455 CalcMVAValues();
456
457 if(!IsSilentFile()) if (fCurrentMethodIdx==0 && fMonitorBoostedMethod) CreateMVAHistorgrams();
458
459 // get ROC integral and overlap integral for single method on
460 // training sample if fMethodWeightType == "ByROC" or the user
461 // wants detailed monitoring
462
463 // boosting (reweight training sample)
464 MonitorBoost(Types::kBeforeBoosting,fCurrentMethodIdx);
465 SingleBoost(fCurrentMethod);
466
467 MonitorBoost(Types::kAfterBoosting,fCurrentMethodIdx);
468 results->GetHist("BoostWeight")->SetBinContent(fCurrentMethodIdx+1,fBoostWeight);
469 results->GetHist("ErrorFraction")->SetBinContent(fCurrentMethodIdx+1,fMethodError);
470
471 if (fDetailedMonitoring) {
472 fROC_training = GetBoostROCIntegral(kTRUE, Types::kTraining, kTRUE);
473 results->GetHist("ROCIntegral_test")->SetBinContent(fCurrentMethodIdx+1, GetBoostROCIntegral(kTRUE, Types::kTesting));
474 results->GetHist("ROCIntegralBoosted_test")->SetBinContent(fCurrentMethodIdx+1, GetBoostROCIntegral(kFALSE, Types::kTesting));
475 results->GetHist("ROCIntegral_train")->SetBinContent(fCurrentMethodIdx+1, fROC_training);
476 results->GetHist("ROCIntegralBoosted_train")->SetBinContent(fCurrentMethodIdx+1, GetBoostROCIntegral(kFALSE, Types::kTraining));
477 results->GetHist("Overlap")->SetBinContent(fCurrentMethodIdx+1, fOverlap_integral);
478 }
479
480
481
482 fMonitorTree->Fill();
483
484 // stop boosting if needed when error has reached 0.5
485 // thought of counting a few steps, but it doesn't seem to be necessary
486 Log() << kDEBUG << "AdaBoost (methodErr) err = " << fMethodError << Endl;
487 if (fMethodError > 0.49999) StopCounter++;
488 if (StopCounter > 0 && fBoostType != "Bagging") {
489 timer.DrawProgressBar( fBoostNum );
490 fBoostNum = fCurrentMethodIdx+1;
491 Log() << kINFO << "Error rate has reached 0.5 ("<< fMethodError<<"), boosting process stopped at #" << fBoostNum << " classifier" << Endl;
492 if (fBoostNum < 5)
493 Log() << kINFO << "The classifier might be too strong to boost with Beta = " << fAdaBoostBeta << ", try reducing it." <<Endl;
494 break;
495 }
496 }
497
498 //as MethodBoost acts not on a private event sample (like MethodBDT does), we need to remember not
499 // to leave "boosted" events to the next classifier in the factory
500
501 ResetBoostWeights();
502
503 Timer* timer1= new Timer( fBoostNum, GetName() );
504 // normalizing the weights of the classifiers
505 for (fCurrentMethodIdx=0;fCurrentMethodIdx<fBoostNum;fCurrentMethodIdx++) {
506 // performing post-boosting actions
507
508 timer1->DrawProgressBar( fCurrentMethodIdx );
509
510 if (fCurrentMethodIdx==fBoostNum) {
511 Log() << kINFO << "Elapsed time: " << timer1->GetElapsedTime()
512 << " " << Endl;
513 }
514
515 TH1F* tmp = dynamic_cast<TH1F*>( results->GetHist("ClassifierWeight") );
516 if (tmp) tmp->SetBinContent(fCurrentMethodIdx+1,fMethodWeight[fCurrentMethodIdx]);
517
518 }
519
520 // Ensure that in case of only 1 boost the method weight equals
521 // 1.0. This avoids unexpected behaviour in case of very bad
522 // classifiers which have fBoostWeight=1 or fMethodError=0.5,
523 // because their weight would be set to zero. This behaviour is
524 // not ok if one boosts just one time.
525 if (fMethods.size()==1) fMethodWeight[0] = 1.0;
526
527 MonitorBoost(Types::kBoostProcEnd);
528
529 delete timer1;
530}
531
532////////////////////////////////////////////////////////////////////////////////
533
535{
536 fBoostedMethodOptions=GetOptions();
537}
538
539////////////////////////////////////////////////////////////////////////////////
540
542{
543 if (fBoostNum <=0) Log() << kFATAL << "CreateHistograms called before fBoostNum is initialized" << Endl;
544 // calculating histograms boundaries and creating histograms..
545 // nrms = number of rms around the average to use for outline (of the 0 classifier)
546 Double_t meanS, meanB, rmsS, rmsB, xmin, xmax, nrms = 10;
547 Int_t signalClass = 0;
548 if (DataInfo().GetClassInfo("Signal") != 0) {
549 signalClass = DataInfo().GetClassInfo("Signal")->GetNumber();
550 }
551 gTools().ComputeStat( GetEventCollection( Types::kMaxTreeType ), fMVAvalues,
552 meanS, meanB, rmsS, rmsB, xmin, xmax, signalClass );
553
555 xmin = TMath::Max( TMath::Min(meanS - nrms*rmsS, meanB - nrms*rmsB ), xmin );
556 xmax = TMath::Min( TMath::Max(meanS + nrms*rmsS, meanB + nrms*rmsB ), xmax ) + 0.00001;
557
558 // creating all the histograms
559 for (UInt_t imtd=0; imtd<fBoostNum; imtd++) {
560 fTrainSigMVAHist .push_back( new TH1F( Form("MVA_Train_S_%04i",imtd), "MVA_Train_S", fNbins, xmin, xmax ) );
561 fTrainBgdMVAHist .push_back( new TH1F( Form("MVA_Train_B%04i", imtd), "MVA_Train_B", fNbins, xmin, xmax ) );
562 fBTrainSigMVAHist.push_back( new TH1F( Form("MVA_BTrain_S%04i",imtd), "MVA_BoostedTrain_S", fNbins, xmin, xmax ) );
563 fBTrainBgdMVAHist.push_back( new TH1F( Form("MVA_BTrain_B%04i",imtd), "MVA_BoostedTrain_B", fNbins, xmin, xmax ) );
564 fTestSigMVAHist .push_back( new TH1F( Form("MVA_Test_S%04i", imtd), "MVA_Test_S", fNbins, xmin, xmax ) );
565 fTestBgdMVAHist .push_back( new TH1F( Form("MVA_Test_B%04i", imtd), "MVA_Test_B", fNbins, xmin, xmax ) );
566 }
567}
568
569////////////////////////////////////////////////////////////////////////////////
570/// resetting back the boosted weights of the events to 1
571
573{
574 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
575 const Event *ev = Data()->GetEvent(ievt);
576 ev->SetBoostWeight( 1.0 );
577 }
578}
579
580////////////////////////////////////////////////////////////////////////////////
581
583{
584 TDirectory* dir=0;
585 if (fMonitorBoostedMethod) {
586 for (UInt_t imtd=0;imtd<fBoostNum;imtd++) {
587
588 //writing the histograms in the specific classifier's directory
589 MethodBase* m = dynamic_cast<MethodBase*>(fMethods[imtd]);
590 if (!m) continue;
591 dir = m->BaseDir();
592 dir->cd();
593 fTrainSigMVAHist[imtd]->SetDirectory(dir);
594 fTrainSigMVAHist[imtd]->Write();
595 fTrainBgdMVAHist[imtd]->SetDirectory(dir);
596 fTrainBgdMVAHist[imtd]->Write();
597 fBTrainSigMVAHist[imtd]->SetDirectory(dir);
598 fBTrainSigMVAHist[imtd]->Write();
599 fBTrainBgdMVAHist[imtd]->SetDirectory(dir);
600 fBTrainBgdMVAHist[imtd]->Write();
601 }
602 }
603
604 // going back to the original folder
605 BaseDir()->cd();
606
607 fMonitorTree->Write();
608}
609
610////////////////////////////////////////////////////////////////////////////////
611
613{
615 if (fMonitorBoostedMethod) {
616 UInt_t nloop = fTestSigMVAHist.size();
617 if (fMethods.size()<nloop) nloop = fMethods.size();
618 //running over all the events and populating the test MVA histograms
619 Data()->SetCurrentType(Types::kTesting);
620 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
621 const Event* ev = GetEvent(ievt);
622 Float_t w = ev->GetWeight();
623 if (DataInfo().IsSignal(ev)) {
624 for (UInt_t imtd=0; imtd<nloop; imtd++) {
625 fTestSigMVAHist[imtd]->Fill(fMethods[imtd]->GetMvaValue(),w);
626 }
627 }
628 else {
629 for (UInt_t imtd=0; imtd<nloop; imtd++) {
630 fTestBgdMVAHist[imtd]->Fill(fMethods[imtd]->GetMvaValue(),w);
631 }
632 }
633 }
634 Data()->SetCurrentType(Types::kTraining);
635 }
636}
637
638////////////////////////////////////////////////////////////////////////////////
639
641{
643 if (treetype==Types::kTraining) return;
644 UInt_t nloop = fTestSigMVAHist.size();
645 if (fMethods.size()<nloop) nloop = fMethods.size();
646 if (fMonitorBoostedMethod) {
647 TDirectory* dir=0;
648 for (UInt_t imtd=0;imtd<nloop;imtd++) {
649 //writing the histograms in the specific classifier's directory
650 MethodBase* mva = dynamic_cast<MethodBase*>(fMethods[imtd]);
651 if (!mva) continue;
652 dir = mva->BaseDir();
653 if (dir==0) continue;
654 dir->cd();
655 fTestSigMVAHist[imtd]->SetDirectory(dir);
656 fTestSigMVAHist[imtd]->Write();
657 fTestBgdMVAHist[imtd]->SetDirectory(dir);
658 fTestBgdMVAHist[imtd]->Write();
659 }
660 }
661}
662
663////////////////////////////////////////////////////////////////////////////////
664/// process user options
665
667{
668}
669
670////////////////////////////////////////////////////////////////////////////////
671/// initialization
672
674{
675 Data()->SetCurrentType(Types::kTraining);
676 MethodBase* meth = dynamic_cast<MethodBase*>(GetLastMethod());
677 if (meth){
678 meth->SetSilentFile(IsSilentFile());
679 if(IsModelPersistence()){
680 TString _fFileDir= DataInfo().GetName();
681 _fFileDir+="/"+gConfig().GetIONames().fWeightFileDir;
682 meth->SetWeightFileDir(_fFileDir);
683 }
684 meth->SetModelPersistence(IsModelPersistence());
685 meth->TrainMethod();
686 }
687}
688
689////////////////////////////////////////////////////////////////////////////////
690/// find the CUT on the individual MVA that defines an event as
691/// correct or misclassified (to be used in the boosting process)
692
694{
695 if (!method || method->GetMethodType() == Types::kDT ){ return;}
696
697 // creating a fine histograms containing the error rate
698 const Int_t nBins=10001;
699 Double_t minMVA=150000;
700 Double_t maxMVA=-150000;
701 for (Long64_t ievt=0; ievt<Data()->GetNEvents(); ievt++) {
702 GetEvent(ievt);
703 Double_t val=method->GetMvaValue();
704 //Helge .. I think one could very well use fMVAValues for that ... -->to do
705 if (val>maxMVA) maxMVA=val;
706 if (val<minMVA) minMVA=val;
707 }
708 maxMVA = maxMVA+(maxMVA-minMVA)/nBins;
709
710 Double_t sum = 0.;
711
712 TH1D *mvaS = new TH1D(Form("MVAS_%d",fCurrentMethodIdx) ,"",nBins,minMVA,maxMVA);
713 TH1D *mvaB = new TH1D(Form("MVAB_%d",fCurrentMethodIdx) ,"",nBins,minMVA,maxMVA);
714 TH1D *mvaSC = new TH1D(Form("MVASC_%d",fCurrentMethodIdx),"",nBins,minMVA,maxMVA);
715 TH1D *mvaBC = new TH1D(Form("MVABC_%d",fCurrentMethodIdx),"",nBins,minMVA,maxMVA);
716
717
718 Results* results = Data()->GetResults(GetMethodName(), Types::kTraining, GetAnalysisType());
719 if (fDetailedMonitoring){
720 results->Store(mvaS, Form("MVAS_%d",fCurrentMethodIdx));
721 results->Store(mvaB, Form("MVAB_%d",fCurrentMethodIdx));
722 results->Store(mvaSC,Form("MVASC_%d",fCurrentMethodIdx));
723 results->Store(mvaBC,Form("MVABC_%d",fCurrentMethodIdx));
724 }
725
726 for (Long64_t ievt=0; ievt<Data()->GetNEvents(); ievt++) {
727
728 Double_t weight = GetEvent(ievt)->GetWeight();
729 Double_t mvaVal=method->GetMvaValue();
730 sum +=weight;
731 if (DataInfo().IsSignal(GetEvent(ievt))){
732 mvaS->Fill(mvaVal,weight);
733 }else {
734 mvaB->Fill(mvaVal,weight);
735 }
736 }
737 SeparationBase *sepGain;
738
739
740 // Boosting should use Misclassification not Gini Index (changed, Helge 31.5.2013)
741 // WARNING! It works with Misclassification only if you fix the signal to
742 // background at every step. Strangely enough, there are better results
743 // ( as seen with BDT ) if you use Gini Index, and accept that sometimes no
744 // sensible cut is found - i.e. the cut is then outside the MVA value range,
745 // all events are classified as background and then according to the Boost
746 // algorithm something is renormed 'automatically' ... so that in the next
747 // step again the result is something sensible.
748 // Strange ... that THIS is supposed to be right?
749
750 // SeparationBase *sepGain2 = new MisClassificationError();
751 //sepGain = new MisClassificationError();
752 sepGain = new GiniIndex();
753 //sepGain = new CrossEntropy();
754
755 Double_t sTot = mvaS->GetSum();
756 Double_t bTot = mvaB->GetSum();
757
758 mvaSC->SetBinContent(1,mvaS->GetBinContent(1));
759 mvaBC->SetBinContent(1,mvaB->GetBinContent(1));
760 Double_t sSel=0;
761 Double_t bSel=0;
762 Double_t separationGain=sepGain->GetSeparationGain(sSel,bSel,sTot,bTot);
763 Double_t mvaCut=mvaSC->GetBinLowEdge(1);
764 Double_t sSelCut=sSel;
765 Double_t bSelCut=bSel;
766 // std::cout << "minMVA =" << minMVA << " maxMVA = " << maxMVA << " width = " << mvaSC->GetBinWidth(1) << std::endl;
767
768 // for (Int_t ibin=1;ibin<=nBins;ibin++) std::cout << " cutvalues[" << ibin<<"]="<<mvaSC->GetBinLowEdge(ibin) << " " << mvaSC->GetBinCenter(ibin) << std::endl;
769 Double_t mvaCutOrientation=1; // 1 if mva > mvaCut --> Signal and -1 if mva < mvaCut (i.e. mva*-1 > mvaCut*-1) --> Signal
770 for (Int_t ibin=1;ibin<=nBins;ibin++){
771 mvaSC->SetBinContent(ibin,mvaS->GetBinContent(ibin)+mvaSC->GetBinContent(ibin-1));
772 mvaBC->SetBinContent(ibin,mvaB->GetBinContent(ibin)+mvaBC->GetBinContent(ibin-1));
773
774 sSel=mvaSC->GetBinContent(ibin);
775 bSel=mvaBC->GetBinContent(ibin);
776
777 // if (ibin==nBins){
778 // std::cout << "Last bin s="<< sSel <<" b="<<bSel << " s="<< sTot-sSel <<" b="<<bTot-bSel << endl;
779 // }
780
781 if (separationGain < sepGain->GetSeparationGain(sSel,bSel,sTot,bTot)
782 // && (mvaSC->GetBinCenter(ibin) >0 || (fCurrentMethodIdx+1)%2 )
783 ){
784 separationGain = sepGain->GetSeparationGain(sSel,bSel,sTot,bTot);
785 // mvaCut=mvaSC->GetBinCenter(ibin);
786 mvaCut=mvaSC->GetBinLowEdge(ibin+1);
787 // if (sSel/bSel > (sTot-sSel)/(bTot-bSel)) mvaCutOrientation=-1;
788 if (sSel*(bTot-bSel) > (sTot-sSel)*bSel) mvaCutOrientation=-1;
789 else mvaCutOrientation=1;
790 sSelCut=sSel;
791 bSelCut=bSel;
792 // std::cout << "new cut at " << mvaCut << "with s="<<sTot-sSel << " b="<<bTot-bSel << std::endl;
793 }
794 /*
795 Double_t ori;
796 if (sSel/bSel > (sTot-sSel)/(bTot-bSel)) ori=-1;
797 else ori=1;
798 std::cout << ibin << " mvacut="<<mvaCut
799 << " sTot=" << sTot
800 << " bTot=" << bTot
801 << " sSel=" << sSel
802 << " bSel=" << bSel
803 << " s/b(1)=" << sSel/bSel
804 << " s/b(2)=" << (sTot-sSel)/(bTot-bSel)
805 << " sepGain="<<sepGain->GetSeparationGain(sSel,bSel,sTot,bTot)
806 << " sepGain2="<<sepGain2->GetSeparationGain(sSel,bSel,sTot,bTot)
807 << " " <<ori
808 << std::endl;
809 */
810
811 }
812
813 if (0){
814 double parentIndex=sepGain->GetSeparationIndex(sTot,bTot);
815 double leftIndex =sepGain->GetSeparationIndex(sSelCut,bSelCut);
816 double rightIndex =sepGain->GetSeparationIndex(sTot-sSelCut,bTot-bSelCut);
817 std::cout
818 << " sTot=" << sTot
819 << " bTot=" << bTot
820 << " s="<<sSelCut
821 << " b="<<bSelCut
822 << " s2="<<(sTot-sSelCut)
823 << " b2="<<(bTot-bSelCut)
824 << " s/b(1)=" << sSelCut/bSelCut
825 << " s/b(2)=" << (sTot-sSelCut)/(bTot-bSelCut)
826 << " index before cut=" << parentIndex
827 << " after: left=" << leftIndex
828 << " after: right=" << rightIndex
829 << " sepGain=" << parentIndex-( (sSelCut+bSelCut) * leftIndex + (sTot-sSelCut+bTot-bSelCut) * rightIndex )/(sTot+bTot)
830 << " sepGain="<<separationGain
831 << " sepGain="<<sepGain->GetSeparationGain(sSelCut,bSelCut,sTot,bTot)
832 << " cut=" << mvaCut
833 << " idx="<<fCurrentMethodIdx
834 << " cutOrientation="<<mvaCutOrientation
835 << std::endl;
836 }
837 method->SetSignalReferenceCut(mvaCut);
838 method->SetSignalReferenceCutOrientation(mvaCutOrientation);
839
840 results->GetHist("SeparationGain")->SetBinContent(fCurrentMethodIdx+1,separationGain);
841
842
843 Log() << kDEBUG << "(old step) Setting method cut to " <<method->GetSignalReferenceCut()<< Endl;
844
845 if(IsSilentFile())
846 {
847 mvaS ->Delete();
848 mvaB ->Delete();
849 mvaSC->Delete();
850 mvaBC->Delete();
851 }
852}
853
854////////////////////////////////////////////////////////////////////////////////
855
857{
858 Double_t returnVal=-1;
859
860
861 if (fBoostType=="AdaBoost") returnVal = this->AdaBoost (method,1);
862 else if (fBoostType=="RealAdaBoost") returnVal = this->AdaBoost (method,0);
863 else if (fBoostType=="Bagging") returnVal = this->Bagging ();
864 else{
865 Log() << kFATAL << "<Boost> unknown boost option " << fBoostType<< " called" << Endl;
866 }
867 fMethodWeight.push_back(returnVal);
868 return returnVal;
869}
870////////////////////////////////////////////////////////////////////////////////
871/// the standard (discrete or real) AdaBoost algorithm
872
874{
875 if (!method) {
876 Log() << kWARNING << " AdaBoost called without classifier reference - needed for calculating AdaBoost " << Endl;
877 return 0;
878 }
879
880 Float_t w,v; Bool_t sig=kTRUE;
881 Double_t sumAll=0, sumWrong=0;
882 Bool_t* WrongDetection=new Bool_t[GetNEvents()];
883 QuickMVAProbEstimator *MVAProb=NULL;
884
885 if (discreteAdaBoost) {
886 FindMVACut(method);
887 Log() << kDEBUG << " individual mva cut value = " << method->GetSignalReferenceCut() << Endl;
888 } else {
889 MVAProb=new TMVA::QuickMVAProbEstimator();
890 // the RealAdaBoost does use a simple "yes (signal)" or "no (background)"
891 // answer from your single MVA, but a "signal probability" instead (in the BDT case,
892 // that would be the 'purity' in the leaf node. For some MLP parameter, the MVA output
893 // can also interpreted as a probability, but here I try a general approach to get this
894 // probability from the MVA distributions...
895
896 for (Long64_t evt=0; evt<GetNEvents(); evt++) {
897 const Event* ev = Data()->GetEvent(evt);
898 MVAProb->AddEvent(fMVAvalues->at(evt),ev->GetWeight(),ev->GetClass());
899 }
900 }
901
902
903 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) WrongDetection[ievt]=kTRUE;
904
905 // finding the wrong events and calculating their total weights
906 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
907 const Event* ev = GetEvent(ievt);
908 sig=DataInfo().IsSignal(ev);
909 v = fMVAvalues->at(ievt);
910 w = ev->GetWeight();
911 sumAll += w;
912 if(!IsSilentFile())
913 {
914 if (fMonitorBoostedMethod) {
915 if (sig) {
916 fBTrainSigMVAHist[fCurrentMethodIdx]->Fill(v,w);
917 fTrainSigMVAHist[fCurrentMethodIdx]->Fill(v,ev->GetOriginalWeight());
918 }
919 else {
920 fBTrainBgdMVAHist[fCurrentMethodIdx]->Fill(v,w);
921 fTrainBgdMVAHist[fCurrentMethodIdx]->Fill(v,ev->GetOriginalWeight());
922 }
923 }
924 }
925
926 if (discreteAdaBoost){
927 if (sig == method->IsSignalLike(fMVAvalues->at(ievt))){
928 WrongDetection[ievt]=kFALSE;
929 }else{
930 WrongDetection[ievt]=kTRUE;
931 sumWrong+=w;
932 }
933 }else{
934 Double_t mvaProb = MVAProb->GetMVAProbAt((Float_t)fMVAvalues->at(ievt));
935 mvaProb = 2*(mvaProb-0.5);
936 Int_t trueType;
937 if (DataInfo().IsSignal(ev)) trueType = 1;
938 else trueType = -1;
939 sumWrong+= w*trueType*mvaProb;
940 }
941 }
942
943 fMethodError=sumWrong/sumAll;
944
945 // calculating the fMethodError and the boostWeight out of it uses the formula
946 // w = ((1-err)/err)^beta
947
948 Double_t boostWeight=0;
949
950 if (fMethodError == 0) { //no misclassification made.. perfect, no boost ;)
951 Log() << kWARNING << "Your classifier worked perfectly on the training sample --> serious overtraining expected and no boosting done " << Endl;
952 }else{
953
954 if (discreteAdaBoost)
955 boostWeight = TMath::Log((1.-fMethodError)/fMethodError)*fAdaBoostBeta;
956 else
957 boostWeight = TMath::Log((1.+fMethodError)/(1-fMethodError))*fAdaBoostBeta;
958
959
960 // std::cout << "boostweight = " << boostWeight << std::endl;
961
962 // ADA boosting, rescaling the weight of the wrong events according to the error level
963 // over the entire test sample rescaling all the weights to have the same sum, but without
964 // touching the original weights (changing only the boosted weight of all the events)
965 // first reweight
966
967 Double_t newSum=0., oldSum=0.;
968
969
970 Double_t boostfactor = TMath::Exp(boostWeight);
971
972
973 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
974 const Event* ev = Data()->GetEvent(ievt);
975 oldSum += ev->GetWeight();
976 if (discreteAdaBoost){
977 // events are classified as Signal OR background .. right or wrong
978 if (WrongDetection[ievt] && boostWeight != 0) {
979 if (ev->GetWeight() > 0) ev->ScaleBoostWeight(boostfactor);
980 else ev->ScaleBoostWeight(1./boostfactor);
981 }
982 // if (ievt<30) std::cout<<ievt<<" var0="<<ev->GetValue(0)<<" var1="<<ev->GetValue(1)<<" weight="<<ev->GetWeight() << " boostby:"<<boostfactor<<std::endl;
983
984 }else{
985 // events are classified by their probability of being signal or background
986 // (eventually you should write this one - i.e. re-use the MVA value that were already
987 // calculated and stored.. however ,for the moment ..
988 Double_t mvaProb = MVAProb->GetMVAProbAt((Float_t)fMVAvalues->at(ievt));
989 mvaProb = 2*(mvaProb-0.5);
990 // mvaProb = (1-mvaProb);
991
992 Int_t trueType=1;
993 if (DataInfo().IsSignal(ev)) trueType = 1;
994 else trueType = -1;
995
996 boostfactor = TMath::Exp(-1*boostWeight*trueType*mvaProb);
997 if (ev->GetWeight() > 0) ev->ScaleBoostWeight(boostfactor);
998 else ev->ScaleBoostWeight(1./boostfactor);
999
1000 }
1001 newSum += ev->GetWeight();
1002 }
1003
1004 Double_t normWeight = oldSum/newSum;
1005 // next normalize the weights
1006 Double_t normSig=0, normBkg=0;
1007 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1008 const Event* ev = Data()->GetEvent(ievt);
1009 ev->ScaleBoostWeight(normWeight);
1010 if (ev->GetClass()) normSig+=ev->GetWeight();
1011 else normBkg+=ev->GetWeight();
1012 }
1013
1014 Results* results = Data()->GetResults(GetMethodName(), Types::kTraining, GetAnalysisType());
1015 results->GetHist("SoverBtotal")->SetBinContent(fCurrentMethodIdx+1, normSig/normBkg);
1016
1017 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1018 const Event* ev = Data()->GetEvent(ievt);
1019
1020 if (ev->GetClass()) ev->ScaleBoostWeight(oldSum/normSig/2);
1021 else ev->ScaleBoostWeight(oldSum/normBkg/2);
1022 }
1023 }
1024
1025 delete[] WrongDetection;
1026 if (MVAProb) delete MVAProb;
1027
1028 fBoostWeight = boostWeight; // used ONLY for the monitoring tree
1029
1030 return boostWeight;
1031}
1032
1033
1034////////////////////////////////////////////////////////////////////////////////
1035/// Bagging or Bootstrap boosting, gives new random poisson weight for every event
1036
1038{
1039 TRandom3 *trandom = new TRandom3(fRandomSeed+fMethods.size());
1040 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1041 const Event* ev = Data()->GetEvent(ievt);
1042 ev->SetBoostWeight(trandom->PoissonD(fBaggedSampleFraction));
1043 }
1044 fBoostWeight = 1; // used ONLY for the monitoring tree
1045 return 1.;
1046}
1047
1048
1049////////////////////////////////////////////////////////////////////////////////
1050/// Get help message text
1051///
1052/// typical length of text line:
1053/// "|--------------------------------------------------------------|"
1054
1056{
1057 Log() << Endl;
1058 Log() << gTools().Color("bold") << "--- Short description:" << gTools().Color("reset") << Endl;
1059 Log() << Endl;
1060 Log() << "This method combines several classifier of one species in a "<<Endl;
1061 Log() << "single multivariate quantity via the boost algorithm." << Endl;
1062 Log() << "the output is a weighted sum over all individual classifiers" <<Endl;
1063 Log() << "By default, the AdaBoost method is employed, which gives " << Endl;
1064 Log() << "events that were misclassified in the previous tree a larger " << Endl;
1065 Log() << "weight in the training of the following classifier."<<Endl;
1066 Log() << "Optionally, Bagged boosting can also be applied." << Endl;
1067 Log() << Endl;
1068 Log() << gTools().Color("bold") << "--- Performance tuning via configuration options:" << gTools().Color("reset") << Endl;
1069 Log() << Endl;
1070 Log() << "The most important parameter in the configuration is the "<<Endl;
1071 Log() << "number of boosts applied (Boost_Num) and the choice of boosting"<<Endl;
1072 Log() << "(Boost_Type), which can be set to either AdaBoost or Bagging." << Endl;
1073 Log() << "AdaBoosting: The most important parameters in this configuration" <<Endl;
1074 Log() << "is the beta parameter (Boost_AdaBoostBeta) " << Endl;
1075 Log() << "When boosting a linear classifier, it is sometimes advantageous"<<Endl;
1076 Log() << "to transform the MVA output non-linearly. The following options" <<Endl;
1077 Log() << "are available: step, log, and minmax, the default is no transform."<<Endl;
1078 Log() <<Endl;
1079 Log() << "Some classifiers are hard to boost and do not improve much in"<<Endl;
1080 Log() << "their performance by boosting them, some even slightly deteriorate"<< Endl;
1081 Log() << "due to the boosting." <<Endl;
1082 Log() << "The booking of the boost method is special since it requires"<<Endl;
1083 Log() << "the booing of the method to be boosted and the boost itself."<<Endl;
1084 Log() << "This is solved by booking the method to be boosted and to add"<<Endl;
1085 Log() << "all Boost parameters, which all begin with \"Boost_\" to the"<<Endl;
1086 Log() << "options string. The factory separates the options and initiates"<<Endl;
1087 Log() << "the boost process. The TMVA macro directory contains the example"<<Endl;
1088 Log() << "macro \"Boost.C\"" <<Endl;
1089}
1090
1091////////////////////////////////////////////////////////////////////////////////
1092
1094{
1095 return 0;
1096}
1097
1098////////////////////////////////////////////////////////////////////////////////
1099/// return boosted MVA response
1100
1102{
1103 Double_t mvaValue = 0;
1104 Double_t norm = 0;
1106 //Double_t fact = TMath::Exp(-1.)+TMath::Exp(1.);
1107 for (UInt_t i=0;i< fMethods.size(); i++){
1108 MethodBase* m = dynamic_cast<MethodBase*>(fMethods[i]);
1109 if (m==0) continue;
1110 Double_t val = fTmpEvent ? m->GetMvaValue(fTmpEvent) : m->GetMvaValue();
1111 Double_t sigcut = m->GetSignalReferenceCut();
1112
1113 // default is no transform
1114 if (fTransformString == "linear"){
1115
1116 }
1117 else if (fTransformString == "log"){
1118 if (val < sigcut) val = sigcut;
1119
1120 val = TMath::Log((val-sigcut)+epsilon);
1121 }
1122 else if (fTransformString == "step" ){
1123 if (m->IsSignalLike(val)) val = 1.;
1124 else val = -1.;
1125 }
1126 else if (fTransformString == "gauss"){
1127 val = TMath::Gaus((val-sigcut),1);
1128 }
1129 else {
1130 Log() << kFATAL << "error unknown transformation " << fTransformString<<Endl;
1131 }
1132 mvaValue+=val*fMethodWeight[i];
1133 norm +=fMethodWeight[i];
1134 // std::cout << "mva("<<i<<") = "<<val<<" " << valx<< " " << mvaValue<<" and sigcut="<<sigcut << std::endl;
1135 }
1136 mvaValue/=norm;
1137 // cannot determine error
1138 NoErrorCalc(err, errUpper);
1139
1140 return mvaValue;
1141}
1142
1143////////////////////////////////////////////////////////////////////////////////
1144/// Calculate the ROC integral of a single classifier or even the
1145/// whole boosted classifier. The tree type (training or testing
1146/// sample) is specified by 'eTT'.
1147///
1148/// If tree type kTraining is set, the original training sample is
1149/// used to compute the ROC integral (original weights).
1150///
1151/// - singleMethod - if kTRUE, return ROC integral of single (last
1152/// trained) classifier; if kFALSE, return ROC
1153/// integral of full classifier
1154///
1155/// - eTT - tree type (Types::kTraining / Types::kTesting)
1156///
1157/// - CalcOverlapIntergral - if kTRUE, the overlap integral of the
1158/// signal/background MVA distributions
1159/// is calculated and stored in
1160/// 'fOverlap_integral'
1161
1163{
1164 // set data sample training / testing
1165 Data()->SetCurrentType(eTT);
1166
1167 MethodBase* method = singleMethod ? dynamic_cast<MethodBase*>(fMethods.back()) : 0; // ToDo CoVerity flags this line as there is no protection against a zero-pointer delivered by dynamic_cast
1168 // to make CoVerity happy (although, OF COURSE, the last method in the committee
1169 // has to be also of type MethodBase as ANY method is... hence the dynamic_cast
1170 // will never by "zero" ...
1171 if (singleMethod && !method) {
1172 Log() << kFATAL << " What do you do? Your method:"
1173 << fMethods.back()->GetName()
1174 << " seems not to be a propper TMVA method"
1175 << Endl;
1176 std::exit(1);
1177 }
1178 Double_t err = 0.0;
1179
1180 // temporary renormalize the method weights in case of evaluation
1181 // of full classifier.
1182 // save the old normalization of the methods
1183 std::vector<Double_t> OldMethodWeight(fMethodWeight);
1184 if (!singleMethod) {
1185 // calculate sum of weights of all methods
1186 Double_t AllMethodsWeight = 0;
1187 for (UInt_t i=0; i<=fCurrentMethodIdx; i++)
1188 AllMethodsWeight += fMethodWeight.at(i);
1189 // normalize the weights of the classifiers
1190 if (AllMethodsWeight != 0.0) {
1191 for (UInt_t i=0; i<=fCurrentMethodIdx; i++)
1192 fMethodWeight[i] /= AllMethodsWeight;
1193 }
1194 }
1195
1196 // calculate MVA values
1197 Double_t meanS, meanB, rmsS, rmsB, xmin, xmax, nrms = 10;
1198 std::vector <Float_t>* mvaRes;
1199 if (singleMethod && eTT==Types::kTraining)
1200 mvaRes = fMVAvalues; // values already calculated
1201 else {
1202 mvaRes = new std::vector <Float_t>(GetNEvents());
1203 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1204 GetEvent(ievt);
1205 (*mvaRes)[ievt] = singleMethod ? method->GetMvaValue(&err) : GetMvaValue(&err);
1206 }
1207 }
1208
1209 // restore the method weights
1210 if (!singleMethod)
1211 fMethodWeight = OldMethodWeight;
1212
1213 // now create histograms for calculation of the ROC integral
1214 Int_t signalClass = 0;
1215 if (DataInfo().GetClassInfo("Signal") != 0) {
1216 signalClass = DataInfo().GetClassInfo("Signal")->GetNumber();
1217 }
1218 gTools().ComputeStat( GetEventCollection(eTT), mvaRes,
1219 meanS, meanB, rmsS, rmsB, xmin, xmax, signalClass );
1220
1222 xmin = TMath::Max( TMath::Min(meanS - nrms*rmsS, meanB - nrms*rmsB ), xmin );
1223 xmax = TMath::Min( TMath::Max(meanS + nrms*rmsS, meanB + nrms*rmsB ), xmax ) + 0.0001;
1224
1225 // calculate ROC integral
1226 TH1* mva_s = new TH1F( "MVA_S", "MVA_S", fNbins, xmin, xmax );
1227 TH1* mva_b = new TH1F( "MVA_B", "MVA_B", fNbins, xmin, xmax );
1228 TH1 *mva_s_overlap=0, *mva_b_overlap=0;
1229 if (CalcOverlapIntergral) {
1230 mva_s_overlap = new TH1F( "MVA_S_OVERLAP", "MVA_S_OVERLAP", fNbins, xmin, xmax );
1231 mva_b_overlap = new TH1F( "MVA_B_OVERLAP", "MVA_B_OVERLAP", fNbins, xmin, xmax );
1232 }
1233 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1234 const Event* ev = GetEvent(ievt);
1235 Float_t w = (eTT==Types::kTesting ? ev->GetWeight() : ev->GetOriginalWeight());
1236 if (DataInfo().IsSignal(ev)) mva_s->Fill( (*mvaRes)[ievt], w );
1237 else mva_b->Fill( (*mvaRes)[ievt], w );
1238
1239 if (CalcOverlapIntergral) {
1240 Float_t w_ov = ev->GetWeight();
1241 if (DataInfo().IsSignal(ev))
1242 mva_s_overlap->Fill( (*mvaRes)[ievt], w_ov );
1243 else
1244 mva_b_overlap->Fill( (*mvaRes)[ievt], w_ov );
1245 }
1246 }
1247 gTools().NormHist( mva_s );
1248 gTools().NormHist( mva_b );
1249 PDF *fS = new PDF( "PDF Sig", mva_s, PDF::kSpline2 );
1250 PDF *fB = new PDF( "PDF Bkg", mva_b, PDF::kSpline2 );
1251
1252 // calculate ROC integral from fS, fB
1254
1255 // calculate overlap integral
1256 if (CalcOverlapIntergral) {
1257 gTools().NormHist( mva_s_overlap );
1258 gTools().NormHist( mva_b_overlap );
1259
1260 fOverlap_integral = 0.0;
1261 for (Int_t bin=1; bin<=mva_s_overlap->GetNbinsX(); bin++){
1262 Double_t bc_s = mva_s_overlap->GetBinContent(bin);
1263 Double_t bc_b = mva_b_overlap->GetBinContent(bin);
1264 if (bc_s > 0.0 && bc_b > 0.0)
1265 fOverlap_integral += TMath::Min(bc_s, bc_b);
1266 }
1267
1268 delete mva_s_overlap;
1269 delete mva_b_overlap;
1270 }
1271
1272 delete mva_s;
1273 delete mva_b;
1274 delete fS;
1275 delete fB;
1276 if (!(singleMethod && eTT==Types::kTraining)) delete mvaRes;
1277
1278 Data()->SetCurrentType(Types::kTraining);
1279
1280 return ROC;
1281}
1282
1284{
1285 // Calculate MVA values of current method fMethods.back() on
1286 // training sample
1287
1288 Data()->SetCurrentType(Types::kTraining);
1289 MethodBase* method = dynamic_cast<MethodBase*>(fMethods.back());
1290 if (!method) {
1291 Log() << kFATAL << "dynamic cast to MethodBase* failed" <<Endl;
1292 return;
1293 }
1294 // calculate MVA values
1295 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1296 GetEvent(ievt);
1297 fMVAvalues->at(ievt) = method->GetMvaValue();
1298 }
1299
1300 // fill cumulative mva distribution
1301
1302
1303}
1304
1305////////////////////////////////////////////////////////////////////////////////
1306/// fill various monitoring histograms from information of the individual classifiers that
1307/// have been boosted.
1308/// of course.... this depends very much on the individual classifiers, and so far, only for
1309/// Decision Trees, this monitoring is actually implemented
1310
1312{
1313 Results* results = Data()->GetResults(GetMethodName(), Types::kTraining, GetAnalysisType());
1314
1315 if (GetCurrentMethod(methodIndex)->GetMethodType() == TMVA::Types::kDT) {
1316 TMVA::MethodDT* currentDT=dynamic_cast<TMVA::MethodDT*>(GetCurrentMethod(methodIndex));
1317 if (currentDT){
1318 if (stage == Types::kBoostProcBegin){
1319 results->Store(new TH1I("NodesBeforePruning","nodes before pruning",this->GetBoostNum(),0,this->GetBoostNum()),"NodesBeforePruning");
1320 results->Store(new TH1I("NodesAfterPruning","nodes after pruning",this->GetBoostNum(),0,this->GetBoostNum()),"NodesAfterPruning");
1321 }
1322
1323 if (stage == Types::kBeforeTraining){
1324 }
1325 else if (stage == Types::kBeforeBoosting){
1326 results->GetHist("NodesBeforePruning")->SetBinContent(methodIndex+1,currentDT->GetNNodesBeforePruning());
1327 results->GetHist("NodesAfterPruning")->SetBinContent(methodIndex+1,currentDT->GetNNodes());
1328 }
1329 else if (stage == Types::kAfterBoosting){
1330
1331 }
1332 else if (stage != Types::kBoostProcEnd){
1333 Log() << kINFO << "<Train> average number of nodes before/after pruning : "
1334 << results->GetHist("NodesBeforePruning")->GetMean() << " / "
1335 << results->GetHist("NodesAfterPruning")->GetMean()
1336 << Endl;
1337 }
1338 }
1339
1340 }else if (GetCurrentMethod(methodIndex)->GetMethodType() == TMVA::Types::kFisher) {
1341 if (stage == Types::kAfterBoosting){
1343 }
1344 }else{
1345 if (methodIndex < 3){
1346 Log() << kDEBUG << "No detailed boost monitoring for "
1347 << GetCurrentMethod(methodIndex)->GetMethodName()
1348 << " yet available " << Endl;
1349 }
1350 }
1351
1352 //boosting plots universal for all classifiers 'typically for debug purposes only as they are not general enough'
1353
1354 if (stage == Types::kBeforeBoosting){
1355 // if you want to display the weighted events for 2D case at each boost step:
1356 if (fDetailedMonitoring){
1357 // the following code is useful only for 2D examples - mainly illustration for debug/educational purposes:
1358 if (DataInfo().GetNVariables() == 2) {
1359 results->Store(new TH2F(Form("EventDistSig_%d",methodIndex),Form("EventDistSig_%d",methodIndex),100,0,7,100,0,7));
1360 results->GetHist(Form("EventDistSig_%d",methodIndex))->SetMarkerColor(4);
1361 results->Store(new TH2F(Form("EventDistBkg_%d",methodIndex),Form("EventDistBkg_%d",methodIndex),100,0,7,100,0,7));
1362 results->GetHist(Form("EventDistBkg_%d",methodIndex))->SetMarkerColor(2);
1363
1364 Data()->SetCurrentType(Types::kTraining);
1365 for (Long64_t ievt=0; ievt<GetNEvents(); ievt++) {
1366 const Event* ev = GetEvent(ievt);
1367 Float_t w = ev->GetWeight();
1368 Float_t v0= ev->GetValue(0);
1369 Float_t v1= ev->GetValue(1);
1370 // if (ievt<3) std::cout<<ievt<<" var0="<<v0<<" var1="<<v1<<" weight="<<w<<std::endl;
1371 TH2* h;
1372 if (DataInfo().IsSignal(ev)) h=results->GetHist2D(Form("EventDistSig_%d",methodIndex));
1373 else h=results->GetHist2D(Form("EventDistBkg_%d",methodIndex));
1374 if (h) h->Fill(v0,v1,w);
1375 }
1376 }
1377 }
1378 }
1379
1380 return;
1381}
1382
1383
#define REGISTER_METHOD(CLASS)
for example
#define h(i)
Definition: RSha256.hxx:106
const Bool_t kFALSE
Definition: RtypesCore.h:90
double Double_t
Definition: RtypesCore.h:57
long long Long64_t
Definition: RtypesCore.h:71
float Float_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:89
#define ClassImp(name)
Definition: Rtypes.h:361
int type
Definition: TGX11.cxx:120
float xmin
Definition: THbookFile.cxx:93
float xmax
Definition: THbookFile.cxx:93
char * Form(const char *fmt,...)
Stat_t GetSum() const
Definition: TArrayD.h:46
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
Describe directory structure in memory.
Definition: TDirectory.h:40
virtual TDirectory * GetDirectory(const char *namecycle, Bool_t printError=false, const char *funcname="GetDirectory")
Find a directory using apath.
Definition: TDirectory.cxx:401
virtual TDirectory * mkdir(const char *name, const char *title="", Bool_t returnExistingDirectory=kFALSE)
Create a sub-directory "a" or a hierarchy of sub-directories "a/b/c/...".
virtual Bool_t cd(const char *path=nullptr)
Change current directory to "this" directory.
Definition: TDirectory.cxx:498
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:614
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:571
1-D histogram with an int per channel (see TH1 documentation)}
Definition: TH1.h:530
The TH1 histogram class.
Definition: TH1.h:56
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
Definition: TH1.cxx:7086
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:316
virtual Int_t GetNbinsX() const
Definition: TH1.h:292
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3275
TAxis * GetYaxis()
Definition: TH1.h:317
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8678
virtual Double_t GetBinLowEdge(Int_t bin) const
Return bin lower edge for 1D histogram.
Definition: TH1.cxx:8608
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4907
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:251
Service class for 2-Dim histogram classes.
Definition: TH2.h:30
IMethod * Create(const std::string &name, const TString &job, const TString &title, DataSetInfo &dsi, const TString &option)
creates the method if needed based on the method name using the creator function the factory has stor...
static ClassifierFactory & Instance()
access to the ClassifierFactory singleton creates the instance if needed
TString fWeightFileDir
Definition: Config.h:124
class TMVA::Config::VariablePlotting fVariablePlotting
IONames & GetIONames()
Definition: Config.h:100
Class that contains all the data information.
Definition: DataSetInfo.h:60
void ScaleBoostWeight(Double_t s) const
Definition: Event.h:112
void SetBoostWeight(Double_t w) const
Definition: Event.h:111
Float_t GetValue(UInt_t ivar) const
return value of i'th variable
Definition: Event.cxx:236
Double_t GetOriginalWeight() const
Definition: Event.h:84
Double_t GetWeight() const
return the event weight - depending on whether the flag IgnoreNegWeightsInTraining is or not.
Definition: Event.cxx:381
UInt_t GetClass() const
Definition: Event.h:86
Implementation of the GiniIndex as separation criterion.
Definition: GiniIndex.h:63
Interface for all concrete MVA method implementations.
Definition: IMethod.h:54
Virtual base Class for all MVA method.
Definition: MethodBase.h:111
void SetSilentFile(Bool_t status)
Definition: MethodBase.h:377
void SetWeightFileDir(TString fileDir)
set directory of weight file
virtual void DeclareCompatibilityOptions()
options that are used ONLY for the READER to ensure backward compatibility they are hence without any...
Definition: MethodBase.cxx:598
virtual void WriteEvaluationHistosToFile(Types::ETreeType treetype)
writes all MVA evaluation histograms to file
TDirectory * BaseDir() const
returns the ROOT directory where info/histograms etc of the corresponding MVA method instance are sto...
virtual Bool_t IsSignalLike()
uses a pre-set cut on the MVA output (SetSignalReferenceCut and SetSignalReferenceCutOrientation) for...
Definition: MethodBase.cxx:856
virtual void TestClassification()
initialization
virtual Double_t GetMvaValue(Double_t *errLower=0, Double_t *errUpper=0)=0
Types::EMVA GetMethodType() const
Definition: MethodBase.h:332
void SetSignalReferenceCut(Double_t cut)
Definition: MethodBase.h:363
void SetSignalReferenceCutOrientation(Double_t cutOrientation)
Definition: MethodBase.h:364
void SetModelPersistence(Bool_t status)
Definition: MethodBase.h:381
Double_t GetSignalReferenceCut() const
Definition: MethodBase.h:359
virtual Double_t GetROCIntegral(TH1D *histS, TH1D *histB) const
calculate the area (integral) under the ROC curve as a overall quality measure of the classification
friend class MethodBoost
Definition: MethodBase.h:116
Class for boosting a TMVA method.
Definition: MethodBoost.h:58
void MonitorBoost(Types::EBoostStage stage, UInt_t methodIdx=0)
fill various monitoring histograms from information of the individual classifiers that have been boos...
void ResetBoostWeights()
resetting back the boosted weights of the events to 1
void SingleTrain()
initialization
DataSetManager * fDataSetManager
Definition: MethodBoost.h:193
virtual void WriteEvaluationHistosToFile(Types::ETreeType treetype)
writes all MVA evaluation histograms to file
void CreateMVAHistorgrams()
Bool_t fHistoricBoolOption
Definition: MethodBoost.h:195
void WriteMonitoringHistosToFile(void) const
write special monitoring histograms to file dummy implementation here --------------—
Double_t AdaBoost(MethodBase *method, Bool_t useYesNoLeaf)
the standard (discrete or real) AdaBoost algorithm
Bool_t BookMethod(Types::EMVA theMethod, TString methodTitle, TString theOption)
just registering the string from which the boosted classifier will be created
void CheckSetup()
check may be overridden by derived class (sometimes, eg, fitters are used which can only be implement...
virtual void TestClassification()
initialization
void InitHistos()
initialisation routine
void ProcessOptions()
process user options
Double_t GetBoostROCIntegral(Bool_t, Types::ETreeType, Bool_t CalcOverlapIntergral=kFALSE)
Calculate the ROC integral of a single classifier or even the whole boosted classifier.
Double_t SingleBoost(MethodBase *method)
Double_t GetMvaValue(Double_t *err=0, Double_t *errUpper=0)
return boosted MVA response
Double_t Bagging()
Bagging or Bootstrap boosting, gives new random poisson weight for every event.
virtual Bool_t HasAnalysisType(Types::EAnalysisType type, UInt_t numberClasses, UInt_t)
Boost can handle classification with 2 classes and regression with one regression-target.
const Ranking * CreateRanking()
void DeclareCompatibilityOptions()
options that are used ONLY for the READER to ensure backward compatibility they are hence without any...
std::vector< Float_t > * fMVAvalues
Definition: MethodBoost.h:191
virtual ~MethodBoost(void)
destructor
void FindMVACut(MethodBase *method)
find the CUT on the individual MVA that defines an event as correct or misclassified (to be used in t...
void GetHelpMessage() const
Get help message text.
Class for categorizing the phase space.
DataSetManager * fDataSetManager
Virtual base class for combining several TMVA method.
std::vector< IMethod * > fMethods
Analysis of Boosted Decision Trees.
Definition: MethodDT.h:49
Int_t GetNNodes()
Definition: MethodDT.h:97
Int_t GetNNodesBeforePruning()
Definition: MethodDT.h:96
static void InhibitOutput()
Definition: MsgLogger.cxx:74
static void EnableOutput()
Definition: MsgLogger.cxx:75
PDF wrapper for histograms; uses user-defined spline interpolation.
Definition: PDF.h:63
@ kSpline2
Definition: PDF.h:70
Double_t GetMVAProbAt(Double_t value)
void AddEvent(Double_t val, Double_t weight, Int_t type)
Ranking for variables in method (implementation)
Definition: Ranking.h:48
Class that is the base-class for a vector of result.
Definition: Results.h:57
TH2 * GetHist2D(const TString &alias) const
Definition: Results.cxx:145
TH1 * GetHist(const TString &alias) const
Definition: Results.cxx:136
void Store(TObject *obj, const char *alias=0)
Definition: Results.cxx:86
An interface to calculate the "SeparationGain" for different separation criteria used in various trai...
virtual Double_t GetSeparationGain(const Double_t nSelS, const Double_t nSelB, const Double_t nTotS, const Double_t nTotB)
Separation Gain: the measure of how the quality of separation of the sample increases by splitting th...
virtual Double_t GetSeparationIndex(const Double_t s, const Double_t b)=0
Timing information for training and evaluation of MVA methods.
Definition: Timer.h:58
TString GetElapsedTime(Bool_t Scientific=kTRUE)
returns pretty string with elapsed time
Definition: Timer.cxx:147
void DrawProgressBar(Int_t, const TString &comment="")
draws progress bar in color or B&W caution:
Definition: Timer.cxx:203
void ComputeStat(const std::vector< TMVA::Event * > &, std::vector< Float_t > *, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Int_t signalClass, Bool_t norm=kFALSE)
sanity check
Definition: Tools.cxx:213
const TString & Color(const TString &)
human readable color strings
Definition: Tools.cxx:839
Double_t NormHist(TH1 *theHist, Double_t norm=1.0)
normalises histogram
Definition: Tools.cxx:394
Singleton class for Global types used by TMVA.
Definition: Types.h:73
@ kBoostProcBegin
Definition: Types.h:152
@ kBeforeBoosting
Definition: Types.h:154
@ kAfterBoosting
Definition: Types.h:155
@ kBoostProcEnd
Definition: Types.h:156
@ kBeforeTraining
Definition: Types.h:153
TString GetMethodName(Types::EMVA method) const
Definition: Types.cxx:136
static Types & Instance()
the the single instance of "Types" if existing already, or create it (Singleton)
Definition: Types.cxx:70
@ kFisher
Definition: Types.h:84
@ kCategory
Definition: Types.h:99
EAnalysisType
Definition: Types.h:127
@ kClassification
Definition: Types.h:128
@ kMaxTreeType
Definition: Types.h:146
@ kTraining
Definition: Types.h:144
@ kTesting
Definition: Types.h:145
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
virtual void Delete(Option_t *option="")
Delete this object.
Definition: TObject.cxx:169
Random number generator class based on M.
Definition: TRandom3.h:27
virtual Double_t PoissonD(Double_t mean)
Generates a random number according to a Poisson law.
Definition: TRandom.cxx:443
Basic string class.
Definition: TString.h:131
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1125
A TTree represents a columnar dataset.
Definition: TTree.h:78
RooCmdArg Timer(Bool_t flag=kTRUE)
void GetMethodName(TString &name, TKey *mkey)
Definition: tmvaglob.cxx:335
create variable transformations
Config & gConfig()
Tools & gTools()
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:158
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculate a gaussian function with mean and sigma.
Definition: TMath.cxx:448
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
Double_t Exp(Double_t x)
Definition: TMath.h:717
Double_t Log(Double_t x)
Definition: TMath.h:750
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
auto * m
Definition: textangle.C:8
static long int sum(long int i)
Definition: Factory.cxx:2275
REAL epsilon
Definition: triangle.c:617