Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
76{
77 fPdf = NULL;
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;
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 TIterator* it = x1.createIterator();
104 RooRealVar* r;
105 while ((r = (RooRealVar*)it->Next()) != NULL)
106 if (r->getVal() != x2.getRealValue(r->GetName())) {
107 delete it;
108 return kFALSE;
109 }
110 delete it;
111 return kTRUE;
112 }
113 return kFALSE;
114}
115
116////////////////////////////////////////////////////////////////////////////////
117/// Populate xPrime with a new proposed point
118
120{
121 if (fLastX.getSize() == 0) {
122 // fLastX not yet initialized
124 // generate initial cache
126 if (fMap.size() > 0) {
127 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
128 fIt->first->setVal(fIt->second->getVal(&x));
129 }
130 fCache = fPdf->generate(xPrime, fCacheSize);
131 }
132
133 Bool_t moved = false;
134 if (fMap.size() > 0) {
135 moved = !Equals(fLastX, x);
136
137 // if we've moved, set the values of the variables in the PDF to the
138 // corresponding values of the variables in x, according to the
139 // mappings (i.e. let the variables in x set the given values for the
140 // PDF that will generate xPrime)
141 if (moved) {
142 // update the pdf parameters
144
145 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
146 fIt->first->setVal(fIt->second->getVal(&x));
147
148 // save the new x in fLastX
150 }
151 }
152
153 // generate new cache if necessary
154 if (moved || fCachePosition >= fCacheSize) {
155 delete fCache;
156 fCache = fPdf->generate(xPrime, fCacheSize);
157 fCachePosition = 0;
158 }
159
160 const RooArgSet* proposal = fCache->get(fCachePosition);
162 RooStats::SetParameters(proposal, &xPrime);
163}
164
165////////////////////////////////////////////////////////////////////////////////
166/// Determine whether or not the proposal density is symmetric for
167/// points x1 and x2 - that is, whether the probabilty of reaching x2
168/// from x1 is equal to the probability of reaching x1 from x2
169
171{
172 // kbelasco: is there a better way to do this?
173 return false;
174}
175
176////////////////////////////////////////////////////////////////////////////////
177/// Return the probability of proposing the point x1 given the starting
178/// point x2
179
181{
183 for (fIt = fMap.begin(); fIt != fMap.end(); fIt++)
184 fIt->first->setVal(fIt->second->getVal(&x2));
187 delete temp;
188 return fPdf->getVal(&x1); // could just as well use x2
189}
190
191////////////////////////////////////////////////////////////////////////////////
192/// specify a mapping between a parameter of the proposal function and
193/// a parameter of interest. this mapping is used to set the value of
194/// proposalParam equal to the value of update to determine the
195/// proposal function.
196/// proposalParam is a parameter of the proposal function that must
197/// be set to the value of update (from the current point) in order to
198/// propose a new point.
199
201{
202 fMaster.add(*update.getParameters((RooAbsData*)NULL));
203 if (update.getParameters((RooAbsData*)NULL)->getSize() == 0)
205 fMap.insert(pair<RooRealVar*, RooAbsReal*>(&proposalParam, &update));
206}
ROOT::R::TRInterface & r
Definition Object.C:4
static void update(gsl_integration_workspace *workspace, double a1, double b1, double area1, double error1, double a2, double b2, double area2, double error2)
static const double x2[5]
static const double x1[5]
const Bool_t kFALSE
Definition RtypesCore.h:101
const Bool_t kTRUE
Definition RtypesCore.h:100
#define ClassImp(name)
Definition Rtypes.h:364
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition RooAbsArg.h:309
Int_t getSize() const
virtual RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE)
Add a clone of the specified argument to list.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:82
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:58
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:64
Double_t 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:35
virtual 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:39
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
virtual Bool_t IsSymmetric(RooArgSet &x1, RooArgSet &x2)
Determine whether or not the proposal density is symmetric for points x1 and x2 - that is,...
RooArgSet fMaster
the cached proposal data set
virtual void AddMapping(RooRealVar &proposalParam, RooAbsReal &update)
specify a mapping between a parameter of the proposal function and a parameter of interest.
virtual Double_t GetProposalDensity(RooArgSet &x1, RooArgSet &x2)
Return the probability of proposing the point x1 given the starting point x2.
virtual Bool_t Equals(RooArgSet &x1, RooArgSet &x2)
whether we own the proposal density function
RooDataSet * fCache
our position in the cached proposal data set
RooArgSet fLastX
pdf iterator
Bool_t fOwnsPdf
pointers to master variables needed for updates
Int_t fCacheSize
the last point we were at
Int_t fCachePosition
how many points to generate each time
std::map< RooRealVar *, RooAbsReal * > fMap
the proposal density function
virtual void Propose(RooArgSet &xPrime, RooArgSet &x)
Populate xPrime with a new proposed point.
PdfProposal()
By default, PdfProposal does NOT own the PDF that serves as the proposal density function.
ProposalFunction is an interface for all proposal functions that would be used with a Markov Chain Mo...
Iterator abstract base class.
Definition TIterator.h:30
virtual TObject * Next()=0
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:429
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)