Logo ROOT   6.08/07
Reference Guide
MetropolisHastings.h
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 #ifndef ROOSTATS_MetropolisHastings
13 #define ROOSTATS_MetropolisHastings
14 
15 #ifndef RooStats_RooStatsUtils
16 #include "RooStats/RooStatsUtils.h"
17 #endif
18 #ifndef ROOT_Rtypes
19 #include "Rtypes.h"
20 #endif
21 #ifndef ROOT_TObject
22 #include "TObject.h"
23 #endif
24 #ifndef ROO_ARG_SET
25 #include "RooArgSet.h"
26 #endif
27 #ifndef ROOSTATS_ProposalFunction
29 #endif
30 #ifndef ROOSTATS_MarkovChain
31 #include "RooStats/MarkovChain.h"
32 #endif
33 
34 namespace RooStats {
35 
36 /**
37 
38  \ingroup Roostats
39 
40  This class uses the Metropolis-Hastings algorithm to construct a Markov Chain
41  of data points using Monte Carlo. In the main algorithm, new points in the
42  parameter space are proposed and then visited based on their relative
43  likelihoods. This class can use any implementation of the ProposalFunction,
44  including non-symmetric proposal functions, to propose parameter points and
45  still maintain detailed balance when constructing the chain.
46 
47 
48 
49  The "Likelihood" function that is sampled when deciding what steps to take in
50  the chain has been given a very generic implementation. The user can create
51  any RooAbsReal based on the parameters and pass it to a MetropolisHastings
52  object with the method SetFunction(RooAbsReal&). Be sure to tell
53  MetropolisHastings whether your RooAbsReal is on a (+/-) regular or log scale,
54  so that it knows what logic to use when sampling your RooAbsReal. For example,
55  a common use is to sample from a -log(Likelihood) distribution (NLL), for which
56  the appropriate configuration calls are SetType(MetropolisHastings::kLog);
57  SetSign(MetropolisHastings::kNegative);
58  If you're using a traditional likelihood function:
59  SetType(MetropolisHastings::kRegular); SetSign(MetropolisHastings::kPositive);
60  You must set these type and sign flags or MetropolisHastings will not construct
61  a MarkovChain.
62 
63  Also note that in ConstructChain(), the values of the variables are randomized
64  uniformly over their intervals before construction of the MarkovChain begins.
65 
66 */
67 
68 
69  class MetropolisHastings : public TObject {
70 
71 
72  public:
75 
76  // default constructor
78 
79  // alternate constructor
80  MetropolisHastings(RooAbsReal& function, const RooArgSet& paramsOfInterest,
81  ProposalFunction& proposalFunction, Int_t numIters);
82 
83  virtual ~MetropolisHastings() {}
84 
85  // main purpose of MetropolisHastings - run Metropolis-Hastings
86  // algorithm to generate Markov Chain of points in the parameter space
87  virtual MarkovChain* ConstructChain();
88 
89  // specify the parameters to store in the chain
90  // if not specified all of them will be stored
91  virtual void SetChainParameters(const RooArgSet& set)
93  // specify all the parameters of interest in the interval
94  virtual void SetParameters(const RooArgSet& set)
96  // set the proposal function for suggesting new points for the MCMC
97  virtual void SetProposalFunction(ProposalFunction& proposalFunction)
98  { fPropFunc = &proposalFunction; }
99  // set the number of iterations to run the metropolis algorithm
100  virtual void SetNumIters(Int_t numIters)
101  { fNumIters = numIters; }
102  // set the number of steps in the chain to discard as burn-in,
103  // starting from the first
104  virtual void SetNumBurnInSteps(Int_t numBurnInSteps)
105  { fNumBurnInSteps = numBurnInSteps; }
106  // set the (likelihood) function
107  virtual void SetFunction(RooAbsReal& function) { fFunction = &function; }
108  // set the sign of the function
109  virtual void SetSign(enum FunctionSign sign) { fSign = sign; }
110  // set the type of the function
111  virtual void SetType(enum FunctionType type) { fType = type; }
112 
113 
114  protected:
115  RooAbsReal* fFunction; // function that will generate likelihood values
116  RooArgSet fParameters; // RooRealVars that define all parameter space
117  RooArgSet fChainParams; // RooRealVars that are stored in the chain
118  ProposalFunction* fPropFunc; // Proposal function for MCMC integration
119  Int_t fNumIters; // number of iterations to run metropolis algorithm
120  Int_t fNumBurnInSteps; // number of iterations to discard as burn-in, starting from the first
121  enum FunctionSign fSign; // whether the likelihood is negative (like NLL) or positive
122  enum FunctionType fType; // whether the likelihood is on a regular, log, (or other) scale
123 
124  // whether we should take the step, based on the value of d, fSign, fType
125  virtual Bool_t ShouldTakeStep(Double_t d);
126  virtual Double_t CalcNLL(Double_t xL);
127 
128  ClassDef(MetropolisHastings,2) // Markov Chain Monte Carlo calculator for Bayesian credible intervals
129  };
130 }
131 
132 
133 #endif
virtual void SetSign(enum FunctionSign sign)
ProposalFunction is an interface for all proposal functions that would be used with a Markov Chain Mo...
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:86
virtual void SetType(enum FunctionType type)
virtual MarkovChain * ConstructChain()
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
This class uses the Metropolis-Hastings algorithm to construct a Markov Chain of data points using Mo...
virtual void SetParameters(const RooArgSet &set)
virtual void SetNumIters(Int_t numIters)
#define ClassDef(name, id)
Definition: Rtypes.h:254
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
virtual void SetFunction(RooAbsReal &function)
virtual void SetProposalFunction(ProposalFunction &proposalFunction)
virtual void SetChainParameters(const RooArgSet &set)
virtual Double_t CalcNLL(Double_t xL)
Stores the steps in a Markov Chain of points.
Definition: MarkovChain.h:53
Namespace for the RooStats classes.
Definition: Asimov.h:20
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
int type
Definition: TGX11.cxx:120
virtual void SetNumBurnInSteps(Int_t numBurnInSteps)
Mother of all ROOT objects.
Definition: TObject.h:37
void RemoveConstantParameters(RooArgSet *set)
Definition: RooStatsUtils.h:73
virtual Bool_t ShouldTakeStep(Double_t d)