Logo ROOT   6.10/09
Reference Guide
BayesianCalculator.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_BayesianCalculator
12 #define ROOSTATS_BayesianCalculator
13 
14 #include "TNamed.h"
15 
16 #include "Math/IFunctionfwd.h"
17 
18 #include "RooArgSet.h"
19 
21 
23 
24 class RooAbsData;
25 class RooAbsPdf;
26 class RooPlot;
27 class RooAbsReal;
28 class TF1;
29 
30 
31 namespace RooStats {
32 
33  class ModelConfig;
34  class SimpleInterval;
35 
36  class BayesianCalculator : public IntervalCalculator, public TNamed {
37 
38  public:
39 
40  // constructor
42 
44  RooAbsPdf& pdf,
45  const RooArgSet& POI,
46  RooAbsPdf& priorPdf,
47  const RooArgSet* nuisanceParameters = 0 );
48 
50  ModelConfig& model );
51 
52  // destructor
53  virtual ~BayesianCalculator();
54 
55  // get the plot with option to get it normalized
56  RooPlot* GetPosteriorPlot(bool norm = false, double precision = 0.01) const;
57 
58  // return posterior pdf (object is managed by the BayesianCalculator class)
59  RooAbsPdf* GetPosteriorPdf() const;
60  // return posterior function (object is managed by the BayesianCalculator class)
62 
63  // compute the interval. By Default a central interval is computed
64  // By using SetLeftTileFraction can control if central/ upper/lower interval
65  // For shortest interval use SetShortestInterval(true)
66  virtual SimpleInterval* GetInterval() const ;
67 
68  virtual void SetData( RooAbsData & data ) {
69  fData = &data;
70  ClearAll();
71  }
72 
73 
74  // set the model via the ModelConfig
75  virtual void SetModel( const ModelConfig& model );
76 
77  // specify the parameters of interest in the interval
78  virtual void SetParameters(const RooArgSet& set) { fPOI.removeAll(); fPOI.add(set); }
79 
80  // specify the nuisance parameters (eg. the rest of the parameters)
82 
83  // Set only the Prior Pdf
84  virtual void SetPriorPdf(RooAbsPdf& pdf) { fPriorPdf = &pdf; }
85 
86  // set the conditional observables which will be used when creating the NLL
87  // so the pdf's will not be normalized on the conditional observables when computing the NLL
89 
90  // set the size of the test (rate of Type I error) ( Eg. 0.05 for a 95% Confidence Interval)
91  virtual void SetTestSize( Double_t size ) {
92  fSize = size;
93  fValidInterval = false;
94  }
95  // set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval)
96  virtual void SetConfidenceLevel( Double_t cl ) { SetTestSize(1.-cl); }
97  // Get the size of the test (eg. rate of Type I error)
98  virtual Double_t Size() const { return fSize; }
99  // Get the Confidence level for the test
100  virtual Double_t ConfidenceLevel() const { return 1.-fSize; }
101 
102  // set the fraction of probability content on the left tail
103  // Central limits use 0.5 (default case)
104  // for upper limits it is 0 and 1 for lower limit
105  // For shortest intervals a negative value (i.e. -1) must be given
106  void SetLeftSideTailFraction(Double_t leftSideFraction ) {fLeftSideFraction = leftSideFraction;}
107 
108  // set the Bayesian calculator to compute the shortest interval (default is central interval)
109  // to switch off SetLeftSideTailFraction to the right value
111 
112  // set the precision of the Root Finder
113  void SetBrfPrecision( double precision ) { fBrfPrecision = precision; }
114 
115  // use directly the approximate posterior function obtained by binning it in nbins
116  // by default the cdf is used by integrating the posterior
117  // if a value of nbin <= 0 the cdf function will be used
118  void SetScanOfPosterior(int nbin = 100) { fNScanBins = nbin; }
119 
120  // set the number of iterations when running a MC integration algorithm
121  // If not set use default algorithmic values
122  // In case of ToyMC sampling of the nuisance the value is 100
123  // In case of using the GSL MCintegrations types the default value is
124  // defined in ROOT::Math::IntegratorMultiDimOptions::DefaultNCalls()
125  virtual void SetNumIters(Int_t numIters) { fNumIterations = numIters; }
126 
127  // set the integration type (possible type are) :
128  void SetIntegrationType(const char * type);
129 
130  // return the mode (most probable value of the posterior function)
131  double GetMode() const;
132 
133  // force the nuisance pdf when using the toy mc sampling
134  void ForceNuisancePdf(RooAbsPdf & pdf) { fNuisancePdf = &pdf; }
135 
136  protected:
137 
138  void ClearAll() const;
139 
140  void ApproximatePosterior() const;
141 
142  void ComputeIntervalFromApproxPosterior(double c1, double c2) const;
143 
144  void ComputeIntervalFromCdf(double c1, double c2) const;
145 
146  void ComputeIntervalUsingRooFit(double c1, double c2) const;
147 
148  void ComputeShortestInterval() const;
149 
150  private:
151 
152  // plan to replace the above: return a SimpleInterval integrating
153  // over all other parameters except the one specified as argument
154  // virtual SimpleInterval* GetInterval( RooRealVar* parameter ) const { return 0; }
155 
156  RooAbsData* fData; // data set
157  RooAbsPdf* fPdf; // model pdf (could contain the nuisance pdf as constraint term)
158  RooArgSet fPOI; // POI
159  RooAbsPdf* fPriorPdf; // prior pdf (typically for the POI)
160  RooAbsPdf* fNuisancePdf; // nuisance pdf (needed when using nuisance sampling technique)
161  RooArgSet fNuisanceParameters; // nuisance parameters
162  RooArgSet fConditionalObs ; // conditional observables
163 
164  mutable RooAbsPdf* fProductPdf; // internal pointer to model * prior
165  mutable RooAbsReal* fLogLike; // internal pointer to log likelihood function
166  mutable RooAbsReal* fLikelihood; // internal pointer to likelihood function
167  mutable RooAbsReal* fIntegratedLikelihood; // integrated likelihood function, i.e - unnormalized posterior function
168  mutable RooAbsPdf* fPosteriorPdf; // normalized (on the poi) posterior pdf
169  mutable ROOT::Math::IGenFunction * fPosteriorFunction; // function representing the posterior
170  mutable TF1 * fApproxPosterior; // TF1 representing the scanned posterior function
171  mutable Double_t fLower; // computer lower interval bound
172  mutable Double_t fUpper; // upper interval bound
173  mutable Double_t fNLLMin; // minimum value of Nll
174  double fSize; // size used for getting the interval
175  double fLeftSideFraction; // fraction of probability content on left side of interval
176  double fBrfPrecision; // root finder precision
177  mutable int fNScanBins; // number of bins to scan, if = -1 no scan is done (default)
178  int fNumIterations; // number of iterations (when using ToyMC)
180 
182 
183  protected:
184 
185  ClassDef(BayesianCalculator,2) // BayesianCalculator class
186 
187  };
188 }
189 
190 #endif
virtual void SetData(RooAbsData &data)
Set the DataSet ( add to the the workspace if not already there ?)
RooAbsPdf * GetPosteriorPdf() const
Build and return the posterior pdf (i.e posterior function normalized to all range of poi) Note that ...
ROOT::Math::IGenFunction * fPosteriorFunction
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:30
Interface (abstract class) for generic functions objects of one-dimension Provides a method to evalua...
Definition: IFunction.h:134
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:86
IntervalCalculator is an interface class for a tools which produce RooStats ConfIntervals.
void SetScanOfPosterior(int nbin=100)
return c1
Definition: legend1.C:41
void SetIntegrationType(const char *type)
set the integration type (possible type are) :
Basic string class.
Definition: TString.h:129
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
double GetMode() const
Returns the value of the parameter for the point in parameter-space that is the most likely...
void SetBrfPrecision(double precision)
virtual void SetNuisanceParameters(const RooArgSet &set)
virtual void SetParameters(const RooArgSet &set)
#define ClassDef(name, id)
Definition: Rtypes.h:297
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual void SetNumIters(Int_t numIters)
void ApproximatePosterior() const
approximate posterior in nbins using a TF1 scan the poi values and evaluate the posterior at each poi...
RooPlot * GetPosteriorPlot(bool norm=false, double precision=0.01) const
return a RooPlot with the posterior and the credibility region NOTE: User takes ownership of the retu...
void ClearAll() const
clear all cached pdf objects
void ComputeIntervalUsingRooFit(double c1, double c2) const
internal function compute the interval using RooFit
virtual void SetConfidenceLevel(Double_t cl)
set the confidence level for the interval (e.g. 0.95 for a 95% Confidence Interval) ...
void SetLeftSideTailFraction(Double_t leftSideFraction)
virtual Double_t ConfidenceLevel() const
Get the Confidence level for the test.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
virtual void SetModel(const ModelConfig &model)
set the model to use The model pdf, prior pdf, parameter of interest and nuisances will be taken acco...
void ComputeShortestInterval() const
compute the shortest interval
Namespace for the RooStats classes.
Definition: Asimov.h:20
RooAbsReal * GetPosteriorFunction() const
Build and return the posterior function (not normalized) as a RooAbsReal the posterior is obtained fr...
return c2
Definition: legend2.C:14
void ComputeIntervalFromCdf(double c1, double c2) const
internal function compute the interval using Cdf integration
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
virtual void SetPriorPdf(RooAbsPdf &pdf)
int type
Definition: TGX11.cxx:120
virtual void SetConditionalObservables(const RooArgSet &set)
virtual void SetTestSize(Double_t size)
set the size of the test (rate of Type I error) ( e.g. 0.05 for a 95% Confidence Interval) ...
BayesianCalculator()
default constructor
SimpleInterval is a concrete implementation of the ConfInterval interface.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
1-Dim function class
Definition: TF1.h:150
virtual Double_t Size() const
Get the size of the test (eg. rate of Type I error)
virtual SimpleInterval * GetInterval() const
Compute the interval.
void ComputeIntervalFromApproxPosterior(double c1, double c2) const
compute the interval using the approximate posterior function
void ForceNuisancePdf(RooAbsPdf &pdf)
double norm(double *x, double *p)
Definition: unuranDistr.cxx:40
BayesianCalculator is a concrete implementation of IntervalCalculator, providing the computation of a...