Logo ROOT   6.21/01
Reference Guide
HybridOriginalDemo.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roostats
3 /// \notebook -js
4 /// Example on how to use the HybridCalculatorOriginal class
5 ///
6 /// With this example, you should get: CL_sb = 0.130 and CL_b = 0.946
7 /// (if data had -2lnQ = -3.0742).
8 ///
9 /// \macro_image
10 /// \macro_output
11 /// \macro_code
12 ///
13 /// \author Gregory Schott
14 
15 #include "RooRandom.h"
16 #include "RooRealVar.h"
17 #include "RooGaussian.h"
18 #include "RooPolynomial.h"
19 #include "RooArgSet.h"
20 #include "RooAddPdf.h"
21 #include "RooDataSet.h"
22 #include "RooExtendPdf.h"
23 #include "RooConstVar.h"
25 #include "RooStats/HybridResult.h"
26 #include "RooStats/HybridPlot.h"
27 
28 void HybridOriginalDemo(int ntoys = 1000)
29 {
30  using namespace RooFit;
31  using namespace RooStats;
32 
33  // set RooFit random seed
35 
36  /// build the models for background and signal+background
37  RooRealVar x("x", "", -3, 3);
38  RooArgList observables(x); // variables to be generated
39 
40  // gaussian signal
41  RooGaussian sig_pdf("sig_pdf", "", x, RooConst(0.0), RooConst(0.8));
42  RooRealVar sig_yield("sig_yield", "", 20, 0, 300);
43 
44  // flat background (extended PDF)
45  RooPolynomial bkg_pdf("bkg_pdf", "", x, RooConst(0));
46  RooRealVar bkg_yield("bkg_yield", "", 40, 0, 300);
47  RooExtendPdf bkg_ext_pdf("bkg_ext_pdf", "", bkg_pdf, bkg_yield);
48 
49  // bkg_yield.setConstant(kTRUE);
50  sig_yield.setConstant(kTRUE);
51 
52  // total sig+bkg (extended PDF)
53  RooAddPdf tot_pdf("tot_pdf", "", RooArgList(sig_pdf, bkg_pdf), RooArgList(sig_yield, bkg_yield));
54 
55  // build the prior PDF on the parameters to be integrated
56  // gaussian contraint on the background yield ( N_B = 40 +/- 10 ie. 25% )
57  RooGaussian bkg_yield_prior("bkg_yield_prior", "", bkg_yield, RooConst(bkg_yield.getVal()), RooConst(10.));
58 
59  RooArgSet nuisance_parameters(bkg_yield); // variables to be integrated
60 
61  /// generate a data sample
62  RooDataSet *data = tot_pdf.generate(observables, RooFit::Extended());
63 
64  // run HybridCalculator on those inputs
65  // use interface from HypoTest calculator by default
66 
67  HybridCalculatorOriginal myHybridCalc(*data, tot_pdf, bkg_ext_pdf, &nuisance_parameters, &bkg_yield_prior);
68 
69  // here I use the default test statistics: 2*lnQ (optional)
70  myHybridCalc.SetTestStatistic(1);
71  // myHybridCalc.SetTestStatistic(3); // profile likelihood ratio
72 
73  myHybridCalc.SetNumberOfToys(ntoys);
74  myHybridCalc.UseNuisance(true);
75 
76  // for speed up generation (do binned data)
77  myHybridCalc.SetGenerateBinned(false);
78 
79  // calculate by running ntoys for the S+B and B hypothesis and retrieve the result
80  HybridResult *myHybridResult = myHybridCalc.GetHypoTest();
81 
82  if (!myHybridResult) {
83  std::cerr << "\nError returned from Hypothesis test" << std::endl;
84  return;
85  }
86 
87  /// nice plot of the results
88  HybridPlot *myHybridPlot =
89  myHybridResult->GetPlot("myHybridPlot", "Plot of results with HybridCalculatorOriginal", 100);
90  myHybridPlot->Draw();
91 
92  /// recover and display the results
93  double clsb_data = myHybridResult->CLsplusb();
94  double clb_data = myHybridResult->CLb();
95  double cls_data = myHybridResult->CLs();
96  double data_significance = myHybridResult->Significance();
97  double min2lnQ_data = myHybridResult->GetTestStat_data();
98 
99  /// compute the mean expected significance from toys
100  double mean_sb_toys_test_stat = myHybridPlot->GetSBmean();
101  myHybridResult->SetDataTestStatistics(mean_sb_toys_test_stat);
102  double toys_significance = myHybridResult->Significance();
103 
104  std::cout << "Completed HybridCalculatorOriginal example:\n";
105  std::cout << " - -2lnQ = " << min2lnQ_data << endl;
106  std::cout << " - CL_sb = " << clsb_data << std::endl;
107  std::cout << " - CL_b = " << clb_data << std::endl;
108  std::cout << " - CL_s = " << cls_data << std::endl;
109  std::cout << " - significance of data = " << data_significance << std::endl;
110  std::cout << " - mean significance of toys = " << toys_significance << std::endl;
111 }
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
HybridPlot * GetPlot(const char *name, const char *title, int n_bins)
prepare a plot showing a result and return a pointer to a HybridPlot object the needed arguments are:...
virtual Double_t CLb() const
Convert NullPValue into a "confidence level".
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
HybridCalculatorOriginal class.
void Draw(const char *options="")
Draw on current pad.
Definition: HybridPlot.cxx:138
RooCmdArg Extended(Bool_t flag=kTRUE)
void SetDataTestStatistics(double testStat_data_val)
set the value of the test statistics on data
virtual Double_t Significance() const
familiar name for the Null p-value in terms of 1-sided Gaussian significance
Double_t x[n]
Definition: legend1.C:17
virtual void SetSeed(ULong_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:597
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:54
This class provides the plots for the result of a study performed with the HybridCalculatorOriginal c...
Definition: HybridPlot.h:36
RooRealVar represents a fundamental (non-derived) real-valued object.
Definition: RooRealVar.h:36
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
virtual Double_t CLs() const
is simply (not a method, but a quantity)
virtual Double_t CLsplusb() const
Convert AlternatePValue into a "confidence level".
double GetSBmean()
Get SB histo mean.
Definition: HybridPlot.h:77
Namespace for the RooStats classes.
Definition: Asimov.h:20
double GetTestStat_data()
Get test statistics value for data.
Definition: HybridResult.h:57
RooExtendPdf is a wrapper around an existing PDF that adds a parameteric extended likelihood term to ...
Definition: RooExtendPdf.h:22
RooConstVar & RooConst(Double_t val)
RooPolynomial implements a polynomial p.d.f of the form By default, the coefficient is chosen to be...
Definition: RooPolynomial.h:28
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
const Bool_t kTRUE
Definition: RtypesCore.h:87
Class encapsulating the result of the HybridCalculatorOriginal.
Definition: HybridResult.h:25