// @(#)root/roostats:$Id$
// Authors: Kevin Belasco        17/06/2009
// Authors: Kyle Cranmer         17/06/2009
/*************************************************************************
 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

#ifndef ROOSTATS_MetropolisHastings
#define ROOSTATS_MetropolisHastings

#ifndef RooStats_RooStatsUtils
#include "RooStats/RooStatsUtils.h"
#endif
#ifndef ROOT_Rtypes
#include "Rtypes.h"
#endif
#ifndef ROOT_TObject
#include "TObject.h"
#endif
#ifndef ROO_ARG_SET
#include "RooArgSet.h"
#endif
#ifndef ROOSTATS_ProposalFunction
#include "RooStats/ProposalFunction.h"
#endif
#ifndef ROOSTATS_MarkovChain
#include "RooStats/MarkovChain.h"
#endif

namespace RooStats {

   class MetropolisHastings :  public TObject {


   public:
      enum FunctionSign {kNegative, kPositive, kSignUnset};
      enum FunctionType {kRegular, kLog, kTypeUnset};

      // default constructor
      MetropolisHastings();

      // alternate constructor
      MetropolisHastings(RooAbsReal& function, const RooArgSet& paramsOfInterest,
            ProposalFunction& proposalFunction, Int_t numIters);

      virtual ~MetropolisHastings() {}

      // main purpose of MetropolisHastings - run Metropolis-Hastings
      // algorithm to generate Markov Chain of points in the parameter space
      virtual MarkovChain* ConstructChain();

      // specify the parameters to store in the chain
      // if not specified all of them will be stored
      virtual void SetChainParameters(const RooArgSet& set)
      { fChainParams.removeAll();  fChainParams.add(set);  RemoveConstantParameters(&fChainParams); }      
      // specify all the parameters of interest in the interval
      virtual void SetParameters(const RooArgSet& set)
      { fParameters.removeAll();  fParameters.add(set);  RemoveConstantParameters(&fParameters); }
      // set the proposal function for suggesting new points for the MCMC
      virtual void SetProposalFunction(ProposalFunction& proposalFunction)
      { fPropFunc = &proposalFunction; }
      // set the number of iterations to run the metropolis algorithm
      virtual void SetNumIters(Int_t numIters)
      { fNumIters = numIters; }
      // set the number of steps in the chain to discard as burn-in,
      // starting from the first
      virtual void SetNumBurnInSteps(Int_t numBurnInSteps)
      { fNumBurnInSteps = numBurnInSteps; }
      // set the (likelihood) function
      virtual void SetFunction(RooAbsReal& function) { fFunction = &function; }
      // set the sign of the function
      virtual void SetSign(enum FunctionSign sign) { fSign = sign; }
      // set the type of the function
      virtual void SetType(enum FunctionType type) { fType = type; }


   protected:
      RooAbsReal* fFunction; // function that will generate likelihood values
      RooArgSet fParameters; // RooRealVars that define all parameter space
      RooArgSet fChainParams; // RooRealVars that are stored in the chain
      ProposalFunction* fPropFunc; // Proposal function for MCMC integration
      Int_t fNumIters; // number of iterations to run metropolis algorithm
      Int_t fNumBurnInSteps; // number of iterations to discard as burn-in, starting from the first
      enum FunctionSign fSign; // whether the likelihood is negative (like NLL) or positive
      enum FunctionType fType; // whether the likelihood is on a regular, log, (or other) scale

      // whether we should take the step, based on the value of d, fSign, fType
      virtual Bool_t ShouldTakeStep(Double_t d);
      virtual Double_t CalcNLL(Double_t xL);

      ClassDef(MetropolisHastings,2) // Markov Chain Monte Carlo calculator for Bayesian credible intervals
   };
}


#endif
 MetropolisHastings.h:1
 MetropolisHastings.h:2
 MetropolisHastings.h:3
 MetropolisHastings.h:4
 MetropolisHastings.h:5
 MetropolisHastings.h:6
 MetropolisHastings.h:7
 MetropolisHastings.h:8
 MetropolisHastings.h:9
 MetropolisHastings.h:10
 MetropolisHastings.h:11
 MetropolisHastings.h:12
 MetropolisHastings.h:13
 MetropolisHastings.h:14
 MetropolisHastings.h:15
 MetropolisHastings.h:16
 MetropolisHastings.h:17
 MetropolisHastings.h:18
 MetropolisHastings.h:19
 MetropolisHastings.h:20
 MetropolisHastings.h:21
 MetropolisHastings.h:22
 MetropolisHastings.h:23
 MetropolisHastings.h:24
 MetropolisHastings.h:25
 MetropolisHastings.h:26
 MetropolisHastings.h:27
 MetropolisHastings.h:28
 MetropolisHastings.h:29
 MetropolisHastings.h:30
 MetropolisHastings.h:31
 MetropolisHastings.h:32
 MetropolisHastings.h:33
 MetropolisHastings.h:34
 MetropolisHastings.h:35
 MetropolisHastings.h:36
 MetropolisHastings.h:37
 MetropolisHastings.h:38
 MetropolisHastings.h:39
 MetropolisHastings.h:40
 MetropolisHastings.h:41
 MetropolisHastings.h:42
 MetropolisHastings.h:43
 MetropolisHastings.h:44
 MetropolisHastings.h:45
 MetropolisHastings.h:46
 MetropolisHastings.h:47
 MetropolisHastings.h:48
 MetropolisHastings.h:49
 MetropolisHastings.h:50
 MetropolisHastings.h:51
 MetropolisHastings.h:52
 MetropolisHastings.h:53
 MetropolisHastings.h:54
 MetropolisHastings.h:55
 MetropolisHastings.h:56
 MetropolisHastings.h:57
 MetropolisHastings.h:58
 MetropolisHastings.h:59
 MetropolisHastings.h:60
 MetropolisHastings.h:61
 MetropolisHastings.h:62
 MetropolisHastings.h:63
 MetropolisHastings.h:64
 MetropolisHastings.h:65
 MetropolisHastings.h:66
 MetropolisHastings.h:67
 MetropolisHastings.h:68
 MetropolisHastings.h:69
 MetropolisHastings.h:70
 MetropolisHastings.h:71
 MetropolisHastings.h:72
 MetropolisHastings.h:73
 MetropolisHastings.h:74
 MetropolisHastings.h:75
 MetropolisHastings.h:76
 MetropolisHastings.h:77
 MetropolisHastings.h:78
 MetropolisHastings.h:79
 MetropolisHastings.h:80
 MetropolisHastings.h:81
 MetropolisHastings.h:82
 MetropolisHastings.h:83
 MetropolisHastings.h:84
 MetropolisHastings.h:85
 MetropolisHastings.h:86
 MetropolisHastings.h:87
 MetropolisHastings.h:88
 MetropolisHastings.h:89
 MetropolisHastings.h:90
 MetropolisHastings.h:91
 MetropolisHastings.h:92
 MetropolisHastings.h:93
 MetropolisHastings.h:94
 MetropolisHastings.h:95
 MetropolisHastings.h:96
 MetropolisHastings.h:97
 MetropolisHastings.h:98
 MetropolisHastings.h:99
 MetropolisHastings.h:100