Logo ROOT   6.08/07
Reference Guide
ResultsRegression.cxx
Go to the documentation of this file.
1 // @(#)root/tmva $Id$
2 // Author: Andreas Hoecker, Peter Speckmayer, Joerg Stelzer, Helge Voss
3 
4 /**********************************************************************************
5  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6  * Package: TMVA *
7  * Class : ResultsRegression *
8  * Web : http://tmva.sourceforge.net *
9  * *
10  * Description: *
11  * Implementation (see header for description) *
12  * *
13  * Authors (alphabetical): *
14  * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15  * Peter Speckmayer <Peter.Speckmayer@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  * *
19  * Copyright (c) 2006: *
20  * CERN, Switzerland *
21  * MPI-K Heidelberg, Germany *
22  * *
23  * Redistribution and use in source and binary forms, with or without *
24  * modification, are permitted according to the terms listed in LICENSE *
25  * (http://tmva.sourceforge.net/LICENSE) *
26  **********************************************************************************/
27 
28 #include "TMVA/ResultsRegression.h"
29 
30 #include "TMVA/DataSet.h"
31 #include "TMVA/DataSetInfo.h"
32 #include "TMVA/MsgLogger.h"
33 #include "TMVA/Results.h"
34 #include "TMVA/Types.h"
35 #include "TMVA/VariableInfo.h"
36 
37 #include "TH1F.h"
38 #include "TH2F.h"
39 #include "TString.h"
40 
41 #include <vector>
42 
43 
44 ////////////////////////////////////////////////////////////////////////////////
45 /// constructor
46 
48  : Results( dsi, resultsName ),
49  fLogger( new MsgLogger(Form("ResultsRegression%s",resultsName.Data()) , kINFO) )
50 {
51 }
52 
53 ////////////////////////////////////////////////////////////////////////////////
54 /// destructor
55 
57 {
58  delete fLogger;
59 }
60 
61 ////////////////////////////////////////////////////////////////////////////////
62 
63 void TMVA::ResultsRegression::SetValue( std::vector<Float_t>& value, Int_t ievt )
64 {
65  if (ievt >= (Int_t)fRegValues.size()) fRegValues.resize( ievt+1 );
66  fRegValues[ievt] = value;
67 }
68 
69 ////////////////////////////////////////////////////////////////////////////////
70 
72 {
73  DataSet* ds = GetDataSet();
74  ds->SetCurrentType( GetTreeType() );
75  const DataSetInfo* dsi = GetDataSetInfo();
76  TString name( Form("tgt_%d",tgtNum) );
77  VariableInfo vinf = dsi->GetTargetInfo(tgtNum);
78  Float_t xmin=0., xmax=0.;
79  if (truncate){
80  xmax = truncvalue;
81  }
82  else{
83  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
84  const Event* ev = ds->GetEvent(ievt);
85  std::vector<Float_t> regVal = fRegValues.at(ievt);
86  Float_t val = regVal.at( tgtNum ) - ev->GetTarget( tgtNum );
87  val *= val;
88  xmax = val> xmax? val: xmax;
89  }
90  }
91  xmax *= 1.1;
92  Int_t nbins = 500;
93  TH1F* h = new TH1F( name, name, nbins, xmin, xmax);
94  h->SetDirectory(0);
95  h->GetXaxis()->SetTitle("Quadratic Deviation");
96  h->GetYaxis()->SetTitle("Weighted Entries");
97 
98  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
99  const Event* ev = ds->GetEvent(ievt);
100  std::vector<Float_t> regVal = fRegValues.at(ievt);
101  Float_t val = regVal.at( tgtNum ) - ev->GetTarget( tgtNum );
102  val *= val;
103  Float_t weight = ev->GetWeight();
104  if (!truncate || val<=truncvalue ) h->Fill( val, weight);
105  }
106  return h;
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 
112 {
113  DataSet* ds = GetDataSet();
114  ds->SetCurrentType( GetTreeType() );
115 
116  TString name( Form("tgt_%d_var_%d",tgtNum, varNum) );
117  const DataSetInfo* dsi = GetDataSetInfo();
118  Float_t xmin, xmax;
119  Bool_t takeTargets = kFALSE;
120  if (varNum >= dsi->GetNVariables()) {
121  takeTargets = kTRUE;
122  varNum -= dsi->GetNVariables();
123  }
124  if (!takeTargets) {
125  VariableInfo vinf = dsi->GetVariableInfo(varNum);
126  xmin = vinf.GetMin();
127  xmax = vinf.GetMax();
128 
129  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
130  const Event* ev = ds->GetEvent(ievt);
131  Float_t val = ev->GetValue(varNum);
132 
133  if (val < xmin ) xmin = val;
134  if (val > xmax ) xmax = val;
135  }
136 
137  }
138  else {
139  VariableInfo vinf = dsi->GetTargetInfo(varNum);
140  xmin = vinf.GetMin();
141  xmax = vinf.GetMax();
142 
143  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
144  const Event* ev = ds->GetEvent(ievt);
145  Float_t val = ev->GetTarget(varNum);
146 
147  if (val < xmin ) xmin = val;
148  if (val > xmax ) xmax = val;
149  }
150  }
151 
152  Float_t ymin = FLT_MAX;
153  Float_t ymax = -FLT_MAX;
154 
155  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
156  const Event* ev = ds->GetEvent(ievt);
157  std::vector<Float_t> regVal = fRegValues.at(ievt);
158 
159  Float_t diff = regVal.at( tgtNum ) - ev->GetTarget( tgtNum );
160  if (diff < ymin ) ymin = diff;
161  if (diff > ymax ) ymax = diff;
162  }
163 
164  Int_t nxbins = 50;
165  Int_t nybins = 50;
166 
167  Float_t epsilon = TMath::Abs(xmax-xmin)/((Float_t)nxbins-1);
168  xmin -= 1.01*epsilon;
169  xmax += 1.01*epsilon;
170 
171  epsilon = (ymax-ymin)/((Float_t)nybins-1);
172  ymin -= 1.01*epsilon;
173  ymax += 1.01*epsilon;
174 
175 
176  TH2F* h = new TH2F( name, name, nxbins, xmin, xmax, nybins, ymin, ymax );
177  h->SetDirectory(0);
178 
179  h->GetXaxis()->SetTitle( (takeTargets ? dsi->GetTargetInfo(varNum).GetTitle() : dsi->GetVariableInfo(varNum).GetTitle() ) );
180  TString varName( dsi->GetTargetInfo(tgtNum).GetTitle() );
181  TString yName( varName+TString("_{regression} - ") + varName+TString("_{true}") );
182  h->GetYaxis()->SetTitle( yName );
183 
184  for (Int_t ievt=0; ievt<ds->GetNEvents(); ievt++) {
185  const Event* ev = ds->GetEvent(ievt);
186  std::vector<Float_t> regVal = fRegValues.at(ievt);
187 
188  Float_t xVal = (takeTargets?ev->GetTarget( varNum ):ev->GetValue( varNum ));
189  Float_t yVal = regVal.at( tgtNum ) - ev->GetTarget( tgtNum );
190 
191  h->Fill( xVal, yVal );
192  }
193 
194  return h;
195 }
196 
197 ////////////////////////////////////////////////////////////////////////////////
198 
200 {
201  Log() << kINFO << "Create variable histograms" << Endl;
202  const DataSetInfo* dsi = GetDataSetInfo();
203 
204  for (UInt_t ivar = 0; ivar < dsi->GetNVariables(); ivar++) {
205  for (UInt_t itgt = 0; itgt < dsi->GetNTargets(); itgt++) {
206  TH2F* h = DeviationAsAFunctionOf( ivar, itgt );
207  TString name( Form("%s_reg_var%d_rtgt%d",prefix.Data(),ivar,itgt) );
208  h->SetName( name );
209  h->SetTitle( name );
210  Store( h );
211  }
212  }
213  Log() << kINFO << "Create regression target histograms" << Endl;
214  for (UInt_t ivar = 0; ivar < dsi->GetNTargets(); ivar++) {
215  for (UInt_t itgt = 0; itgt < dsi->GetNTargets(); itgt++) {
216  TH2F* h = DeviationAsAFunctionOf( dsi->GetNVariables()+ivar, itgt );
217  TString name( Form("%s_reg_tgt%d_rtgt%d",prefix.Data(),ivar,itgt) );
218  h->SetName( name );
219  h->SetTitle( name );
220  Store( h );
221  }
222  }
223 
224  Log() << kINFO << "Create regression average deviation" << Endl;
225  for (UInt_t itgt = 0; itgt < dsi->GetNTargets(); itgt++) {
226  TH1F* h = QuadraticDeviation(itgt);
227  TString name( Form("%s_Quadr_Deviation_target_%d_",prefix.Data(),itgt) );
228  h->SetName( name );
229  h->SetTitle( name );
230  Double_t yq[1], xq[]={0.9};
231  h->GetQuantiles(1,yq,xq);
232  Store( h );
233 
234  TH1F* htrunc = QuadraticDeviation(itgt, true, yq[0]);
235  TString name2( Form("%s_Quadr_Dev_best90perc_target_%d_",prefix.Data(),itgt) );
236  htrunc->SetName( name2 );
237  htrunc->SetTitle( name2 );
238  Store( htrunc );
239  }
240  Log() << kINFO << "Results created" << Endl;
241 }
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3125
UInt_t GetNVariables() const
Definition: DataSetInfo.h:128
std::vector< std::vector< Float_t > > fRegValues
float xmin
Definition: THbookFile.cxx:93
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:162
float Float_t
Definition: RtypesCore.h:53
virtual void SetDirectory(TDirectory *dir)
By default when an histogram is created, it is added to the list of histogram objects in the current ...
Definition: TH1.cxx:8051
Double_t GetMin() const
Definition: VariableInfo.h:71
float ymin
Definition: THbookFile.cxx:93
THist< 1, float, THistStatContent, THistStatUncertainty > TH1F
Definition: THist.hxx:302
TH1 * h
Definition: legend2.C:5
virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum=0)
Compute Quantiles for this histogram Quantile x_q of a probability distribution Function F is defined...
Definition: TH1.cxx:4190
const DataSetInfo * GetDataSetInfo() const
Definition: Results.h:77
DataSet * GetDataSet() const
Definition: Results.h:78
Basic string class.
Definition: TString.h:137
tomato 1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:575
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
const Bool_t kFALSE
Definition: Rtypes.h:92
TH2F * DeviationAsAFunctionOf(UInt_t varNum, UInt_t tgtNum)
int nbins[3]
void SetValue(std::vector< Float_t > &value, Int_t ievt)
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
std::vector< std::vector< double > > Data
Double_t GetWeight() const
return the event weight - depending on whether the flag IgnoreNegWeightsInTraining is or not...
Definition: Event.cxx:378
Double_t GetMax() const
Definition: VariableInfo.h:72
float ymax
Definition: THbookFile.cxx:93
Float_t GetTarget(UInt_t itgt) const
Definition: Event.h:104
UInt_t GetNTargets() const
Definition: DataSetInfo.h:129
tomato 2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:255
VariableInfo & GetTargetInfo(Int_t i)
Definition: DataSetInfo.h:119
unsigned int UInt_t
Definition: RtypesCore.h:42
char * Form(const char *fmt,...)
TAxis * GetYaxis()
Definition: TH1.h:325
float xmax
Definition: THbookFile.cxx:93
ResultsRegression(const DataSetInfo *dsi, TString resultsName)
constructor
MsgLogger & Log() const
message logger
REAL epsilon
Definition: triangle.c:617
Float_t GetValue(UInt_t ivar) const
return value of i&#39;th variable
Definition: Event.cxx:233
virtual void SetName(const char *name)
Change the name of this histogram.
Definition: TH1.cxx:8071
double Double_t
Definition: RtypesCore.h:55
void SetCurrentType(Types::ETreeType type) const
Definition: DataSet.h:114
VariableInfo & GetVariableInfo(Int_t i)
Definition: DataSetInfo.h:114
void CreateDeviationHistograms(TString prefix)
Types::ETreeType GetTreeType() const
Definition: Results.h:76
Long64_t GetNEvents(Types::ETreeType type=Types::kMaxTreeType) const
Definition: DataSet.h:229
virtual void SetTitle(const char *title)
Change (i.e.
Definition: TH1.cxx:6027
void Store(TObject *obj, const char *alias=0)
Definition: Results.cxx:83
const Bool_t kTRUE
Definition: Rtypes.h:91
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:292
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:155
THist< 2, float, THistStatContent, THistStatUncertainty > TH2F
Definition: THist.hxx:308
const Event * GetEvent() const
Definition: DataSet.cxx:211
char name[80]
Definition: TGX11.cxx:109
TH1F * QuadraticDeviation(UInt_t tgtNum, Bool_t truncate=false, Double_t truncvalue=0.)
TAxis * GetXaxis()
Definition: TH1.h:324
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:52
const char * Data() const
Definition: TString.h:349