Logo ROOT   6.08/07
Reference Guide
HybridCalculatorOriginal.h
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 
3 /*************************************************************************
4  * Project: RooStats *
5  * Package: RooFit/RooStats *
6  * Authors: *
7  * Danilo Piparo, Gregory Schott *
8  *************************************************************************
9  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
10  * All rights reserved. *
11  * *
12  * For the licensing terms see $ROOTSYS/LICENSE. *
13  * For the list of contributors see $ROOTSYS/README/CREDITS. *
14  *************************************************************************/
15 
16 #ifndef ROOSTATS_HybridCalculatorOriginal
17 #define ROOSTATS_HybridCalculatorOriginal
18 
19 #ifndef ROOSTATS_HypoTestCalculator
21 #endif
22 
23 #include <vector>
24 
25 
26 #ifndef ROOSTATS_HybridResult
27 #include "RooStats/HybridResult.h"
28 #endif
29 
30 #ifndef ROOSTATS_ModelConfig
31 #include "RooStats/ModelConfig.h"
32 #endif
33 
34 class TH1;
35 
36 namespace RooStats {
37 
38  class HybridResult;
39 
40  /**
41 
42 
43 HybridCalculatorOriginal class. This class is depracated and it is replaced by the HybridCalculator.
44 This is a fresh rewrite in RooStats of
45  RooStatsCms/LimitCalculator developped by D. Piparo and G. Schott
46 Authors: D. Piparo, G. Schott - Universitaet Karlsruhe
47 
48 The class is born from the need to have an implementation of the CLs
49 method that could take advantage from the RooFit Package.
50 The basic idea is the following:
51 - Instantiate an object specifying a signal+background model, a background model and a dataset.
52 - Perform toy MC experiments to know the distributions of -2lnQ
53 - Calculate the CLsb and CLs values as "integrals" of these distributions.
54 
55 The class allows the user to input models as RooAbsPdf ( TH1 object could be used
56 by using the RooHistPdf class)
57 The pdfs must be "extended": for more information please refer to
58 http://roofit.sourceforge.net). The dataset can be entered as a
59 RooAbsData objects.
60 
61 Unlike the TLimit Class a complete MC generation is performed at each step
62 and not a simple Poisson fluctuation of the contents of the bins.
63 Another innovation is the treatment of the nuisance parameters. The user
64 can input in the constructor nuisance parameters.
65 To include the information that we have about the nuisance parameters a prior
66 PDF (RooAbsPdf) should be specified
67 
68 Different test statistic can be used (likelihood ratio, number of events or
69 profile likelihood ratio. The default is the likelihood ratio.
70 See the method SetTestStatistic.
71 
72 The number of toys to be generated is controlled by SetNumberOfToys(n).
73 
74 The result of the calculations is returned as a HybridResult object pointer.
75 
76 see also the following interesting references:
77 - Alex Read, "Presentation of search results: the CLs technique",
78  Journal of Physics G: Nucl. Part. Phys. 28 2693-2704 (2002).
79  see http://www.iop.org/EJ/abstract/0954-3899/28/10/313/
80 
81 - Alex Read, "Modified Frequentist Analysis of Search Results (The CLs Method)" CERN 2000-005 (30 May 2000)
82 
83 - V. Bartsch, G.Quast, "Expected signal observability at future experiments" CMS NOTE 2005/004
84 
85 - http://root.cern.ch/root/html/src/TLimit.html
86 */
87 
88 
90 
91  public:
92 
93 
94  /// Dummy Constructor with only name
95  explicit HybridCalculatorOriginal(const char *name = 0);
96 
97  /// Constructor for HybridCalculator from pdf instances but without a data-set
99  RooAbsPdf& b_model,
100  RooArgList& observables,
101  const RooArgSet* nuisance_parameters = 0,
102  RooAbsPdf* prior_pdf = 0,
103  bool GenerateBinned = false, int testStatistics = 1, int ntoys = 1000 );
104 
105  /// Constructor for HybridCalculator using a data set and pdf instances
107  RooAbsPdf& sb_model,
108  RooAbsPdf& b_model,
109  const RooArgSet* nuisance_parameters = 0,
110  RooAbsPdf* prior_pdf = 0,
111  bool GenerateBinned = false, int testStatistics = 1, int ntoys = 1000 );
112 
113 
114  /// Constructor passing a ModelConfig for the SBmodel and a ModelConfig for the B Model
116  const ModelConfig& sb_model,
117  const ModelConfig& b_model,
118  bool GenerateBinned = false, int testStatistics = 1, int ntoys = 1000 );
119 
120 
121  public:
122 
123  /// Destructor of HybridCalculator
124  virtual ~HybridCalculatorOriginal();
125 
126  /// inherited methods from HypoTestCalculator interface
127  virtual HybridResult* GetHypoTest() const;
128 
129  // inherited setter methods from HypoTestCalculator
130 
131 
132  // set the model for the null hypothesis (only B)
133  virtual void SetNullModel(const ModelConfig & );
134  // set the model for the alternate hypothesis (S+B)
135  virtual void SetAlternateModel(const ModelConfig & );
136 
137 
138  // Set a common PDF for both the null and alternate
139  virtual void SetCommonPdf(RooAbsPdf & pdf) { fSbModel = &pdf; }
140  // Set the PDF for the null (only B)
141  virtual void SetNullPdf(RooAbsPdf& pdf) { fBModel = &pdf; }
142  // Set the PDF for the alternate hypothesis ( i.e. S+B)
143  virtual void SetAlternatePdf(RooAbsPdf& pdf) { fSbModel = &pdf; }
144 
145  // Set the DataSet
146  virtual void SetData(RooAbsData& data) { fData = &data; }
147 
148  // set parameter values for the null if using a common PDF
149  virtual void SetNullParameters(const RooArgSet& ) { } // not needed
150  // set parameter values for the alternate if using a common PDF
151  virtual void SetAlternateParameters(const RooArgSet&) {} // not needed
152 
153  // additional methods specific for HybridCalculator
154  // set a prior pdf for the nuisance parameters
155  void SetNuisancePdf(RooAbsPdf & prior_pdf) {
156  fPriorPdf = &prior_pdf;
157  fUsePriorPdf = true; // if set by default turn it on
158  }
159 
160  // set the nuisance parameters to be marginalized
161  void SetNuisanceParameters(const RooArgSet & params) { fNuisanceParameters = &params; }
162 
163  // set number of toy MC (Default is 1000)
164  void SetNumberOfToys(unsigned int ntoys) { fNToys = ntoys; }
165 
166  // return number of toys used
167  unsigned int GetNumberOfToys() const { return fNToys; }
168 
169  // control use of the pdf for the nuisance parameter and marginalize them
170  void UseNuisance(bool on = true) { fUsePriorPdf = on; }
171 
172  // control to use bin data generation
173  void SetGenerateBinned(bool on = true) { fGenerateBinned = on; }
174 
175  /// set the desired test statistics:
176  /// index=1 : 2 * log( L_sb / L_b ) (DEFAULT)
177  /// index=2 : number of generated events
178  /// index=3 : profiled likelihood ratio
179  /// if the index is different to any of those values, the default is used
180  void SetTestStatistic(int index);
181 
182  HybridResult* Calculate(TH1& data, unsigned int nToys, bool usePriors) const;
183  HybridResult* Calculate(RooAbsData& data, unsigned int nToys, bool usePriors) const;
184  HybridResult* Calculate(unsigned int nToys, bool usePriors) const;
185  void PrintMore(const char* options) const;
186 
187  void PatchSetExtended(bool on = true) { fTmpDoExtended = on; std::cout << "extended patch set to " << on << std::endl; } // patch to test with RooPoisson (or other non-extended models)
188 
189  private:
190 
191  void RunToys(std::vector<double>& bVals, std::vector<double>& sbVals, unsigned int nToys, bool usePriors) const;
192 
193  // check input parameters before performing the calculation
194  bool DoCheckInputs() const;
195 
196  unsigned int fTestStatisticsIdx; // Index of the test statistics to use
197  unsigned int fNToys; // number of Toys MC
198  RooAbsPdf* fSbModel; // The pdf of the signal+background model
199  RooAbsPdf* fBModel; // The pdf of the background model
200  mutable RooArgList* fObservables; // Collection of the observables of the model
201  const RooArgSet* fNuisanceParameters; // Collection of the nuisance parameters in the model
202  RooAbsPdf* fPriorPdf; // Prior PDF of the nuisance parameters
203  RooAbsData * fData; // pointer to the data sets
204  bool fGenerateBinned; //Flag to control binned generation
205  bool fUsePriorPdf; // use a prior for nuisance parameters
207 
208 // TString fSbModelName; // name of pdf of the signal+background model
209 // TString fBModelName; // name of pdf of the background model
210 // TString fPriorPdfName; // name of pdf of the background model
211 // TString fDataName; // name of the dataset in the workspace
212 
213  protected:
214  ClassDef(HybridCalculatorOriginal,1) // Hypothesis test calculator using a Bayesian-frequentist hybrid method
215  };
216 
217 }
218 
219 #endif
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:52
virtual void SetCommonPdf(RooAbsPdf &pdf)
HybridCalculatorOriginal class.
virtual void SetNullPdf(RooAbsPdf &pdf)
HybridResult * Calculate(TH1 &data, unsigned int nToys, bool usePriors) const
#define ClassDef(name, id)
Definition: Rtypes.h:254
virtual void SetNullParameters(const RooArgSet &)
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
virtual HybridResult * GetHypoTest() const
inherited methods from HypoTestCalculator interface
virtual void SetAlternateModel(const ModelConfig &)
virtual void SetNullModel(const ModelConfig &)
void RunToys(std::vector< double > &bVals, std::vector< double > &sbVals, unsigned int nToys, bool usePriors) const
HypoTestCalculator is an interface class for a tools which produce RooStats HypoTestResults.
HybridCalculatorOriginal(const char *name=0)
Dummy Constructor with only name.
virtual void SetData(RooAbsData &data)
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
void SetNuisanceParameters(const RooArgSet &params)
Namespace for the RooStats classes.
Definition: Asimov.h:20
The TH1 histogram class.
Definition: TH1.h:80
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual void SetAlternateParameters(const RooArgSet &)
void PrintMore(const char *options) const
char name[80]
Definition: TGX11.cxx:109
virtual void SetAlternatePdf(RooAbsPdf &pdf)
void SetTestStatistic(int index)
set the desired test statistics: index=1 : 2 * log( L_sb / L_b ) (DEFAULT) index=2 : number of genera...
virtual ~HybridCalculatorOriginal()
Destructor of HybridCalculator.