Logo ROOT  
Reference Guide
rs701_BayesianCalculator.C File Reference

Detailed Description

Bayesian calculator: basic exmple.

View in nbviewer Open in SWAN

RooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby
Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
All rights reserved, please read http://roofit.sourceforge.net/license.txt
[#1] INFO:ObjectHandling -- RooWorkspace::exportToCint(w) INFO: references to all objects in this workspace will be created in CINT in 'namespace w'
Bayesian Result using a Flat prior
0.9% CL central interval: [ 0.50606 - 6.89326 ] or 0.95% CL limits
Bayesian Result using a 1/sqrt(s) prior
0.9% CL central interval: [ 0.0746294 - 5.85425 ] or 0.95% CL limits
#include "RooRealVar.h"
#include "RooWorkspace.h"
#include "RooDataSet.h"
#include "RooPlot.h"
#include "RooMsgService.h"
#include "TCanvas.h"
using namespace RooFit;
using namespace RooStats;
void rs701_BayesianCalculator(bool useBkg = true, double confLevel = 0.90)
{
RooWorkspace *w = new RooWorkspace("w", true);
w->factory("SUM::pdf(s[0.001,15]*Uniform(x[0,1]),b[1,0,2]*Uniform(x))");
w->factory("Gaussian::prior_b(b,1,1)");
w->factory("PROD::model(pdf,prior_b)");
RooAbsPdf *model = w->pdf("model"); // pdf*priorNuisance
RooArgSet nuisanceParameters(*(w->var("b")));
RooAbsRealLValue *POI = w->var("s");
RooAbsPdf *priorPOI = (RooAbsPdf *)w->factory("Uniform::priorPOI(s)");
RooAbsPdf *priorPOI2 = (RooAbsPdf *)w->factory("GenericPdf::priorPOI2('1/sqrt(@0)',s)");
w->factory("n[3]"); // observed number of events
// create a data set with n observed events
RooDataSet data("data", "", RooArgSet(*(w->var("x")), *(w->var("n"))), "n");
data.add(RooArgSet(*(w->var("x"))), w->var("n")->getVal());
// to suppress messages when pdf goes to zero
RooArgSet *nuisPar = 0;
if (useBkg)
nuisPar = &nuisanceParameters;
// if (!useBkg) ((RooRealVar *)w->var("b"))->setVal(0);
double size = 1. - confLevel;
std::cout << "\nBayesian Result using a Flat prior " << std::endl;
BayesianCalculator bcalc(data, *model, RooArgSet(*POI), *priorPOI, nuisPar);
bcalc.SetTestSize(size);
SimpleInterval *interval = bcalc.GetInterval();
double cl = bcalc.ConfidenceLevel();
std::cout << cl << "% CL central interval: [ " << interval->LowerLimit() << " - " << interval->UpperLimit()
<< " ] or " << cl + (1. - cl) / 2 << "% CL limits\n";
RooPlot *plot = bcalc.GetPosteriorPlot();
TCanvas *c1 = new TCanvas("c1", "Bayesian Calculator Result");
c1->Divide(1, 2);
c1->cd(1);
plot->Draw();
c1->Update();
std::cout << "\nBayesian Result using a 1/sqrt(s) prior " << std::endl;
BayesianCalculator bcalc2(data, *model, RooArgSet(*POI), *priorPOI2, nuisPar);
bcalc2.SetTestSize(size);
SimpleInterval *interval2 = bcalc2.GetInterval();
cl = bcalc2.ConfidenceLevel();
std::cout << cl << "% CL central interval: [ " << interval2->LowerLimit() << " - " << interval2->UpperLimit()
<< " ] or " << cl + (1. - cl) / 2 << "% CL limits\n";
RooPlot *plot2 = bcalc2.GetPosteriorPlot();
c1->cd(2);
plot2->Draw();
gPad->SetLogy();
c1->Update();
// observe one event while expecting one background event -> the 95% CL upper limit on s is 4.10
// observe one event while expecting zero background event -> the 95% CL upper limit on s is 4.74
}
Author
Gregory Schott

Definition in file rs701_BayesianCalculator.C.

RooWorkspace.h
RooMsgService.h
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooStats::SimpleInterval::LowerLimit
virtual Double_t LowerLimit()
Definition: SimpleInterval.h:47
TCanvas.h
RooDataSet.h
RooStats::SimpleInterval::UpperLimit
virtual Double_t UpperLimit()
Definition: SimpleInterval.h:49
RooWorkspace::factory
RooFactoryWSTool & factory()
Return instance to factory tool.
Definition: RooWorkspace.cxx:2166
RooStats::SimpleInterval
Definition: SimpleInterval.h:20
RooDataSet::add
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0) override
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
Definition: RooDataSet.cxx:1164
RooFit
Definition: RooCFunction1Binding.h:29
RooFit::FATAL
@ FATAL
Definition: RooGlobalFunc.h:65
BayesianCalculator.h
RooPlot.h
RooMsgService::setGlobalKillBelow
void setGlobalKillBelow(RooFit::MsgLevel level)
Definition: RooMsgService.h:160
SimpleInterval.h
RooPlot
Definition: RooPlot.h:44
RooWorkspace::pdf
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
Definition: RooWorkspace.cxx:1277
RooRealVar.h
RooStats::SimpleInterval::ConfidenceLevel
virtual Double_t ConfidenceLevel() const
return confidence level
Definition: SimpleInterval.h:44
RooWorkspace
Definition: RooWorkspace.h:43
TCanvas
Definition: TCanvas.h:23
RooStats
Definition: Asimov.h:19
RooStats::BayesianCalculator
Definition: BayesianCalculator.h:37
gPad
#define gPad
Definition: TVirtualPad.h:287
RooWorkspace::var
RooRealVar * var(const char *name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found.
Definition: RooWorkspace.cxx:1295
RooDataSet
Definition: RooDataSet.h:33
make_cnn_model.model
model
Definition: make_cnn_model.py:6
RooAbsPdf
Definition: RooAbsPdf.h:40
RooAbsRealLValue
Definition: RooAbsRealLValue.h:31
RooMsgService::instance
static RooMsgService & instance()
Return reference to singleton instance.
Definition: RooMsgService.cxx:363
RooArgSet
Definition: RooArgSet.h:28
c1
return c1
Definition: legend1.C:41