Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
MethodLikelihood.cxx
Go to the documentation of this file.
1// @(#)root/tmva $Id$
2// Author: Andreas Hoecker, Joerg Stelzer, Helge Voss, Kai Voss, Eckhard von Toerne, Jan Therhaag
3
4/**********************************************************************************
5 * Project: TMVA - a Root-integrated toolkit for multivariate Data analysis *
6 * Package: TMVA *
7 * Class : TMVA::MethodLikelihood *
8 * *
9 * *
10 * Description: *
11 * Implementation (see header for description) *
12 * *
13 * Authors (alphabetical): *
14 * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15 * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
16 * Kai Voss <Kai.Voss@cern.ch> - U. of Victoria, Canada *
17 * Jan Therhaag <Jan.Therhaag@cern.ch> - U of Bonn, Germany *
18 * Eckhard v. Toerne <evt@uni-bonn.de> - U. of Bonn, Germany *
19 * *
20 * Copyright (c) 2005-2011: *
21 * CERN, Switzerland *
22 * U. of Victoria, Canada *
23 * MPI-K Heidelberg, Germany *
24 * U. of Bonn, Germany *
25 * *
26 * Redistribution and use in source and binary forms, with or without *
27 * modification, are permitted according to the terms listed in LICENSE *
28 * (see tmva/doc/LICENSE) *
29 **********************************************************************************/
30
31/*! \class TMVA::MethodLikelihood
32\ingroup TMVA
33
34Likelihood analysis ("non-parametric approach")
35
36
37Also implemented is a "diagonalized likelihood approach",
38which improves over the uncorrelated likelihood approach by
39transforming linearly the input variables into a diagonal space,
40using the square-root of the covariance matrix
41
42
43The method of maximum likelihood is the most straightforward, and
44certainly among the most elegant multivariate analyser approaches.
45We define the likelihood ratio, \f$ R_L \f$, for event
46\f$ i \f$, by:
47
48\f[
49R_L(i) = \frac{L_S(i)}{L_B(i) + L_B(i)}
50\f]
51
52Here the signal and background likelihoods, \f$ L_S \f$,
53\f$ L_B \f$, are products of the corresponding probability
54densities, \f$ p_S \f$, \f$ p_B \f$, of the
55\f$ N_{var} \f$ discriminating variables used in the MVA:
56
57\f[
58L_S(i) \ \prod_{j=1}^{N_{var}} p_{Sj} (i)
59\f]
60
61and accordingly for \f$ L_B \f$.
62In practise, TMVA uses polynomial splines to estimate the probability
63density functions (PDF) obtained from the distributions of the
64training variables.
65
66
67Note that in TMVA the output of the likelihood ratio is transformed by:
68
69\f[
70R_L(i) \to R'_L(i) = -\frac{1}{\tau} ln(R_L^{-1}(i) -1)
71\f]
72
73to avoid the occurrence of heavy peaks at \f$ R_L = 0.1 \f$ .
74
75#### Decorrelated (or "diagonalized") Likelihood
76
77The biggest drawback of the Likelihood approach is that it assumes
78that the discriminant variables are uncorrelated. If it were the case,
79it can be proven that the discrimination obtained by the above likelihood
80ratio is optimal, ie, no other method can beat it. However, in most
81practical applications of MVAs correlations are present. </p>
82
83
84Linear correlations, measured from the training sample, can be taken
85into account in a straightforward manner through the square-root
86of the covariance matrix. The square-root of a matrix
87\f$ C \f$ is the matrix \f$ C&prime; \f$ that multiplied with itself
88yields \f$ C \f$: \f$ C \f$=\f$ C&prime;C&prime; \f$. We compute the
89square-root matrix (SQM) by means of diagonalising (\f$ D \f$) the
90covariance matrix:
91
92\f[
93D = S^TCS \Rightarrow C' = S \sqrt{DS^T}
94\f]
95
96and the linear transformation of the linearly correlated into the
97uncorrelated variables space is then given by multiplying the measured
98variable tuple by the inverse of the SQM. Note that these transformations
99are performed for both signal and background separately, since the
100correlation pattern is not the same in the two samples.
101
102
103The above diagonalisation is complete for linearly correlated,
104Gaussian distributed variables only. In real-world examples this
105is not often the case, so that only little additional information
106may be recovered by the diagonalisation procedure. In these cases,
107non-linear methods must be applied.
108*/
109
111
112#include "TMVA/Configurable.h"
114#include "TMVA/DataSet.h"
115#include "TMVA/DataSetInfo.h"
116#include "TMVA/IMethod.h"
117#include "TMVA/MethodBase.h"
118#include "TMVA/MsgLogger.h"
119#include "TMVA/PDF.h"
120#include "TMVA/Ranking.h"
121#include "TMVA/Tools.h"
122#include "TMVA/Types.h"
123#include "TMVA/VariableInfo.h"
124
125#include "TVector.h"
126#include "TMath.h"
127#include "TFile.h"
128#include "TH1.h"
129
130#include <iostream>
131#include <iomanip>
132#include <vector>
133#include <cstdlib>
134
135REGISTER_METHOD(Likelihood)
136
137
138////////////////////////////////////////////////////////////////////////////////
139/// standard constructor
140
142 const TString& methodTitle,
145 TMVA::MethodBase( jobName, Types::kLikelihood, methodTitle, theData, theOption),
146 fEpsilon ( 1.e3 * DBL_MIN ),
147 fTransformLikelihoodOutput( kFALSE ),
148 fDropVariable ( 0 ),
149 fHistSig ( 0 ),
150 fHistBgd ( 0 ),
151 fHistSig_smooth( 0 ),
152 fHistBgd_smooth( 0 ),
153 fDefaultPDFLik ( 0 ),
154 fPDFSig ( 0 ),
155 fPDFBgd ( 0 ),
156 fNsmooth ( 2 ),
157 fNsmoothVarS ( 0 ),
158 fNsmoothVarB ( 0 ),
159 fAverageEvtPerBin( 0 ),
160 fAverageEvtPerBinVarS (0),
161 fAverageEvtPerBinVarB (0),
162 fKDEfineFactor ( 0 ),
163 fInterpolateString(0)
164{
165}
166
167////////////////////////////////////////////////////////////////////////////////
168/// construct likelihood references from file
169
171 const TString& theWeightFile) :
172 TMVA::MethodBase( Types::kLikelihood, theData, theWeightFile),
173 fEpsilon ( 1.e3 * DBL_MIN ),
174 fTransformLikelihoodOutput( kFALSE ),
175 fDropVariable ( 0 ),
176 fHistSig ( 0 ),
177 fHistBgd ( 0 ),
178 fHistSig_smooth( 0 ),
179 fHistBgd_smooth( 0 ),
180 fDefaultPDFLik ( 0 ),
181 fPDFSig ( 0 ),
182 fPDFBgd ( 0 ),
183 fNsmooth ( 2 ),
184 fNsmoothVarS ( 0 ),
185 fNsmoothVarB ( 0 ),
186 fAverageEvtPerBin( 0 ),
187 fAverageEvtPerBinVarS (0),
188 fAverageEvtPerBinVarB (0),
189 fKDEfineFactor ( 0 ),
190 fInterpolateString(0)
191{
192}
193
194////////////////////////////////////////////////////////////////////////////////
195/// destructor
196
198{
199 if (NULL != fDefaultPDFLik) delete fDefaultPDFLik;
200 if (NULL != fHistSig) delete fHistSig;
201 if (NULL != fHistBgd) delete fHistBgd;
202 if (NULL != fHistSig_smooth) delete fHistSig_smooth;
203 if (NULL != fHistBgd_smooth) delete fHistBgd_smooth;
204 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
205 if ((*fPDFSig)[ivar] !=0) delete (*fPDFSig)[ivar];
206 if ((*fPDFBgd)[ivar] !=0) delete (*fPDFBgd)[ivar];
207 }
208 if (NULL != fPDFSig) delete fPDFSig;
209 if (NULL != fPDFBgd) delete fPDFBgd;
210}
211
212////////////////////////////////////////////////////////////////////////////////
213/// FDA can handle classification with 2 classes
214
221
222////////////////////////////////////////////////////////////////////////////////
223/// default initialisation called by all constructors
224
226{
227 // no ranking test
228 fDropVariable = -1;
229 fHistSig = new std::vector<TH1*> ( GetNvar(), (TH1*)0 );
230 fHistBgd = new std::vector<TH1*> ( GetNvar(), (TH1*)0 );
231 fHistSig_smooth = new std::vector<TH1*> ( GetNvar(), (TH1*)0 );
232 fHistBgd_smooth = new std::vector<TH1*> ( GetNvar(), (TH1*)0 );
233 fPDFSig = new std::vector<TMVA::PDF*>( GetNvar(), (TMVA::PDF*)0 );
234 fPDFBgd = new std::vector<TMVA::PDF*>( GetNvar(), (TMVA::PDF*)0 );
235}
236
237////////////////////////////////////////////////////////////////////////////////
238/// define the options (their key words) that can be set in the option string
239///
240/// TransformOutput `<bool>` transform (often strongly peaked) likelihood output through sigmoid inversion
241
243{
244 DeclareOptionRef( fTransformLikelihoodOutput = kFALSE, "TransformOutput",
245 "Transform likelihood output by inverse sigmoid function" );
246
247 // initialize
248
249 // reading every PDF's definition and passing the option string to the next one to be read and marked
250 TString updatedOptions = GetOptions();
251 fDefaultPDFLik = new PDF( TString(GetName()) + " PDF", updatedOptions );
252 fDefaultPDFLik->DeclareOptions();
253 fDefaultPDFLik->ParseOptions();
254 updatedOptions = fDefaultPDFLik->GetOptions();
255 for (UInt_t ivar = 0; ivar< DataInfo().GetNVariables(); ivar++) {
256 (*fPDFSig)[ivar] = new PDF( TString::Format("%s PDF Sig[%d]", GetName(), ivar), updatedOptions,
257 TString::Format("Sig[%d]",ivar), fDefaultPDFLik );
258 (*fPDFSig)[ivar]->DeclareOptions();
259 (*fPDFSig)[ivar]->ParseOptions();
260 updatedOptions = (*fPDFSig)[ivar]->GetOptions();
261 (*fPDFBgd)[ivar] = new PDF( TString::Format("%s PDF Bkg[%d]", GetName(), ivar), updatedOptions,
262 TString::Format("Bkg[%d]",ivar), fDefaultPDFLik );
263 (*fPDFBgd)[ivar]->DeclareOptions();
264 (*fPDFBgd)[ivar]->ParseOptions();
265 updatedOptions = (*fPDFBgd)[ivar]->GetOptions();
266 }
267
268 // the final marked option string is written back to the original likelihood
269 SetOptions( updatedOptions );
270}
271
272
274{
275 // options that are used ONLY for the READER to ensure backward compatibility
276
278 DeclareOptionRef( fNsmooth = 1, "NSmooth",
279 "Number of smoothing iterations for the input histograms");
280 DeclareOptionRef( fAverageEvtPerBin = 50, "NAvEvtPerBin",
281 "Average number of events per PDF bin");
282 DeclareOptionRef( fKDEfineFactor =1. , "KDEFineFactor",
283 "Fine tuning factor for Adaptive KDE: Factor to multiply the width of the kernel");
284 DeclareOptionRef( fBorderMethodString = "None", "KDEborder",
285 "Border effects treatment (1=no treatment , 2=kernel renormalization, 3=sample mirroring)" );
286 DeclareOptionRef( fKDEiterString = "Nonadaptive", "KDEiter",
287 "Number of iterations (1=non-adaptive, 2=adaptive)" );
288 DeclareOptionRef( fKDEtypeString = "Gauss", "KDEtype",
289 "KDE kernel type (1=Gauss)" );
290 fAverageEvtPerBinVarS = new Int_t[GetNvar()];
291 fAverageEvtPerBinVarB = new Int_t[GetNvar()];
292 fNsmoothVarS = new Int_t[GetNvar()];
293 fNsmoothVarB = new Int_t[GetNvar()];
294 fInterpolateString = new TString[GetNvar()];
295 for(UInt_t i=0; i<GetNvar(); ++i) {
296 fAverageEvtPerBinVarS[i] = fAverageEvtPerBinVarB[i] = 0;
297 fNsmoothVarS[i] = fNsmoothVarB[i] = 0;
298 fInterpolateString[i] = "";
299 }
300 DeclareOptionRef( fAverageEvtPerBinVarS, GetNvar(), "NAvEvtPerBinSig",
301 "Average num of events per PDF bin and variable (signal)");
302 DeclareOptionRef( fAverageEvtPerBinVarB, GetNvar(), "NAvEvtPerBinBkg",
303 "Average num of events per PDF bin and variable (background)");
304 DeclareOptionRef(fNsmoothVarS, GetNvar(), "NSmoothSig",
305 "Number of smoothing iterations for the input histograms");
306 DeclareOptionRef(fNsmoothVarB, GetNvar(), "NSmoothBkg",
307 "Number of smoothing iterations for the input histograms");
308 DeclareOptionRef(fInterpolateString, GetNvar(), "PDFInterpol", "Method of interpolating reference histograms (e.g. Spline2 or KDE)");
309}
310
311////////////////////////////////////////////////////////////////////////////////
312/// process user options
313/// reference cut value to distinguish signal-like from background-like events
314
316{
317 SetSignalReferenceCut( TransformLikelihoodOutput( 0.5, 0.5 ) );
318
319 fDefaultPDFLik->ProcessOptions();
320 for (UInt_t ivar = 0; ivar< DataInfo().GetNVariables(); ivar++) {
321 (*fPDFBgd)[ivar]->ProcessOptions();
322 (*fPDFSig)[ivar]->ProcessOptions();
323 }
324}
325
326////////////////////////////////////////////////////////////////////////////////
327/// create reference distributions (PDFs) from signal and background events:
328/// fill histograms and smooth them; if decorrelation is required, compute
329/// corresponding square-root matrices
330/// the reference histograms require the correct boundaries. Since in Likelihood classification
331/// the transformations are applied using both classes, also the corresponding boundaries
332/// need to take this into account
333
335{
336 UInt_t nvar=GetNvar();
337 std::vector<Double_t> xmin(nvar), xmax(nvar);
338 for (UInt_t ivar=0; ivar<nvar; ivar++) {xmin[ivar]=1e30; xmax[ivar]=-1e30;}
339
340 UInt_t nevents=Data()->GetNEvents();
341 for (UInt_t ievt=0; ievt<nevents; ievt++) {
342 // use the true-event-type's transformation
343 // set the event true event types transformation
344 const Event* origEv = Data()->GetEvent(ievt);
345 if (IgnoreEventsWithNegWeightsInTraining() && origEv->GetWeight()<=0) continue;
346 // loop over classes
347 for (int cls=0;cls<2;cls++){
348 GetTransformationHandler().SetTransformationReferenceClass(cls);
349 const Event* ev = GetTransformationHandler().Transform( origEv );
350 for (UInt_t ivar=0; ivar<nvar; ivar++) {
351 Float_t value = ev->GetValue(ivar);
352 if (value < xmin[ivar]) xmin[ivar] = value;
353 if (value > xmax[ivar]) xmax[ivar] = value;
354 }
355 }
356 }
357
358 // create reference histograms
359 // (KDE smoothing requires very finely binned reference histograms)
360 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
361 TString var = (*fInputVars)[ivar];
362
363 // the reference histograms require the correct boundaries. Since in Likelihood classification
364 // the transformations are applied using both classes, also the corresponding boundaries
365 // need to take this into account
366
367 // special treatment for discrete variables
368 if (DataInfo().GetVariableInfo(ivar).GetVarType() == 'I') {
369 // special treatment for integer variables
371 xmax[ivar]=xmax[ivar]+1; // make sure that all entries are included in histogram
373 Int_t nbins = ixmax - ixmin;
374 (*fHistSig)[ivar] = new TH1F(GetMethodName()+"_"+var + "_sig", var + " signal training", nbins, ixmin, ixmax );
375 (*fHistBgd)[ivar] = new TH1F(GetMethodName()+"_"+var + "_bgd", var + " background training", nbins, ixmin, ixmax );
376 } else {
377
378 UInt_t minNEvt = TMath::Min(Data()->GetNEvtSigTrain(),Data()->GetNEvtBkgdTrain());
379 Int_t nbinsS = (*fPDFSig)[ivar]->GetHistNBins( minNEvt );
380 Int_t nbinsB = (*fPDFBgd)[ivar]->GetHistNBins( minNEvt );
381
382 (*fHistSig)[ivar] = new TH1F( TString::Format("%s_%s_%s_sig",DataInfo().GetName(),GetMethodName().Data(),var.Data()),
383 TString::Format("%s_%s_%s signal training",DataInfo().GetName(),GetMethodName().Data(),var.Data()),
384 nbinsS, xmin[ivar], xmax[ivar] );
385 (*fHistBgd)[ivar] = new TH1F( TString::Format("%s_%s_%s_bgd",DataInfo().GetName(),GetMethodName().Data(),var.Data()),
386 TString::Format("%s_%s_%s background training",DataInfo().GetName(),GetMethodName().Data(),var.Data()),
387 nbinsB, xmin[ivar], xmax[ivar] );
388 }
389 }
390
391 // ----- fill the reference histograms
392 Log() << kINFO << "Filling reference histograms" << Endl;
393
394 // event loop
395 for (Int_t ievt=0; ievt<Data()->GetNEvents(); ievt++) {
396
397 // use the true-event-type's transformation
398 // set the event true event types transformation
399 const Event* origEv = Data()->GetEvent(ievt);
400 if (IgnoreEventsWithNegWeightsInTraining() && origEv->GetWeight()<=0) continue;
401 GetTransformationHandler().SetTransformationReferenceClass( origEv->GetClass() );
402 const Event* ev = GetTransformationHandler().Transform( origEv );
403
404 // the event weight
405 Float_t weight = ev->GetWeight();
406
407 // fill variable vector
408 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
409 Double_t value = ev->GetValue(ivar);
410 // verify limits
411 if (value >= xmax[ivar]) value = xmax[ivar] - 1.0e-10;
412 else if (value < xmin[ivar]) value = xmin[ivar] + 1.0e-10;
413 // inserting check if there are events in overflow or underflow
414 if (value >=(*fHistSig)[ivar]->GetXaxis()->GetXmax() ||
415 value <(*fHistSig)[ivar]->GetXaxis()->GetXmin()){
416 Log()<<kWARNING
417 <<"error in filling likelihood reference histograms var="
418 <<(*fInputVars)[ivar]
419 << ", xmin="<<(*fHistSig)[ivar]->GetXaxis()->GetXmin()
420 << ", value="<<value
421 << ", xmax="<<(*fHistSig)[ivar]->GetXaxis()->GetXmax()
422 << Endl;
423 }
424 if (DataInfo().IsSignal(ev)) (*fHistSig)[ivar]->Fill( value, weight );
425 else (*fHistBgd)[ivar]->Fill( value, weight );
426 }
427 }
428
429 // building the pdfs
430 Log() << kINFO << "Building PDF out of reference histograms" << Endl;
431 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
432
433 // the PDF is built from (binned) reference histograms
434 // in case of KDE, this has a large number of bins, which makes it quasi-unbinned
435 (*fPDFSig)[ivar]->BuildPDF( (*fHistSig)[ivar] );
436 (*fPDFBgd)[ivar]->BuildPDF( (*fHistBgd)[ivar] );
437
438 (*fPDFSig)[ivar]->ValidatePDF( (*fHistSig)[ivar] );
439 (*fPDFBgd)[ivar]->ValidatePDF( (*fHistBgd)[ivar] );
440
441 // saving the smoothed histograms
442 if ((*fPDFSig)[ivar]->GetSmoothedHist() != 0) (*fHistSig_smooth)[ivar] = (*fPDFSig)[ivar]->GetSmoothedHist();
443 if ((*fPDFBgd)[ivar]->GetSmoothedHist() != 0) (*fHistBgd_smooth)[ivar] = (*fPDFBgd)[ivar]->GetSmoothedHist();
444 }
445 ExitFromTraining();
446}
447
448////////////////////////////////////////////////////////////////////////////////
449/// returns the likelihood estimator for signal
450/// fill a new Likelihood branch into the testTree
451
453{
454 UInt_t ivar;
455
456 // cannot determine error
457 NoErrorCalc(err, errUpper);
458
459 // retrieve variables, and transform, if required
460 TVector vs( GetNvar() );
461 TVector vb( GetNvar() );
462
463 // need to distinguish signal and background in case of variable transformation
464 // signal first
465
466 GetTransformationHandler().SetTransformationReferenceClass( fSignalClass );
467 // temporary: JS --> FIX
468 //GetTransformationHandler().SetTransformationReferenceClass( 0 );
469 const Event* ev = GetEvent();
470 for (ivar=0; ivar<GetNvar(); ivar++) vs(ivar) = ev->GetValue(ivar);
471
472 GetTransformationHandler().SetTransformationReferenceClass( fBackgroundClass );
473 // temporary: JS --> FIX
474 //GetTransformationHandler().SetTransformationReferenceClass( 1 );
475 ev = GetEvent();
476 for (ivar=0; ivar<GetNvar(); ivar++) vb(ivar) = ev->GetValue(ivar);
477
478 // compute the likelihood (signal)
479 Double_t ps(1), pb(1), p(0);
480 for (ivar=0; ivar<GetNvar(); ivar++) {
481
482 // drop one variable (this is ONLY used for internal variable ranking !)
483 if ((Int_t)ivar == fDropVariable) continue;
484
485 Double_t x[2] = { vs(ivar), vb(ivar) };
486
487 for (UInt_t itype=0; itype < 2; itype++) {
488
489 // verify limits
490 if (x[itype] >= (*fPDFSig)[ivar]->GetXmax()) x[itype] = (*fPDFSig)[ivar]->GetXmax() - 1.0e-10;
491 else if (x[itype] < (*fPDFSig)[ivar]->GetXmin()) x[itype] = (*fPDFSig)[ivar]->GetXmin();
492
493 // find corresponding histogram from cached indices
494 PDF* pdf = (itype == 0) ? (*fPDFSig)[ivar] : (*fPDFBgd)[ivar];
495 if (pdf == 0) Log() << kFATAL << "<GetMvaValue> Reference histograms don't exist" << Endl;
496 TH1* hist = pdf->GetPDFHist();
497
498 // interpolate linearly between adjacent bins
499 // this is not useful for discrete variables
500 Int_t bin = hist->FindBin(x[itype]);
501
502 // **** POTENTIAL BUG: PREFORMANCE IS WORSE WHEN USING TRUE TYPE ***
503 // ==> commented out at present
504 if ((*fPDFSig)[ivar]->GetInterpolMethod() == TMVA::PDF::kSpline0 ||
505 DataInfo().GetVariableInfo(ivar).GetVarType() == 'N') {
506 p = TMath::Max( hist->GetBinContent(bin), fEpsilon );
507 } else { // splined PDF
508 Int_t nextbin = bin;
509 if ((x[itype] > hist->GetBinCenter(bin) && bin != hist->GetNbinsX()) || bin == 1)
510 nextbin++;
511 else
512 nextbin--;
513
514
515 Double_t dx = hist->GetBinCenter(bin) - hist->GetBinCenter(nextbin);
516 Double_t dy = hist->GetBinContent(bin) - hist->GetBinContent(nextbin);
517 Double_t like = hist->GetBinContent(bin) + (x[itype] - hist->GetBinCenter(bin)) * dy/dx;
518
519 p = TMath::Max( like, fEpsilon );
520 }
521
522 if (itype == 0) ps *= p;
523 else pb *= p;
524 }
525 }
526
527 // the likelihood ratio (transform it ?)
528 return TransformLikelihoodOutput( ps, pb );
529}
530
531////////////////////////////////////////////////////////////////////////////////
532/// returns transformed or non-transformed output
533
535{
536 if (ps < fEpsilon) ps = fEpsilon;
537 if (pb < fEpsilon) pb = fEpsilon;
538 Double_t r = ps/(ps + pb);
539 if (r >= 1.0) r = 1. - 1.e-15;
540
541 if (fTransformLikelihoodOutput) {
542 // inverse Fermi function
543
544 // sanity check
545 if (r <= 0.0) r = fEpsilon;
546 else if (r >= 1.0) r = 1. - 1.e-15;
547
548 Double_t tau = 15.0;
549 r = - TMath::Log(1.0/r - 1.0)/tau;
550 }
551
552 return r;
553}
554
555////////////////////////////////////////////////////////////////////////////////
556/// write options to stream
557
558void TMVA::MethodLikelihood::WriteOptionsToStream( std::ostream& o, const TString& prefix ) const
559{
561
562 // writing the options defined for the different pdfs
563 if (fDefaultPDFLik != 0) {
564 o << prefix << std::endl << prefix << "#Default Likelihood PDF Options:" << std::endl << prefix << std::endl;
565 fDefaultPDFLik->WriteOptionsToStream( o, prefix );
566 }
567 for (UInt_t ivar = 0; ivar < fPDFSig->size(); ivar++) {
568 if ((*fPDFSig)[ivar] != 0) {
569 o << prefix << std::endl << prefix << TString::Format("#Signal[%d] Likelihood PDF Options:",ivar) << std::endl << prefix << std::endl;
570 (*fPDFSig)[ivar]->WriteOptionsToStream( o, prefix );
571 }
572 if ((*fPDFBgd)[ivar] != 0) {
573 o << prefix << std::endl << prefix << "#Background[%d] Likelihood PDF Options:" << std::endl << prefix << std::endl;
574 (*fPDFBgd)[ivar]->WriteOptionsToStream( o, prefix );
575 }
576 }
577}
578
579////////////////////////////////////////////////////////////////////////////////
580/// write weights to XML
581
583{
584 void* wght = gTools().AddChild(parent, "Weights");
585 gTools().AddAttr(wght, "NVariables", GetNvar());
586 gTools().AddAttr(wght, "NClasses", 2);
587 void* pdfwrap;
588 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
589 if ( (*fPDFSig)[ivar]==0 || (*fPDFBgd)[ivar]==0 )
590 Log() << kFATAL << "Reference histograms for variable " << ivar
591 << " don't exist, can't write it to weight file" << Endl;
592 pdfwrap = gTools().AddChild(wght, "PDFDescriptor");
593 gTools().AddAttr(pdfwrap, "VarIndex", ivar);
594 gTools().AddAttr(pdfwrap, "ClassIndex", 0);
595 (*fPDFSig)[ivar]->AddXMLTo(pdfwrap);
596 pdfwrap = gTools().AddChild(wght, "PDFDescriptor");
597 gTools().AddAttr(pdfwrap, "VarIndex", ivar);
598 gTools().AddAttr(pdfwrap, "ClassIndex", 1);
599 (*fPDFBgd)[ivar]->AddXMLTo(pdfwrap);
600 }
601}
602
603////////////////////////////////////////////////////////////////////////////////
604/// computes ranking of input variables
605
607{
608 // create the ranking object
609 if (fRanking) delete fRanking;
610 fRanking = new Ranking( GetName(), "Delta Separation" );
611
612 Double_t sepRef = -1, sep = -1;
613 for (Int_t ivar=-1; ivar<(Int_t)GetNvar(); ivar++) {
614
615 // this variable should not be used
616 fDropVariable = ivar;
617
618 TString nameS = TString::Format( "rS_%i", ivar+1 );
619 TString nameB = TString::Format( "rB_%i", ivar+1 );
620 TH1* rS = new TH1F( nameS, nameS, 80, 0, 1 );
621 TH1* rB = new TH1F( nameB, nameB, 80, 0, 1 );
622
623 // the event loop
624 for (Int_t ievt=0; ievt<Data()->GetNTrainingEvents(); ievt++) {
625
626 const Event* origEv = Data()->GetEvent(ievt);
627 GetTransformationHandler().SetTransformationReferenceClass( origEv->GetClass() );
628 const Event* ev = GetTransformationHandler().Transform(Data()->GetEvent(ievt));
629
630 Double_t lk = this->GetMvaValue();
631 Double_t w = ev->GetWeight();
632 if (DataInfo().IsSignal(ev)) rS->Fill( lk, w );
633 else rB->Fill( lk, w );
634 }
635
636 // compute separation
637 sep = TMVA::gTools().GetSeparation( rS, rB );
638 if (ivar == -1) sepRef = sep;
639 sep = sepRef - sep;
640
641 // don't need these histograms anymore
642 delete rS;
643 delete rB;
644
645 if (ivar >= 0) fRanking->AddRank( Rank( DataInfo().GetVariableInfo(ivar).GetInternalName(), sep ) );
646 }
647
648 fDropVariable = -1;
649
650 return fRanking;
651}
652
653////////////////////////////////////////////////////////////////////////////////
654/// write reference PDFs to ROOT file
655
657{
658 TString pname = "PDF_";
659 for (UInt_t ivar=0; ivar<GetNvar(); ivar++){
660 (*fPDFSig)[ivar]->Write( pname + GetInputVar( ivar ) + "_S" );
661 (*fPDFBgd)[ivar]->Write( pname + GetInputVar( ivar ) + "_B" );
662 }
663}
664
665////////////////////////////////////////////////////////////////////////////////
666/// read weights from XML
667
669{
670 TString pname = "PDF_";
672 TH1::AddDirectory(0); // this avoids the binding of the hists in TMVA::PDF to the current ROOT file
673 UInt_t nvars=0;
674 gTools().ReadAttr(wghtnode, "NVariables",nvars);
675 void* descnode = gTools().GetChild(wghtnode);
676 for (UInt_t ivar=0; ivar<nvars; ivar++){
677 void* pdfnode = gTools().GetChild(descnode);
678 Log() << kDEBUG << "Reading signal and background PDF for variable: " << GetInputVar( ivar ) << Endl;
679 if ((*fPDFSig)[ivar] !=0) delete (*fPDFSig)[ivar];
680 if ((*fPDFBgd)[ivar] !=0) delete (*fPDFBgd)[ivar];
681 (*fPDFSig)[ivar] = new PDF( GetInputVar( ivar ) + " PDF Sig" );
682 (*fPDFBgd)[ivar] = new PDF( GetInputVar( ivar ) + " PDF Bkg" );
683 (*fPDFSig)[ivar]->SetReadingVersion( GetTrainingTMVAVersionCode() );
684 (*fPDFBgd)[ivar]->SetReadingVersion( GetTrainingTMVAVersionCode() );
685 (*(*fPDFSig)[ivar]).ReadXML(pdfnode);
688 (*(*fPDFBgd)[ivar]).ReadXML(pdfnode);
690 }
692}
693
694////////////////////////////////////////////////////////////////////////////////
695/// read weight info from file
696/// nothing to do for this method
697
699{
700 TString pname = "PDF_";
702 TH1::AddDirectory(0); // this avoids the binding of the hists in TMVA::PDF to the current ROOT file
703 for (UInt_t ivar=0; ivar<GetNvar(); ivar++){
704 Log() << kDEBUG << "Reading signal and background PDF for variable: " << GetInputVar( ivar ) << Endl;
705 if ((*fPDFSig)[ivar] !=0) delete (*fPDFSig)[ivar];
706 if ((*fPDFBgd)[ivar] !=0) delete (*fPDFBgd)[ivar];
707 (*fPDFSig)[ivar] = new PDF(GetInputVar( ivar ) + " PDF Sig" );
708 (*fPDFBgd)[ivar] = new PDF(GetInputVar( ivar ) + " PDF Bkg");
709 (*fPDFSig)[ivar]->SetReadingVersion( GetTrainingTMVAVersionCode() );
710 (*fPDFBgd)[ivar]->SetReadingVersion( GetTrainingTMVAVersionCode() );
711 istr >> *(*fPDFSig)[ivar];
712 istr >> *(*fPDFBgd)[ivar];
713 }
715}
716
717////////////////////////////////////////////////////////////////////////////////
718/// read reference PDF from ROOT file
719
721{
722 TString pname = "PDF_";
724 TH1::AddDirectory(0); // this avoids the binding of the hists in TMVA::PDF to the current ROOT file
725 for (UInt_t ivar=0; ivar<GetNvar(); ivar++){
726 (*fPDFSig)[ivar] = (TMVA::PDF*)rf.Get( TString::Format( "PDF_%s_S", GetInputVar( ivar ).Data() ) );
727 (*fPDFBgd)[ivar] = (TMVA::PDF*)rf.Get( TString::Format( "PDF_%s_B", GetInputVar( ivar ).Data() ) );
728 }
730}
731
732////////////////////////////////////////////////////////////////////////////////
733/// write histograms and PDFs to file for monitoring purposes
734
736{
737 Log() << kINFO << "Write monitoring histograms to file: " << BaseDir()->GetPath() << Endl;
738 BaseDir()->cd();
739
740 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
741 (*fHistSig)[ivar]->Write();
742 (*fHistBgd)[ivar]->Write();
743 if ((*fHistSig_smooth)[ivar] != 0) (*fHistSig_smooth)[ivar]->Write();
744 if ((*fHistBgd_smooth)[ivar] != 0) (*fHistBgd_smooth)[ivar]->Write();
745 (*fPDFSig)[ivar]->GetPDFHist()->Write();
746 (*fPDFBgd)[ivar]->GetPDFHist()->Write();
747
748 if ((*fPDFSig)[ivar]->GetNSmoothHist() != 0) (*fPDFSig)[ivar]->GetNSmoothHist()->Write();
749 if ((*fPDFBgd)[ivar]->GetNSmoothHist() != 0) (*fPDFBgd)[ivar]->GetNSmoothHist()->Write();
750
751 // add special plots to check the smoothing in the GetVal method
752 Float_t xmin=((*fPDFSig)[ivar]->GetPDFHist()->GetXaxis())->GetXmin();
753 Float_t xmax=((*fPDFSig)[ivar]->GetPDFHist()->GetXaxis())->GetXmax();
754 TH1F* mm = new TH1F( (*fInputVars)[ivar]+"_additional_check",
755 (*fInputVars)[ivar]+"_additional_check", 15000, xmin, xmax );
756 Double_t intBin = (xmax-xmin)/15000;
757 for (Int_t bin=0; bin < 15000; bin++) {
758 Double_t x = (bin + 0.5)*intBin + xmin;
759 mm->SetBinContent(bin+1 ,(*fPDFSig)[ivar]->GetVal(x));
760 }
761 mm->Write();
762
763 // ---------- create cloned low-binned histogram for comparison in macros (mainly necessary for KDE)
764 TH1* h[2] = { (*fHistSig)[ivar], (*fHistBgd)[ivar] };
765 for (UInt_t i=0; i<2; i++) {
766 TH1* hclone = (TH1F*)h[i]->Clone( TString(h[i]->GetName()) + "_nice" );
767 hclone->SetName ( TString(h[i]->GetName()) + "_nice" );
768 hclone->SetTitle( TString(h[i]->GetTitle()) + "" );
769 if (hclone->GetNbinsX() > 100) {
770 Int_t resFactor = 5;
771 hclone->Rebin( resFactor );
772 hclone->Scale( 1.0/resFactor );
773 }
774 hclone->Write();
775 }
776 // ----------
777 }
778}
779
780////////////////////////////////////////////////////////////////////////////////
781/// write specific header of the classifier (mostly include files)
782
784{
785 fout << "#include <cmath>" << std::endl;
786 fout << "#include <cstdlib>" << std::endl;
787}
788
789////////////////////////////////////////////////////////////////////////////////
790/// write specific classifier response
791
792void TMVA::MethodLikelihood::MakeClassSpecific( std::ostream& fout, const TString& className ) const
793{
794 Int_t dp = fout.precision();
795 fout << " double fEpsilon;" << std::endl;
796
797 Int_t * nbin = new Int_t[GetNvar()];
798
799 Int_t nbinMax=-1;
800 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
801 nbin[ivar]=(*fPDFSig)[ivar]->GetPDFHist()->GetNbinsX();
802 if (nbin[ivar] > nbinMax) nbinMax=nbin[ivar];
803 }
804
805 fout << " static float fRefS[][" << nbinMax << "]; "
806 << "// signal reference vector [nvars][max_nbins]" << std::endl;
807 fout << " static float fRefB[][" << nbinMax << "]; "
808 << "// backgr reference vector [nvars][max_nbins]" << std::endl << std::endl;
809 fout << "// if a variable has its PDF encoded as a spline0 --> treat it like an Integer valued one" <<std::endl;
810 fout << " bool fHasDiscretPDF[" << GetNvar() <<"]; "<< std::endl;
811 fout << " int fNbin[" << GetNvar() << "]; "
812 << "// number of bins (discrete variables may have less bins)" << std::endl;
813 fout << " double fHistMin[" << GetNvar() << "]; " << std::endl;
814 fout << " double fHistMax[" << GetNvar() << "]; " << std::endl;
815
816 fout << " double TransformLikelihoodOutput( double, double ) const;" << std::endl;
817 fout << "};" << std::endl;
818 fout << "" << std::endl;
819 fout << "inline void " << className << "::Initialize() " << std::endl;
820 fout << "{" << std::endl;
821 fout << " fEpsilon = " << fEpsilon << ";" << std::endl;
822 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
823 fout << " fNbin[" << ivar << "] = " << (*fPDFSig)[ivar]->GetPDFHist()->GetNbinsX() << ";" << std::endl;
824 fout << " fHistMin[" << ivar << "] = " << (*fPDFSig)[ivar]->GetPDFHist()->GetXaxis()->GetXmin() << ";" << std::endl;
825 fout << " fHistMax[" << ivar << "] = " << (*fPDFSig)[ivar]->GetPDFHist()->GetXaxis()->GetXmax() << ";" << std::endl;
826 // sanity check (for previous code lines)
827 if ((((*fPDFSig)[ivar]->GetPDFHist()->GetNbinsX() != nbin[ivar] ||
828 (*fPDFBgd)[ivar]->GetPDFHist()->GetNbinsX() != nbin[ivar])
829 ) ||
830 (*fPDFSig)[ivar]->GetPDFHist()->GetNbinsX() != (*fPDFBgd)[ivar]->GetPDFHist()->GetNbinsX()) {
831 Log() << kFATAL << "<MakeClassSpecific> Mismatch in binning of variable "
832 << "\"" << GetOriginalVarName(ivar) << "\" of type: \'" << DataInfo().GetVariableInfo(ivar).GetVarType()
833 << "\' : "
834 << "nxS = " << (*fPDFSig)[ivar]->GetPDFHist()->GetNbinsX() << ", "
835 << "nxB = " << (*fPDFBgd)[ivar]->GetPDFHist()->GetNbinsX()
836 << " while we expect " << nbin[ivar]
837 << Endl;
838 }
839 }
840 for (UInt_t ivar=0; ivar<GetNvar(); ivar++){
841 if ((*fPDFSig)[ivar]->GetInterpolMethod() == TMVA::PDF::kSpline0)
842 fout << " fHasDiscretPDF[" << ivar <<"] = true; " << std::endl;
843 else
844 fout << " fHasDiscretPDF[" << ivar <<"] = false; " << std::endl;
845 }
846
847 fout << "}" << std::endl << std::endl;
848
849 fout << "inline double " << className
850 << "::GetMvaValue__( const std::vector<double>& inputValues ) const" << std::endl;
851 fout << "{" << std::endl;
852 fout << " double ps(1), pb(1);" << std::endl;
853 fout << " std::vector<double> inputValuesSig = inputValues;" << std::endl;
854 fout << " std::vector<double> inputValuesBgd = inputValues;" << std::endl;
855 if (GetTransformationHandler().GetTransformationList().GetSize() != 0) {
856 fout << " Transform(inputValuesSig,0);" << std::endl;
857 fout << " Transform(inputValuesBgd,1);" << std::endl;
858 }
859 fout << " for (size_t ivar = 0; ivar < GetNvar(); ivar++) {" << std::endl;
860 fout << std::endl;
861 fout << " // dummy at present... will be used for variable transforms" << std::endl;
862 fout << " double x[2] = { inputValuesSig[ivar], inputValuesBgd[ivar] };" << std::endl;
863 fout << std::endl;
864 fout << " for (int itype=0; itype < 2; itype++) {" << std::endl;
865 fout << std::endl;
866 fout << " // interpolate linearly between adjacent bins" << std::endl;
867 fout << " // this is not useful for discrete variables (or forced Spline0)" << std::endl;
868 fout << " int bin = int((x[itype] - fHistMin[ivar])/(fHistMax[ivar] - fHistMin[ivar])*fNbin[ivar]) + 0;" << std::endl;
869 fout << std::endl;
870 fout << " // since the test data sample is in general different from the training sample" << std::endl;
871 fout << " // it can happen that the min/max of the training sample are trespassed --> correct this" << std::endl;
872 fout << " if (bin < 0) {" << std::endl;
873 fout << " bin = 0;" << std::endl;
874 fout << " x[itype] = fHistMin[ivar];" << std::endl;
875 fout << " }" << std::endl;
876 fout << " else if (bin >= fNbin[ivar]) {" << std::endl;
877 fout << " bin = fNbin[ivar]-1;" << std::endl;
878 fout << " x[itype] = fHistMax[ivar];" << std::endl;
879 fout << " }" << std::endl;
880 fout << std::endl;
881 fout << " // find corresponding histogram from cached indices" << std::endl;
882 fout << " float ref = (itype == 0) ? fRefS[ivar][bin] : fRefB[ivar][bin];" << std::endl;
883 fout << std::endl;
884 fout << " // sanity check" << std::endl;
885 fout << " if (ref < 0) {" << std::endl;
886 fout << " std::cout << \"Fatal error in " << className
887 << ": bin entry < 0 ==> abort\" << std::endl;" << std::endl;
888 fout << " std::exit(1);" << std::endl;
889 fout << " }" << std::endl;
890 fout << std::endl;
891 fout << " double p = ref;" << std::endl;
892 fout << std::endl;
893 fout << " if (GetType(ivar) != 'I' && !fHasDiscretPDF[ivar]) {" << std::endl;
894 fout << " float bincenter = (bin + 0.5)/fNbin[ivar]*(fHistMax[ivar] - fHistMin[ivar]) + fHistMin[ivar];" << std::endl;
895 fout << " int nextbin = bin;" << std::endl;
896 fout << " if ((x[itype] > bincenter && bin != fNbin[ivar]-1) || bin == 0) " << std::endl;
897 fout << " nextbin++;" << std::endl;
898 fout << " else" << std::endl;
899 fout << " nextbin--; " << std::endl;
900 fout << std::endl;
901 fout << " double refnext = (itype == 0) ? fRefS[ivar][nextbin] : fRefB[ivar][nextbin];" << std::endl;
902 fout << " float nextbincenter = (nextbin + 0.5)/fNbin[ivar]*(fHistMax[ivar] - fHistMin[ivar]) + fHistMin[ivar];" << std::endl;
903 fout << std::endl;
904 fout << " double dx = bincenter - nextbincenter;" << std::endl;
905 fout << " double dy = ref - refnext;" << std::endl;
906 fout << " p += (x[itype] - bincenter) * dy/dx;" << std::endl;
907 fout << " }" << std::endl;
908 fout << std::endl;
909 fout << " if (p < fEpsilon) p = fEpsilon; // avoid zero response" << std::endl;
910 fout << std::endl;
911 fout << " if (itype == 0) ps *= p;" << std::endl;
912 fout << " else pb *= p;" << std::endl;
913 fout << " } " << std::endl;
914 fout << " } " << std::endl;
915 fout << std::endl;
916 fout << " // the likelihood ratio (transform it ?)" << std::endl;
917 fout << " return TransformLikelihoodOutput( ps, pb ); " << std::endl;
918 fout << "}" << std::endl << std::endl;
919
920 fout << "inline double " << className << "::TransformLikelihoodOutput( double ps, double pb ) const" << std::endl;
921 fout << "{" << std::endl;
922 fout << " // returns transformed or non-transformed output" << std::endl;
923 fout << " if (ps < fEpsilon) ps = fEpsilon;" << std::endl;
924 fout << " if (pb < fEpsilon) pb = fEpsilon;" << std::endl;
925 fout << " double r = ps/(ps + pb);" << std::endl;
926 fout << " if (r >= 1.0) r = 1. - 1.e-15;" << std::endl;
927 fout << std::endl;
928 fout << " if (" << (fTransformLikelihoodOutput ? "true" : "false") << ") {" << std::endl;
929 fout << " // inverse Fermi function" << std::endl;
930 fout << std::endl;
931 fout << " // sanity check" << std::endl;
932 fout << " if (r <= 0.0) r = fEpsilon;" << std::endl;
933 fout << " else if (r >= 1.0) r = 1. - 1.e-15;" << std::endl;
934 fout << std::endl;
935 fout << " double tau = 15.0;" << std::endl;
936 fout << " r = - log(1.0/r - 1.0)/tau;" << std::endl;
937 fout << " }" << std::endl;
938 fout << std::endl;
939 fout << " return r;" << std::endl;
940 fout << "}" << std::endl;
941 fout << std::endl;
942
943 fout << "// Clean up" << std::endl;
944 fout << "inline void " << className << "::Clear() " << std::endl;
945 fout << "{" << std::endl;
946 fout << " // nothing to clear" << std::endl;
947 fout << "}" << std::endl << std::endl;
948
949 fout << "// signal map" << std::endl;
950 fout << "float " << className << "::fRefS[][" << nbinMax << "] = " << std::endl;
951 fout << "{ " << std::endl;
952 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
953 fout << " { ";
954 for (Int_t ibin=1; ibin<=nbinMax; ibin++) {
955 if (ibin-1 < nbin[ivar])
956 fout << (*fPDFSig)[ivar]->GetPDFHist()->GetBinContent(ibin);
957 else
958 fout << -1;
959
960 if (ibin < nbinMax) fout << ", ";
961 }
962 fout << " }, " << std::endl;
963 }
964 fout << "}; " << std::endl;
965 fout << std::endl;
966
967 fout << "// background map" << std::endl;
968 fout << "float " << className << "::fRefB[][" << nbinMax << "] = " << std::endl;
969 fout << "{ " << std::endl;
970 for (UInt_t ivar=0; ivar<GetNvar(); ivar++) {
971 fout << " { ";
972 fout << std::setprecision(8);
973 for (Int_t ibin=1; ibin<=nbinMax; ibin++) {
974 if (ibin-1 < nbin[ivar])
975 fout << (*fPDFBgd)[ivar]->GetPDFHist()->GetBinContent(ibin);
976 else
977 fout << -1;
978
979 if (ibin < nbinMax) fout << ", ";
980 }
981 fout << " }, " << std::endl;
982 }
983 fout << "}; " << std::endl;
984 fout << std::endl;
985 fout << std::setprecision(dp);
986
987 delete[] nbin;
988}
989
990////////////////////////////////////////////////////////////////////////////////
991/// get help message text
992///
993/// typical length of text line:
994/// "|--------------------------------------------------------------|"
995
997{
998 Log() << Endl;
999 Log() << gTools().Color("bold") << "--- Short description:" << gTools().Color("reset") << Endl;
1000 Log() << Endl;
1001 Log() << "The maximum-likelihood classifier models the data with probability " << Endl;
1002 Log() << "density functions (PDF) reproducing the signal and background" << Endl;
1003 Log() << "distributions of the input variables. Correlations among the " << Endl;
1004 Log() << "variables are ignored." << Endl;
1005 Log() << Endl;
1006 Log() << gTools().Color("bold") << "--- Performance optimisation:" << gTools().Color("reset") << Endl;
1007 Log() << Endl;
1008 Log() << "Required for good performance are decorrelated input variables" << Endl;
1009 Log() << "(PCA transformation via the option \"VarTransform=Decorrelate\"" << Endl;
1010 Log() << "may be tried). Irreducible non-linear correlations may be reduced" << Endl;
1011 Log() << "by precombining strongly correlated input variables, or by simply" << Endl;
1012 Log() << "removing one of the variables." << Endl;
1013 Log() << Endl;
1014 Log() << gTools().Color("bold") << "--- Performance tuning via configuration options:" << gTools().Color("reset") << Endl;
1015 Log() << Endl;
1016 Log() << "High fidelity PDF estimates are mandatory, i.e., sufficient training " << Endl;
1017 Log() << "statistics is required to populate the tails of the distributions" << Endl;
1018 Log() << "It would be a surprise if the default Spline or KDE kernel parameters" << Endl;
1019 Log() << "provide a satisfying fit to the data. The user is advised to properly" << Endl;
1020 Log() << "tune the events per bin and smooth options in the spline cases" << Endl;
1021 Log() << "individually per variable. If the KDE kernel is used, the adaptive" << Endl;
1022 Log() << "Gaussian kernel may lead to artefacts, so please always also try" << Endl;
1023 Log() << "the non-adaptive one." << Endl;
1024 Log() << "" << Endl;
1025 Log() << "All tuning parameters must be adjusted individually for each input" << Endl;
1026 Log() << "variable!" << Endl;
1027}
1028
#define REGISTER_METHOD(CLASS)
for example
#define h(i)
Definition RSha256.hxx:106
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
float xmin
float xmax
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:131
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:879
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:109
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition TH1.cxx:9166
static void AddDirectory(Bool_t add=kTRUE)
Sets the flag controlling the automatic add of histograms in memory.
Definition TH1.cxx:1263
virtual Int_t GetNbinsX() const
Definition TH1.h:542
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:9247
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition TH1.cxx:5063
virtual Int_t FindBin(Double_t x, Double_t y=0, Double_t z=0)
Return Global bin number corresponding to x,y,z.
Definition TH1.cxx:3649
static Bool_t AddDirectoryStatus()
Static function: cannot be inlined on Windows/NT.
Definition TH1.cxx:741
void WriteOptionsToStream(std::ostream &o, const TString &prefix) const
write options to output stream (e.g. in writing the MVA weight files
Class that contains all the data information.
Definition DataSetInfo.h:62
Virtual base Class for all MVA method.
Definition MethodBase.h:111
virtual void DeclareCompatibilityOptions()
options that are used ONLY for the READER to ensure backward compatibility they are hence without any...
Likelihood analysis ("non-parametric approach")
Bool_t HasAnalysisType(Types::EAnalysisType type, UInt_t numberClasses, UInt_t numberTargets) override
FDA can handle classification with 2 classes.
void ProcessOptions() override
process user options reference cut value to distinguish signal-like from background-like events
virtual void WriteOptionsToStream(std::ostream &o, const TString &prefix) const
write options to stream
void MakeClassSpecific(std::ostream &, const TString &) const override
write specific classifier response
const Ranking * CreateRanking() override
computes ranking of input variables
void ReadWeightsFromStream(std::istream &istr) override
read weight info from file nothing to do for this method
void Train() override
create reference distributions (PDFs) from signal and background events: fill histograms and smooth t...
void GetHelpMessage() const override
get help message text
virtual ~MethodLikelihood()
destructor
void DeclareOptions() override
define the options (their key words) that can be set in the option string
void ReadWeightsFromXML(void *wghtnode) override
read weights from XML
MethodLikelihood(const TString &jobName, const TString &methodTitle, DataSetInfo &theData, const TString &theOption="")
standard constructor
void WriteMonitoringHistosToFile() const override
write histograms and PDFs to file for monitoring purposes
void Init() override
default initialisation called by all constructors
void DeclareCompatibilityOptions() override
options that are used ONLY for the READER to ensure backward compatibility they are hence without any...
void AddWeightsXMLTo(void *parent) const override
write weights to XML
Double_t TransformLikelihoodOutput(Double_t ps, Double_t pb) const
returns transformed or non-transformed output
Double_t GetMvaValue(Double_t *err=nullptr, Double_t *errUpper=nullptr) override
returns the likelihood estimator for signal fill a new Likelihood branch into the testTree
void WriteWeightsToStream(TFile &rf) const
write reference PDFs to ROOT file
void MakeClassSpecificHeader(std::ostream &, const TString &="") const override
write specific header of the classifier (mostly include files)
PDF wrapper for histograms; uses user-defined spline interpolation.
Definition PDF.h:63
TH1 * GetPDFHist() const
Definition PDF.h:92
@ kSpline0
Definition PDF.h:70
Ranking for variables in method (implementation)
Definition Ranking.h:48
Double_t GetSeparation(TH1 *S, TH1 *B) const
compute "separation" defined as
Definition Tools.cxx:121
const TString & Color(const TString &)
human readable color strings
Definition Tools.cxx:828
void ReadAttr(void *node, const char *, T &value)
read attribute from xml
Definition Tools.h:329
void * GetChild(void *parent, const char *childname=nullptr)
get child node
Definition Tools.cxx:1150
void AddAttr(void *node, const char *, const T &value, Int_t precision=16)
add attribute to xml
Definition Tools.h:347
void * AddChild(void *parent, const char *childname, const char *content=nullptr, bool isRootNode=false)
add child node
Definition Tools.cxx:1124
void * GetNextChild(void *prevchild, const char *childname=nullptr)
XML helpers.
Definition Tools.cxx:1162
Singleton class for Global types used by TMVA.
Definition Types.h:71
@ kClassification
Definition Types.h:127
virtual Int_t Write(const char *name=nullptr, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
Definition TObject.cxx:964
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2384
Double_t x[n]
Definition legend1.C:17
create variable transformations
Tools & gTools()
MsgLogger & Endl(MsgLogger &ml)
Definition MsgLogger.h:148
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition TMath.h:704
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition TMath.h:767
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:199