TMVAClassification.C
Go to the documentation of this file.
1 // @(#)root/tmva $Id$
2 /**********************************************************************************
3  * Project : TMVA - a ROOT-integrated toolkit for multivariate data analysis *
4  * Package : TMVA *
5  * Root Macro: TMVAClassification *
6  * *
7  * This macro provides examples for the training and testing of the *
8  * TMVA classifiers. *
9  * *
10  * As input data is used a toy-MC sample consisting of four Gaussian-distributed *
11  * and linearly correlated input variables. *
12  * *
13  * The methods to be used can be switched on and off by means of booleans, or *
14  * via the prompt command, for example: *
15  * *
16  * root -l ./TMVAClassification.C\(\"Fisher,Likelihood\"\) *
17  * *
18  * (note that the backslashes are mandatory) *
19  * If no method given, a default set of classifiers is used. *
20  * *
21  * The output file "TMVA.root" can be analysed with the use of dedicated *
22  * macros (simply say: root -l <macro.C>), which can be conveniently *
23  * invoked through a GUI that will appear at the end of the run of this macro. *
24  * Launch the GUI via the command: *
25  * *
26  * root -l ./TMVAGui.C *
27  * *
28  **********************************************************************************/
29 
30 #include <cstdlib>
31 #include <iostream>
32 #include <map>
33 #include <string>
34 
35 #include "TChain.h"
36 #include "TFile.h"
37 #include "TTree.h"
38 #include "TString.h"
39 #include "TObjString.h"
40 #include "TSystem.h"
41 #include "TROOT.h"
42 
43 
44 #if not defined(__CINT__) || defined(__MAKECINT__)
45 // needs to be included when makecint runs (ACLIC)
46 #include "TMVA/Factory.h"
47 #include "TMVA/Tools.h"
48 #endif
49 
50 void TMVAClassification( TString myMethodList = "" )
51 {
52  // The explicit loading of the shared libTMVA is done in TMVAlogon.C, defined in .rootrc
53  // if you use your private .rootrc, or run from a different directory, please copy the
54  // corresponding lines from .rootrc
55 
56  // methods to be processed can be given as an argument; use format:
57  //
58  // mylinux~> root -l TMVAClassification.C\(\"myMethod1,myMethod2,myMethod3\"\)
59  //
60  // if you like to use a method via the plugin mechanism, we recommend using
61  //
62  // mylinux~> root -l TMVAClassification.C\(\"P_myMethod\"\)
63  // (an example is given for using the BDT as plugin (see below),
64  // but of course the real application is when you write your own
65  // method based)
66 
67  //---------------------------------------------------------------
68  // This loads the library
69  TMVA::Tools::Instance();
70 
71  // to get access to the GUI and all tmva macros
72  TString tmva_dir(TString(gRootDir) + "/tmva");
73  if(gSystem->Getenv("TMVASYS"))
74  tmva_dir = TString(gSystem->Getenv("TMVASYS"));
75  gROOT->SetMacroPath(tmva_dir + "/test/:" + gROOT->GetMacroPath() );
76  gROOT->ProcessLine(".L TMVAGui.C");
77 
78 
79 
80  // Default MVA methods to be trained + tested
81  std::map<std::string,int> Use;
82 
83  // --- Cut optimisation
84  Use["Cuts"] = 0;
85  Use["CutsD"] = 0;
86  Use["CutsPCA"] = 0;
87  Use["CutsGA"] = 0;
88  Use["CutsSA"] = 0;
89  //
90  // --- 1-dimensional likelihood ("naive Bayes estimator")
91  Use["Likelihood"] = 0;
92  Use["LikelihoodD"] = 0; // the "D" extension indicates decorrelated input variables (see option strings)
93  Use["LikelihoodPCA"] = 0; // the "PCA" extension indicates PCA-transformed input variables (see option strings)
94  Use["LikelihoodKDE"] = 0;
95  Use["LikelihoodMIX"] = 0;
96  //
97  // --- Mutidimensional likelihood and Nearest-Neighbour methods
98  Use["PDERS"] = 0;
99  Use["PDERSD"] = 0;
100  Use["PDERSPCA"] = 0;
101  Use["PDEFoam"] = 0;
102  Use["PDEFoamBoost"] = 0; // uses generalised MVA method boosting
103  Use["KNN"] = 0; // k-nearest neighbour method
104  //
105  // --- Linear Discriminant Analysis
106  Use["LD"] = 0; // Linear Discriminant identical to Fisher
107  Use["Fisher"] = 0;
108  Use["FisherG"] = 0;
109  Use["BoostedFisher"] = 0; // uses generalised MVA method boosting
110  Use["HMatrix"] = 0;
111  //
112  // --- Function Discriminant analysis
113  Use["FDA_GA"] = 0; // minimisation of user-defined function using Genetics Algorithm
114  Use["FDA_SA"] = 0;
115  Use["FDA_MC"] = 0;
116  Use["FDA_MT"] = 0;
117  Use["FDA_GAMT"] = 0;
118  Use["FDA_MCMT"] = 0;
119  //
120  // --- Neural Networks (all are feed-forward Multilayer Perceptrons)
121  Use["MLP"] = 0; // Recommended ANN
122  Use["MLPBFGS"] = 0; // Recommended ANN with optional training method
123  Use["MLPBNN"] = 0; // Recommended ANN with BFGS training method and bayesian regulator
124  Use["CFMlpANN"] = 0; // Depreciated ANN from ALEPH
125  Use["TMlpANN"] = 0; // ROOT's own ANN
126  //
127  // --- Support Vector Machine
128  Use["SVM"] = 0;
129  //
130  // --- Boosted Decision Trees
131  Use["BDT"] = 0; // uses Adaptive Boost
132  Use["BDTG"] = 1; // uses Gradient Boost
133  Use["BDTB"] = 0; // uses Bagging
134  Use["BDTD"] = 0; // decorrelation + Adaptive Boost
135  Use["BDTF"] = 0; // allow usage of fisher discriminant for node splitting
136  //
137  // --- Friedman's RuleFit method, ie, an optimised series of cuts ("rules")
138  Use["RuleFit"] = 0;
139  // ---------------------------------------------------------------
140 
141  std::cout << std::endl;
142  std::cout << "==> Start TMVAClassification" << std::endl;
143 
144  // Select methods (don't look at this code - not of interest)
145  if (myMethodList != "") {
146  for (std::map<std::string,int>::iterator it = Use.begin(); it != Use.end(); it++) it->second = 0;
147 
148  std::vector<TString> mlist = TMVA::gTools().SplitString( myMethodList, ',' );
149  for (UInt_t i=0; i<mlist.size(); i++) {
150  std::string regMethod(mlist[i]);
151 
152  if (Use.find(regMethod) == Use.end()) {
153  std::cout << "Method \"" << regMethod << "\" not known in TMVA under this name. Choose among the following:" << std::endl;
154  for (std::map<std::string,int>::iterator it = Use.begin(); it != Use.end(); it++) std::cout << it->first << " ";
155  std::cout << std::endl;
156  return;
157  }
158  Use[regMethod] = 1;
159  }
160  }
161 
162  // --------------------------------------------------------------------------------------------------
163 
164  // --- Here the preparation phase begins
165 
166  // Create a ROOT output file where TMVA will store ntuples, histograms, etc.
167  TString outfileName( "TMVA.root" );
168  TFile* outputFile = TFile::Open( outfileName, "RECREATE" );
169 
170  // Create the factory object. Later you can choose the methods
171  // whose performance you'd like to investigate. The factory is
172  // the only TMVA object you have to interact with
173  //
174  // The first argument is the base of the name of all the
175  // weightfiles in the directory weight/
176  //
177  // The second argument is the output file for the training results
178  // All TMVA output can be suppressed by removing the "!" (not) in
179  // front of the "Silent" argument in the option string
180  TMVA::Factory *factory = new TMVA::Factory( "TMVAClassification", outputFile,
181  "!V:!Silent:Color:DrawProgressBar:Transformations=I;D;P;G,D:AnalysisType=Classification" );
182 
183  // If you wish to modify default settings
184  // (please check "src/Config.h" to see all available global options)
185  // (TMVA::gConfig().GetVariablePlotting()).fTimesRMS = 8.0;
186  // (TMVA::gConfig().GetIONames()).fWeightFileDir = "myWeightDirectory";
187 
188  // Define the input variables that shall be used for the MVA training
189  // note that you may also use variable expressions, such as: "3*var1/var2*abs(var3)"
190  // [all types of expressions that can also be parsed by TTree::Draw( "expression" )]
191  factory->AddVariable( "TMath::ACos(trk.kalman.tracks[0].dir.z)*180./TMath::Pi()","Theta","",'F');
192  factory->AddVariable( "energy.numu.hadclust.nhit > 0? energy.numu.hadcalE / energy.numu.hadclust.nhit : 0", "hadEhit", "GeV", 'F' );
193  factory->AddVariable( "energy.numu.trkccE > 0? energy.numu.recotrkcchadE/energy.numu.trkccE : 0", "hadEFrac", "", 'F' );
194  factory->AddVariable( "energy.numu.recomuonE > 0? energy.numu.recomuonE : 0", "muonE", "GeV", 'F' );
195  factory->AddVariable("trk.kalman.tracks[0].nplanegap","nGaps","",'F');
196 
197  // You can add so-called "Spectator variables", which are not used in the MVA training,
198  // but will appear in the final "TestTree" produced by TMVA. This TestTree will contain the
199  // input variables, the response values of all trained MVAs, and the spectator variables
200  factory->AddSpectator( "mc.nu.pdg", "PDG code", "", 'F' );
201 
202  // Register signal and background trees
203  TChain *signal = new TChain("recTree","signal");
204  TChain *background = new TChain("recTree","background");
205 
206  // signal->Add("RHC_Files/nd_nonswap_rhc.1?_of_100.root");
207  // background->Add("RHC_Files/nd_nonswap_rhc.1?_of_100.root");
208 
209  signal->Add("/nova/ana/users/bzamoran/RHC_decaf/nd_nonswap_rhc.vol1.root");
210  background->Add("/nova/ana/users/bzamoran/RHC_decaf/nd_nonswap_rhc.vol1.root");
211 
212  // global event weights per tree (see below for setting event-wise weights)
213  Double_t signalWeight = 1.0;
214  Double_t backgroundWeight = 7.0;
215 
216  // You can add an arbitrary number of signal or background trees
217  factory->AddSignalTree ( signal, signalWeight );
218  factory->AddBackgroundTree( background, backgroundWeight );
219 
220  //factory->SetBackgroundWeightExpression( "weight" );
221 
222  // Apply additional cuts on the signal and background samples (can be different)
223  TCut mycuts = "mc.nu.pdg == -14";
224  TCut mycutb = "mc.nu.pdg == 14";
225 
226  // Tell the factory how to use the training and testing events
227  //
228  // If no numbers of events are given, half of the events in the tree are used
229  // for training, and the other half for testing:
230  // factory->PrepareTrainingAndTestTree( mycut, "SplitMode=random:!V" );
231  // To also specify the number of testing events, use:
232  // factory->PrepareTrainingAndTestTree( mycut,
233  // "NSigTrain=3000:NBkgTrain=3000:NSigTest=3000:NBkgTest=3000:SplitMode=Random:!V" );
234  factory->PrepareTrainingAndTestTree( mycuts, mycutb,
235  "nTrain_Signal=0:nTrain_Background=0:SplitMode=Random:NormMode=NumEvents:!V" );
236 
237  // ---- Book MVA methods
238  //
239  // Please lookup the various method configuration options in the corresponding cxx files, eg:
240  // src/MethoCuts.cxx, etc, or here: http://tmva.sourceforge.net/optionRef.html
241  // it is possible to preset ranges in the option string in which the cut optimisation should be done:
242  // "...:CutRangeMin[2]=-1:CutRangeMax[2]=1"...", where [2] is the third input variable
243 
244  // Cut optimisation
245  if (Use["Cuts"])
246  factory->BookMethod( TMVA::Types::kCuts, "Cuts",
247  "!H:!V:FitMethod=MC:EffSel:SampleSize=200000:VarProp=FSmart" );
248 
249  if (Use["CutsD"])
250  factory->BookMethod( TMVA::Types::kCuts, "CutsD",
251  "!H:!V:FitMethod=MC:EffSel:SampleSize=200000:VarProp=FSmart:VarTransform=Decorrelate" );
252 
253  if (Use["CutsPCA"])
254  factory->BookMethod( TMVA::Types::kCuts, "CutsPCA",
255  "!H:!V:FitMethod=MC:EffSel:SampleSize=200000:VarProp=FSmart:VarTransform=PCA" );
256 
257  if (Use["CutsGA"])
258  factory->BookMethod( TMVA::Types::kCuts, "CutsGA",
259  "H:!V:FitMethod=GA:CutRangeMin[0]=-10:CutRangeMax[0]=10:VarProp[1]=FMax:EffSel:Steps=30:Cycles=3:PopSize=400:SC_steps=10:SC_rate=5:SC_factor=0.95" );
260 
261  if (Use["CutsSA"])
262  factory->BookMethod( TMVA::Types::kCuts, "CutsSA",
263  "!H:!V:FitMethod=SA:EffSel:MaxCalls=150000:KernelTemp=IncAdaptive:InitialTemp=1e+6:MinTemp=1e-6:Eps=1e-10:UseDefaultScale" );
264 
265  // Likelihood ("naive Bayes estimator")
266  if (Use["Likelihood"])
267  factory->BookMethod( TMVA::Types::kLikelihood, "Likelihood",
268  "H:!V:TransformOutput:PDFInterpol=Spline2:NSmoothSig[0]=20:NSmoothBkg[0]=20:NSmoothBkg[1]=10:NSmooth=1:NAvEvtPerBin=50" );
269 
270  // Decorrelated likelihood
271  if (Use["LikelihoodD"])
272  factory->BookMethod( TMVA::Types::kLikelihood, "LikelihoodD",
273  "!H:!V:TransformOutput:PDFInterpol=Spline2:NSmoothSig[0]=20:NSmoothBkg[0]=20:NSmooth=5:NAvEvtPerBin=50:VarTransform=Decorrelate" );
274 
275  // PCA-transformed likelihood
276  if (Use["LikelihoodPCA"])
277  factory->BookMethod( TMVA::Types::kLikelihood, "LikelihoodPCA",
278  "!H:!V:!TransformOutput:PDFInterpol=Spline2:NSmoothSig[0]=20:NSmoothBkg[0]=20:NSmooth=5:NAvEvtPerBin=50:VarTransform=PCA" );
279 
280  // Use a kernel density estimator to approximate the PDFs
281  if (Use["LikelihoodKDE"])
282  factory->BookMethod( TMVA::Types::kLikelihood, "LikelihoodKDE",
283  "!H:!V:!TransformOutput:PDFInterpol=KDE:KDEtype=Gauss:KDEiter=Adaptive:KDEFineFactor=0.3:KDEborder=None:NAvEvtPerBin=50" );
284 
285  // Use a variable-dependent mix of splines and kernel density estimator
286  if (Use["LikelihoodMIX"])
287  factory->BookMethod( TMVA::Types::kLikelihood, "LikelihoodMIX",
288  "!H:!V:!TransformOutput:PDFInterpolSig[0]=KDE:PDFInterpolBkg[0]=KDE:PDFInterpolSig[1]=KDE:PDFInterpolBkg[1]=KDE:PDFInterpolSig[2]=Spline2:PDFInterpolBkg[2]=Spline2:PDFInterpolSig[3]=Spline2:PDFInterpolBkg[3]=Spline2:KDEtype=Gauss:KDEiter=Nonadaptive:KDEborder=None:NAvEvtPerBin=50" );
289 
290  // Test the multi-dimensional probability density estimator
291  // here are the options strings for the MinMax and RMS methods, respectively:
292  // "!H:!V:VolumeRangeMode=MinMax:DeltaFrac=0.2:KernelEstimator=Gauss:GaussSigma=0.3" );
293  // "!H:!V:VolumeRangeMode=RMS:DeltaFrac=3:KernelEstimator=Gauss:GaussSigma=0.3" );
294  if (Use["PDERS"])
295  factory->BookMethod( TMVA::Types::kPDERS, "PDERS",
296  "!H:!V:NormTree=T:VolumeRangeMode=Adaptive:KernelEstimator=Gauss:GaussSigma=0.3:NEventsMin=400:NEventsMax=600" );
297 
298  if (Use["PDERSD"])
299  factory->BookMethod( TMVA::Types::kPDERS, "PDERSD",
300  "!H:!V:VolumeRangeMode=Adaptive:KernelEstimator=Gauss:GaussSigma=0.3:NEventsMin=400:NEventsMax=600:VarTransform=Decorrelate" );
301 
302  if (Use["PDERSPCA"])
303  factory->BookMethod( TMVA::Types::kPDERS, "PDERSPCA",
304  "!H:!V:VolumeRangeMode=Adaptive:KernelEstimator=Gauss:GaussSigma=0.3:NEventsMin=400:NEventsMax=600:VarTransform=PCA" );
305 
306  // Multi-dimensional likelihood estimator using self-adapting phase-space binning
307  if (Use["PDEFoam"])
308  factory->BookMethod( TMVA::Types::kPDEFoam, "PDEFoam",
309  "!H:!V:SigBgSeparate=F:TailCut=0.001:VolFrac=0.0666:nActiveCells=500:nSampl=2000:nBin=5:Nmin=100:Kernel=None:Compress=T" );
310 
311  if (Use["PDEFoamBoost"])
312  factory->BookMethod( TMVA::Types::kPDEFoam, "PDEFoamBoost",
313  "!H:!V:Boost_Num=30:Boost_Transform=linear:SigBgSeparate=F:MaxDepth=4:UseYesNoCell=T:DTLogic=MisClassificationError:FillFoamWithOrigWeights=F:TailCut=0:nActiveCells=500:nBin=20:Nmin=400:Kernel=None:Compress=T" );
314 
315  // K-Nearest Neighbour classifier (KNN)
316  if (Use["KNN"])
317  factory->BookMethod( TMVA::Types::kKNN, "KNN",
318  "H:nkNN=20:ScaleFrac=0.8:SigmaFact=1.0:Kernel=Gaus:UseKernel=F:UseWeight=T:!Trim" );
319 
320  // H-Matrix (chi2-squared) method
321  if (Use["HMatrix"])
322  factory->BookMethod( TMVA::Types::kHMatrix, "HMatrix", "!H:!V:VarTransform=None" );
323 
324  // Linear discriminant (same as Fisher discriminant)
325  if (Use["LD"])
326  factory->BookMethod( TMVA::Types::kLD, "LD", "H:!V:VarTransform=None:CreateMVAPdfs:PDFInterpolMVAPdf=Spline2:NbinsMVAPdf=50:NsmoothMVAPdf=10" );
327 
328  // Fisher discriminant (same as LD)
329  if (Use["Fisher"])
330  factory->BookMethod( TMVA::Types::kFisher, "Fisher", "H:!V:Fisher:VarTransform=None:CreateMVAPdfs:PDFInterpolMVAPdf=Spline2:NbinsMVAPdf=50:NsmoothMVAPdf=10" );
331 
332  // Fisher with Gauss-transformed input variables
333  if (Use["FisherG"])
334  factory->BookMethod( TMVA::Types::kFisher, "FisherG", "H:!V:VarTransform=Gauss" );
335 
336  // Composite classifier: ensemble (tree) of boosted Fisher classifiers
337  if (Use["BoostedFisher"])
338  factory->BookMethod( TMVA::Types::kFisher, "BoostedFisher",
339  "H:!V:Boost_Num=20:Boost_Transform=log:Boost_Type=AdaBoost:Boost_AdaBoostBeta=0.2:!Boost_DetailedMonitoring" );
340 
341  // Function discrimination analysis (FDA) -- test of various fitters - the recommended one is Minuit (or GA or SA)
342  if (Use["FDA_MC"])
343  factory->BookMethod( TMVA::Types::kFDA, "FDA_MC",
344  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=MC:SampleSize=100000:Sigma=0.1" );
345 
346  if (Use["FDA_GA"]) // can also use Simulated Annealing (SA) algorithm (see Cuts_SA options])
347  factory->BookMethod( TMVA::Types::kFDA, "FDA_GA",
348  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=GA:PopSize=300:Cycles=3:Steps=20:Trim=True:SaveBestGen=1" );
349 
350  if (Use["FDA_SA"]) // can also use Simulated Annealing (SA) algorithm (see Cuts_SA options])
351  factory->BookMethod( TMVA::Types::kFDA, "FDA_SA",
352  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=SA:MaxCalls=15000:KernelTemp=IncAdaptive:InitialTemp=1e+6:MinTemp=1e-6:Eps=1e-10:UseDefaultScale" );
353 
354  if (Use["FDA_MT"])
355  factory->BookMethod( TMVA::Types::kFDA, "FDA_MT",
356  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=MINUIT:ErrorLevel=1:PrintLevel=-1:FitStrategy=2:UseImprove:UseMinos:SetBatch" );
357 
358  if (Use["FDA_GAMT"])
359  factory->BookMethod( TMVA::Types::kFDA, "FDA_GAMT",
360  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=GA:Converger=MINUIT:ErrorLevel=1:PrintLevel=-1:FitStrategy=0:!UseImprove:!UseMinos:SetBatch:Cycles=1:PopSize=5:Steps=5:Trim" );
361 
362  if (Use["FDA_MCMT"])
363  factory->BookMethod( TMVA::Types::kFDA, "FDA_MCMT",
364  "H:!V:Formula=(0)+(1)*x0+(2)*x1+(3)*x2+(4)*x3:ParRanges=(-1,1);(-10,10);(-10,10);(-10,10);(-10,10):FitMethod=MC:Converger=MINUIT:ErrorLevel=1:PrintLevel=-1:FitStrategy=0:!UseImprove:!UseMinos:SetBatch:SampleSize=20" );
365 
366  // TMVA ANN: MLP (recommended ANN) -- all ANNs in TMVA are Multilayer Perceptrons
367  if (Use["MLP"])
368  factory->BookMethod( TMVA::Types::kMLP, "MLP", "H:!V:NeuronType=tanh:VarTransform=N:NCycles=600:HiddenLayers=N+5:TestRate=5:!UseRegulator" );
369 
370  if (Use["MLPBFGS"])
371  factory->BookMethod( TMVA::Types::kMLP, "MLPBFGS", "H:!V:NeuronType=tanh:VarTransform=N:NCycles=600:HiddenLayers=N+5:TestRate=5:TrainingMethod=BFGS:!UseRegulator" );
372 
373  if (Use["MLPBNN"])
374  factory->BookMethod( TMVA::Types::kMLP, "MLPBNN", "H:!V:NeuronType=tanh:VarTransform=N:NCycles=600:HiddenLayers=N+5:TestRate=5:TrainingMethod=BFGS:UseRegulator" ); // BFGS training with bayesian regulators
375 
376  // CF(Clermont-Ferrand)ANN
377  if (Use["CFMlpANN"])
378  factory->BookMethod( TMVA::Types::kCFMlpANN, "CFMlpANN", "!H:!V:NCycles=2000:HiddenLayers=N+1,N" ); // n_cycles:#nodes:#nodes:...
379 
380  // Tmlp(Root)ANN
381  if (Use["TMlpANN"])
382  factory->BookMethod( TMVA::Types::kTMlpANN, "TMlpANN", "!H:!V:NCycles=200:HiddenLayers=N+1,N:LearningMethod=BFGS:ValidationFraction=0.3" ); // n_cycles:#nodes:#nodes:...
383 
384  // Support Vector Machine
385  if (Use["SVM"])
386  factory->BookMethod( TMVA::Types::kSVM, "SVM", "Gamma=0.25:Tol=0.001:VarTransform=Norm" );
387 
388  // Boosted Decision Trees
389  if (Use["BDTG"]) // Gradient Boost
390  factory->BookMethod( TMVA::Types::kBDT, "BDTG",
391  "!H:!V:NTrees=1000:MinNodeSize=2.5%:BoostType=Grad:Shrinkage=0.10:UseBaggedBoost:BaggedSampleFraction=0.5:nCuts=20:MaxDepth=2" );
392 
393  if (Use["BDT"]) // Adaptive Boost
394  factory->BookMethod( TMVA::Types::kBDT, "BDT",
395  "!H:!V:NTrees=850:MinNodeSize=2.5%:MaxDepth=3:BoostType=AdaBoost:AdaBoostBeta=0.5:UseBaggedBoost:BaggedSampleFraction=0.5:SeparationType=GiniIndex:nCuts=20" );
396 
397  if (Use["BDTB"]) // Bagging
398  factory->BookMethod( TMVA::Types::kBDT, "BDTB",
399  "!H:!V:NTrees=400:BoostType=Bagging:SeparationType=GiniIndex:nCuts=20" );
400 
401  if (Use["BDTD"]) // Decorrelation + Adaptive Boost
402  factory->BookMethod( TMVA::Types::kBDT, "BDTD",
403  "!H:!V:NTrees=400:MinNodeSize=5%:MaxDepth=3:BoostType=AdaBoost:SeparationType=GiniIndex:nCuts=20:VarTransform=Decorrelate" );
404 
405  if (Use["BDTF"]) // Allow Using Fisher discriminant in node splitting for (strong) linearly correlated variables
406  factory->BookMethod( TMVA::Types::kBDT, "BDTMitFisher",
407  "!H:!V:NTrees=50:MinNodeSize=2.5%:UseFisherCuts:MaxDepth=3:BoostType=AdaBoost:AdaBoostBeta=0.5:SeparationType=GiniIndex:nCuts=20" );
408 
409  // RuleFit -- TMVA implementation of Friedman's method
410  if (Use["RuleFit"])
411  factory->BookMethod( TMVA::Types::kRuleFit, "RuleFit",
412  "H:!V:RuleFitModule=RFTMVA:Model=ModRuleLinear:MinImp=0.001:RuleMinDist=0.001:NTrees=20:fEventsMin=0.01:fEventsMax=0.5:GDTau=-1.0:GDTauPrec=0.01:GDStep=0.01:GDNSteps=10000:GDErrScale=1.02" );
413 
414  // For an example of the category classifier usage, see: TMVAClassificationCategory
415 
416  // --------------------------------------------------------------------------------------------------
417 
418  // ---- Now you can optimize the setting (configuration) of the MVAs using the set of training events
419 
420  // ---- STILL EXPERIMENTAL and only implemented for BDT's !
421  // factory->OptimizeAllMethods("SigEffAt001","Scan");
422  // factory->OptimizeAllMethods("ROCIntegral","FitGA");
423 
424  // --------------------------------------------------------------------------------------------------
425 
426  // ---- Now you can tell the factory to train, test, and evaluate the MVAs
427 
428  // Train MVAs using the set of training events
429  factory->TrainAllMethods();
430 
431  // ---- Evaluate all MVAs using the set of test events
432  factory->TestAllMethods();
433 
434  // ----- Evaluate and compare performance of all configured MVAs
435  factory->EvaluateAllMethods();
436 
437  // --------------------------------------------------------------
438 
439  // Save the output
440  outputFile->Close();
441 
442  std::cout << "==> Wrote root file: " << outputFile->GetName() << std::endl;
443  std::cout << "==> TMVAClassification is done!" << std::endl;
444 
445  delete factory;
446 
447  // Launch the GUI for the root macros
448  if (!gROOT->IsBatch()) TMVAGui( outfileName );
449 }
set< int >::iterator it
void TMVAClassification(TString myMethodList="")
OStream cout
Definition: OStream.cxx:6
enum BeamMode string