Logo ROOT   6.10/09
Reference Guide
ConfidenceBelt.h
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 // Author: Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke
3 /*************************************************************************
4  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5  * All rights reserved. *
6  * *
7  * For the licensing terms see $ROOTSYS/LICENSE. *
8  * For the list of contributors see $ROOTSYS/README/CREDITS. *
9  *************************************************************************/
10 
11 #ifndef RooStats_ConfidenceBelt
12 #define RooStats_ConfidenceBelt
13 
14 #include "RooArgSet.h"
15 #include "RooAbsData.h"
16 #include "RooStats/ConfInterval.h"
17 
19 
20 #include "TRef.h"
21 
22 #include <vector>
23 #include <map>
24 
25 
26 namespace RooStats {
27 
28 ////////////////////////////////////////////////////////////////////////////////
29 
30  class SamplingSummaryLookup : public TObject {
31 
32  typedef std::pair<Double_t, Double_t> AcceptanceCriteria; // defined by Confidence level, leftside tail probability
33  typedef std::map<Int_t, AcceptanceCriteria> LookupTable; // map ( Index, ( CL, leftside tail prob) )
34 
35  public:
38 
39  void Add(Double_t cl, Double_t leftside){
40  // add cl,leftside pair to lookup table
41  AcceptanceCriteria tmp(cl, leftside);
42 
43  // should check to see if this is already in the map
44  if(GetLookupIndex(cl,leftside) >=0 ){
45  std::cout<< "SamplingSummaryLookup::Add, already in lookup table" << std::endl;
46  } else
47  fLookupTable[fLookupTable.size()]= tmp;
48  }
49 
51  // get index for cl,leftside pair
52  AcceptanceCriteria tmp(cl, leftside);
53 
54  Double_t tolerance = 1E-6; // some small number to protect floating point comparison. What is better way?
55  LookupTable::iterator it = fLookupTable.begin();
56  Int_t index = -1;
57  for(; it!= fLookupTable.end(); ++it) {
58  index++;
59  if( TMath::Abs( (*it).second.first - cl ) < tolerance &&
60  TMath::Abs( (*it).second.second - leftside ) < tolerance )
61  break; // exit loop, found
62  }
63 
64  // check that it was found
65  if(index == (Int_t)fLookupTable.size())
66  index = -1;
67 
68  return index;
69  }
70 
72  if(index<0 || index>(Int_t)fLookupTable.size()) {
73  std::cout << "SamplingSummaryLookup::GetConfidenceLevel, index not in lookup table" << std::endl;
74  return -1;
75  }
76  return fLookupTable[index].first;
77  }
78 
80  if(index<0 || index>(Int_t)fLookupTable.size()) {
81  std::cout << "SamplingSummaryLookup::GetLeftSideTailFraction, index not in lookup table" << std::endl;
82  return -1;
83  }
84  return fLookupTable[index].second;
85  }
86 
87  private:
88  LookupTable fLookupTable; // map ( Index, ( CL, leftside tail prob) )
89 
90  protected:
91  ClassDef(SamplingSummaryLookup,1) // A simple class used by ConfidenceBelt
92  };
93 
94 
95  ///////////////////////////
96  class AcceptanceRegion : public TObject{
97  public:
98  AcceptanceRegion() : fLookupIndex(0), fLowerLimit(0), fUpperLimit(0) {}
99  virtual ~AcceptanceRegion() {}
100 
102  fLookupIndex = lu;
103  fLowerLimit = ll;
104  fUpperLimit = ul;
105  }
106  Int_t GetLookupIndex(){return fLookupIndex;}
107  Double_t GetLowerLimit(){return fLowerLimit;}
108  Double_t GetUpperLimit(){return fUpperLimit;}
109 
110  private:
111  Int_t fLookupIndex; // want a small footprint reference to the RooArgSet for particular parameter point
112  Double_t fLowerLimit; // lower limit on test statistic
113  Double_t fUpperLimit; // upper limit on test statistic
114 
115  protected:
116  ClassDef(AcceptanceRegion,1) // A simple class for acceptance regions used for ConfidenceBelt
117 
118  };
119 
120 ////////////////////////////////////////////////////////////////////////////////
121 
122  class SamplingSummary : public TObject {
123  public:
124  SamplingSummary() : fParameterPointIndex(0) {}
125  virtual ~SamplingSummary() {}
126  SamplingSummary(AcceptanceRegion& ar) : fParameterPointIndex(0) {
127  AddAcceptanceRegion(ar);
128  }
129  Int_t GetParameterPointIndex(){return fParameterPointIndex;}
131  return (SamplingDistribution*) fSamplingDistribution.GetObject(); // dereference TRef
132  }
133  AcceptanceRegion& GetAcceptanceRegion(Int_t index=0){return fAcceptanceRegions[index];}
134 
136  Int_t index = ar.GetLookupIndex();
137  if( fAcceptanceRegions.count(index) !=0) {
138  std::cout << "SamplingSummary::AddAcceptanceRegion, need to implement merging protocol" << std::endl;
139  } else {
140  fAcceptanceRegions[index]=ar;
141  }
142  }
143 
144  private:
145  Int_t fParameterPointIndex; // want a small footprint reference to the RooArgSet for particular parameter point
146  TRef fSamplingDistribution; // persistent pointer to a SamplingDistribution
147  std::map<Int_t, AcceptanceRegion> fAcceptanceRegions;
148 
149  protected:
150  ClassDef(SamplingSummary,1) // A summary of acceptance regions for confidence belt
151 
152  };
153 
154 ////////////////////////////////////////////////////////////////////////////////
155 
156  class ConfidenceBelt : public TNamed {
157 
158  private:
160  std::vector<SamplingSummary> fSamplingSummaries; // composite of several AcceptanceRegions
161  RooAbsData* fParameterPoints; // either a histogram (RooDataHist) or a tree (RooDataSet)
162 
163 
164  public:
165  // constructors,destructors
166  ConfidenceBelt();
167  ConfidenceBelt(const char* name);
168  ConfidenceBelt(const char* name, const char* title);
169  ConfidenceBelt(const char* name, RooAbsData&);
170  ConfidenceBelt(const char* name, const char* title, RooAbsData&);
171  virtual ~ConfidenceBelt();
172 
173  // add after creating a region
174  void AddAcceptanceRegion(RooArgSet&, AcceptanceRegion region, Double_t cl=-1., Double_t leftside=-1.);
175 
176  // add without creating a region, more useful for clients
177  void AddAcceptanceRegion(RooArgSet& point, Int_t dataSetIndex, Double_t lower, Double_t upper, Double_t cl=-1., Double_t leftside=-1.);
178 
179  AcceptanceRegion* GetAcceptanceRegion(RooArgSet&, Double_t cl=-1., Double_t leftside=-1.);
180  Double_t GetAcceptanceRegionMin(RooArgSet&, Double_t cl=-1., Double_t leftside=-1.);
181  Double_t GetAcceptanceRegionMax(RooArgSet&, Double_t cl=-1., Double_t leftside=-1.);
182  std::vector<Double_t> ConfidenceLevels() const ;
183 
184  // Method to return lower limit on a given parameter
185  // Double_t LowerLimit(RooRealVar& param) ; // could provide, but misleading?
186  // Double_t UpperLimit(RooRealVar& param) ; // could provide, but misleading?
187 
188  // do we want it to return list of parameters
189  virtual RooArgSet* GetParameters() const;
190 
191  // check if parameters are correct. (dummy implementation to start)
192  Bool_t CheckParameters(RooArgSet&) const ;
193 
194  protected:
195  ClassDef(ConfidenceBelt,1) // A confidence belt for the Neyman Construction
196 
197  };
198 }
199 
200 #endif
SamplingSummary(AcceptanceRegion &ar)
void Add(Double_t cl, Double_t leftside)
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Short_t Abs(Short_t d)
Definition: TMathBase.h:108
std::pair< Double_t, Double_t > AcceptanceCriteria
Persistent Reference link to a TObject A TRef is a lightweight object pointing to any TObject...
Definition: TRef.h:32
Double_t GetConfidenceLevel(Int_t index)
#define ClassDef(name, id)
Definition: Rtypes.h:297
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
ConfidenceBelt is a concrete implementation of the ConfInterval interface.
std::map< Int_t, AcceptanceCriteria > LookupTable
SamplingSummaryLookup fSamplingSummaryLookup
AcceptanceRegion(Int_t lu, Double_t ll, Double_t ul)
Int_t GetLookupIndex(Double_t cl, Double_t leftside)
std::map< Int_t, AcceptanceRegion > fAcceptanceRegions
AcceptanceRegion & GetAcceptanceRegion(Int_t index=0)
void GetParameters(TFitEditor::FuncParams_t &pars, TF1 *func)
Stores the parameters of the given function into pars.
Definition: TFitEditor.cxx:270
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
This class simply holds a sampling distribution of some test statistic.
constexpr Double_t E()
Definition: TMath.h:74
Namespace for the RooStats classes.
Definition: Asimov.h:20
double Double_t
Definition: RtypesCore.h:55
Double_t GetLeftSideTailFraction(Int_t index)
Mother of all ROOT objects.
Definition: TObject.h:37
std::vector< SamplingSummary > fSamplingSummaries
void AddAcceptanceRegion(AcceptanceRegion &ar)
SamplingDistribution * GetSamplingDistribution()