00001 // @(#)root/roostats:$Id: MetropolisHastings.h 31276 2009-11-18 15:06:42Z moneta $ 00002 // Authors: Kevin Belasco 17/06/2009 00003 // Authors: Kyle Cranmer 17/06/2009 00004 /************************************************************************* 00005 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. * 00006 * All rights reserved. * 00007 * * 00008 * For the licensing terms see $ROOTSYS/LICENSE. * 00009 * For the list of contributors see $ROOTSYS/README/CREDITS. * 00010 *************************************************************************/ 00011 00012 #ifndef ROOSTATS_MetropolisHastings 00013 #define ROOSTATS_MetropolisHastings 00014 00015 #ifndef RooStats_RooStatsUtils 00016 #include "RooStats/RooStatsUtils.h" 00017 #endif 00018 #ifndef ROOT_Rtypes 00019 #include "Rtypes.h" 00020 #endif 00021 #ifndef ROOT_TObject 00022 #include "TObject.h" 00023 #endif 00024 #ifndef ROO_ARG_SET 00025 #include "RooArgSet.h" 00026 #endif 00027 #ifndef ROOSTATS_ProposalFunction 00028 #include "RooStats/ProposalFunction.h" 00029 #endif 00030 #ifndef ROOSTATS_MarkovChain 00031 #include "RooStats/MarkovChain.h" 00032 #endif 00033 00034 namespace RooStats { 00035 00036 class MetropolisHastings : public TObject { 00037 00038 00039 public: 00040 enum FunctionSign {kNegative, kPositive, kSignUnset}; 00041 enum FunctionType {kRegular, kLog, kTypeUnset}; 00042 00043 // default constructor 00044 MetropolisHastings(); 00045 00046 // alternate constructor 00047 MetropolisHastings(RooAbsReal& function, RooArgSet& paramsOfInterest, 00048 ProposalFunction& proposalFunction, Int_t numIters); 00049 00050 virtual ~MetropolisHastings() {} 00051 00052 // main purpose of MetropolisHastings - run Metropolis-Hastings 00053 // algorithm to generate Markov Chain of points in the parameter space 00054 virtual MarkovChain* ConstructChain(); 00055 00056 // specify the parameters of interest in the interval 00057 // kbelasco: should clone before removing constant parameters? 00058 virtual void SetParameters(RooArgSet& set) 00059 { fParameters = &set; RemoveConstantParameters(fParameters); } 00060 // set the proposal function for suggesting new points for the MCMC 00061 virtual void SetProposalFunction(ProposalFunction& proposalFunction) 00062 { fPropFunc = &proposalFunction; } 00063 // set the number of iterations to run the metropolis algorithm 00064 virtual void SetNumIters(Int_t numIters) 00065 { fNumIters = numIters; } 00066 // set the number of steps in the chain to discard as burn-in, 00067 // starting from the first 00068 virtual void SetNumBurnInSteps(Int_t numBurnInSteps) 00069 { fNumBurnInSteps = numBurnInSteps; } 00070 // set the (likelihood) function 00071 virtual void SetFunction(RooAbsReal& function) { fFunction = &function; } 00072 // set the sign of the function 00073 virtual void SetSign(enum FunctionSign sign) { fSign = sign; } 00074 // set the type of the function 00075 virtual void SetType(enum FunctionType type) { fType = type; } 00076 00077 00078 protected: 00079 RooAbsReal* fFunction; // function that will generate likelihood values 00080 RooArgSet* fParameters; // RooRealVars that define parameter space 00081 ProposalFunction* fPropFunc; // Proposal function for MCMC integration 00082 Int_t fNumIters; // number of iterations to run metropolis algorithm 00083 Int_t fNumBurnInSteps; // number of iterations to discard as burn-in, starting from the first 00084 enum FunctionSign fSign; // whether the likelihood is negative (like NLL) or positive 00085 enum FunctionType fType; // whether the likelihood is on a regular, log, (or other) scale 00086 00087 // whether we should take the step, based on the value of d, fSign, fType 00088 virtual Bool_t ShouldTakeStep(Double_t d); 00089 virtual Double_t CalcNLL(Double_t xL); 00090 00091 ClassDef(MetropolisHastings,1) // Markov Chain Monte Carlo calculator for Bayesian credible intervals 00092 }; 00093 } 00094 00095 00096 #endif