ROOT   Reference Guide
Zbi_Zgamma.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roostats
3/// \notebook -js
4/// Demonstrate Z_Bi = Z_Gamma
5///
6/// \macro_image
7/// \macro_output
8/// \macro_code
9///
10/// \authors Kyle Cranmer & Wouter Verkerke
11
12#include "RooRealVar.h"
13#include "RooProdPdf.h"
14#include "RooWorkspace.h"
15#include "RooDataSet.h"
16#include "TCanvas.h"
17#include "TH1.h"
18
19using namespace RooFit;
20using namespace RooStats;
21
22void Zbi_Zgamma()
23{
24
25 // Make model for prototype on/off problem
26 // Pois(x | s+b) * Pois(y | tau b )
27 // for Z_Gamma, use uniform prior on b.
28 RooWorkspace *w1 = new RooWorkspace("w", true);
29 w1->factory("Poisson::px(x[150,0,500],sum::splusb(s[0,0,100],b[100,0,300]))");
30 w1->factory("Poisson::py(y[100,0,500],prod::taub(tau[1.],b))");
31 w1->factory("Uniform::prior_b(b)");
32
33 // construct the Bayesian-averaged model (eg. a projection pdf)
34 // p'(x|s) = \int db p(x|s+b) * [ p(y|b) * prior(b) ]
35 w1->factory("PROJ::averagedModel(PROD::foo(px|b,py,prior_b),b)");
36
37 // plot it, blue is averaged model, red is b known exactly
38 RooPlot *frame = w1->var("x")->frame();
39 w1->pdf("averagedModel")->plotOn(frame);
40 w1->pdf("px")->plotOn(frame, LineColor(kRed));
41 frame->Draw();
42
43 // compare analytic calculation of Z_Bi
44 // with the numerical RooFit implementation of Z_Gamma
45 // for an example with x = 150, y = 100
46
47 // numeric RooFit Z_Gamma
48 w1->var("y")->setVal(100);
49 w1->var("x")->setVal(150);
50 RooAbsReal *cdf = w1->pdf("averagedModel")->createCdf(*w1->var("x"));
51 cdf->getVal(); // get ugly print messages out of the way
52
53 cout << "Hybrid p-value = " << cdf->getVal() << endl;
54 cout << "Z_Gamma Significance = " << PValueToSignificance(1 - cdf->getVal()) << endl;
55
56 // analytic Z_Bi
57 double Z_Bi = NumberCountingUtils::BinomialWithTauObsZ(150, 100, 1);
58 std::cout << "Z_Bi significance estimation: " << Z_Bi << std::endl;
59
60 // OUTPUT
61 // Hybrid p-value = 0.999058
62 // Z_Gamma Significance = 3.10804
63 // Z_Bi significance estimation: 3.10804
64}
@ kRed
Definition: Rtypes.h:66
RooAbsReal * createCdf(const RooArgSet &iset, const RooArgSet &nset=RooArgSet())
Create a cumulative distribution function of this p.d.f in terms of the observables listed in iset.
Definition: RooAbsPdf.cxx:3413
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none(), const RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Definition: RooAbsPdf.h:121
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create a new RooPlot on the heap with a drawing frame initialized for this object,...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:282
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43
RooRealVar * var(const char *name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found.
RooFactoryWSTool & factory()
Return instance to factory tool.
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
RooCmdArg LineColor(Color_t color)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Double_t BinomialWithTauObsZ(Double_t nObs, Double_t bExp, Double_t tau)
Namespace for the RooStats classes.
Definition: Asimov.h:19
Double_t PValueToSignificance(Double_t pvalue)
returns one-sided significance corresponding to a p-value
Definition: RooStatsUtils.h:51