Logo ROOT  
Reference Guide
PdfProposal.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Authors: Kevin Belasco 17/06/2009
3// Authors: Kyle Cranmer 17/06/2009
4/*************************************************************************
5 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12////////////////////////////////////////////////////////////////////////////////
13
14/** \class RooStats::PdfProposal
15 \ingroup Roostats
16
17PdfProposal is a concrete implementation of the ProposalFunction interface.
18It proposes points across the parameter space in the distribution of the
19given PDF.
20
21To make Propose(xPrime, x) dependent on x, configure with
22PdfProposal::AddMapping(varToUpdate, valueToUse). For example, suppose we have:
23
24~~~{.cpp}
25// our parameter
26RooRealVar p("p", "p", 5, 0, 10);
27
28// create mean and sigma for gaussian proposal function
29RooRealVar meanP("meanP", "meanP", 0, 10);
30RooRealVar sigma("sigma", "sigma", 1, 0, 5);
31RooGaussian pGaussian("pGaussian", "pGaussian", p, meanP, sigma);
32
33// configure proposal function
34PdfProposal pdfProposal(pGaussian);
35pdfProposal.AddMapping(meanP, p); // each call of Propose(xPrime, x), meanP in
36 // the proposal function will be updated to
37 // the value of p in x. this will center the
38 // proposal function about x's p when
39 // proposing for xPrime
40
41// To improve performance, PdfProposal has the ability to cache a specified
42// number of proposals. If you don't call this function, the default cache size
43// is 1, which can be slow.
44pdfProposal.SetCacheSize(desiredCacheSize);
45~~~
46
47PdfProposal currently uses a fixed cache size. Adaptive caching methods are in the works
48for future versions.
49*/
50
51
52#include "Rtypes.h"
53
56#include "RooArgSet.h"
57#include "RooDataSet.h"
58#include "RooAbsPdf.h"
59#include "RooMsgService.h"
60#include "RooRealVar.h"
61#include "TIterator.h"
62
63#include <map>
64
66
67using namespace RooFit;
68using namespace RooStats;
69using namespace std;
70
71////////////////////////////////////////////////////////////////////////////////
72/// By default, PdfProposal does NOT own the PDF that serves as the
73/// proposal density function
74
75PdfProposal::PdfProposal() : ProposalFunction()
76{
77 fPdf = NULL;
78 fOwnsPdf = false;
79 fCacheSize = 1;
81 fCache = NULL;
82}
83
84////////////////////////////////////////////////////////////////////////////////
85/// By default, PdfProposal does NOT own the PDF that serves as the
86/// proposal density function
87
89{
90 fPdf = &pdf;
91 fOwnsPdf = false;
92 fCacheSize = 1;
94 fCache = NULL;
95}
96
97////////////////////////////////////////////////////////////////////////////////
98/// determine whether these two RooArgSets represent the same point
99
101{
102 if (x1.equals(x2)) {
103 for (auto const *r : static_range_cast<RooRealVar*>(x1))
104 if (r->getVal() != x2.getRealValue(r->GetName())) {
105 return false;
106 }
107 return true;
108 }
109 return false;
110}
111
112////////////////////////////////////////////////////////////////////////////////
113/// Populate xPrime with a new proposed point
114
116{
117 if (fLastX.getSize() == 0) {
118 // fLastX not yet initialized
120 // generate initial cache
122 if (fMap.size() > 0) {
123 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
124 fIt->first->setVal(fIt->second->getVal(&x));
125 }
126 fCache = fPdf->generate(xPrime, fCacheSize);
127 }
128
129 bool moved = false;
130 if (fMap.size() > 0) {
131 moved = !Equals(fLastX, x);
132
133 // if we've moved, set the values of the variables in the PDF to the
134 // corresponding values of the variables in x, according to the
135 // mappings (i.e. let the variables in x set the given values for the
136 // PDF that will generate xPrime)
137 if (moved) {
138 // update the pdf parameters
140
141 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
142 fIt->first->setVal(fIt->second->getVal(&x));
143
144 // save the new x in fLastX
146 }
147 }
148
149 // generate new cache if necessary
150 if (moved || fCachePosition >= fCacheSize) {
151 delete fCache;
152 fCache = fPdf->generate(xPrime, fCacheSize);
153 fCachePosition = 0;
154 }
155
156 const RooArgSet* proposal = fCache->get(fCachePosition);
158 RooStats::SetParameters(proposal, &xPrime);
159}
160
161////////////////////////////////////////////////////////////////////////////////
162/// Determine whether or not the proposal density is symmetric for
163/// points x1 and x2 - that is, whether the probabilty of reaching x2
164/// from x1 is equal to the probability of reaching x1 from x2
165
167{
168 // kbelasco: is there a better way to do this?
169 return false;
170}
171
172////////////////////////////////////////////////////////////////////////////////
173/// Return the probability of proposing the point x1 given the starting
174/// point x2
175
177{
179 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
180 fIt->first->setVal(fIt->second->getVal(&x2));
183 delete temp;
184 return fPdf->getVal(&x1); // could just as well use x2
185}
186
187////////////////////////////////////////////////////////////////////////////////
188/// specify a mapping between a parameter of the proposal function and
189/// a parameter of interest. this mapping is used to set the value of
190/// proposalParam equal to the value of update to determine the
191/// proposal function.
192/// proposalParam is a parameter of the proposal function that must
193/// be set to the value of update (from the current point) in order to
194/// propose a new point.
195
197{
198 fMaster.add(*update.getParameters((RooAbsData*)NULL));
199 if (update.getParameters((RooAbsData*)NULL)->getSize() == 0)
201 fMap.insert(pair<RooRealVar*, RooAbsReal*>(&proposalParam, &update));
202}
static void update(gsl_integration_workspace *workspace, double a1, double b1, double area1, double error1, double a2, double b2, double area2, double error2)
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:317
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:61
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,...
Definition: RooAbsPdf.h:60
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
PdfProposal is a concrete implementation of the ProposalFunction interface.
Definition: PdfProposal.h:30
std::map< RooRealVar *, RooAbsReal * >::iterator fIt
map of values in pdf to update
Definition: PdfProposal.h:106
RooArgSet fMaster
the cached proposal data set
Definition: PdfProposal.h:111
bool IsSymmetric(RooArgSet &x1, RooArgSet &x2) override
Determine whether or not the proposal density is symmetric for points x1 and x2 - that is,...
virtual void AddMapping(RooRealVar &proposalParam, RooAbsReal &update)
specify a mapping between a parameter of the proposal function and a parameter of interest.
bool fOwnsPdf
pointers to master variables needed for updates
Definition: PdfProposal.h:112
RooDataSet * fCache
our position in the cached proposal data set
Definition: PdfProposal.h:110
virtual bool Equals(RooArgSet &x1, RooArgSet &x2)
whether we own the proposal density function
RooArgSet fLastX
pdf iterator
Definition: PdfProposal.h:107
void Propose(RooArgSet &xPrime, RooArgSet &x) override
Populate xPrime with a new proposed point.
Int_t fCacheSize
the last point we were at
Definition: PdfProposal.h:108
Int_t fCachePosition
how many points to generate each time
Definition: PdfProposal.h:109
std::map< RooRealVar *, RooAbsReal * > fMap
the proposal density function
Definition: PdfProposal.h:105
PdfProposal()
By default, PdfProposal does NOT own the PDF that serves as the proposal density function.
Definition: PdfProposal.cxx:75
double GetProposalDensity(RooArgSet &x1, RooArgSet &x2) override
Return the probability of proposing the point x1 given the starting point x2.
ProposalFunction is an interface for all proposal functions that would be used with a Markov Chain Mo...
Double_t x[n]
Definition: legend1.C:17
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition: Common.h:18
Namespace for the RooStats classes.
Definition: Asimov.h:19
void SetParameters(const RooArgSet *desiredVals, RooArgSet *paramsToChange)
Definition: RooStatsUtils.h:65