Base class for MCMC fitting algorithms. More...

#include <Fitters/MCMCBase.h>

Inheritance diagram for MCMCBase:

Collaboration diagram for MCMCBase:

Public Member Functions
	MCMCBase (Manager *const fitMan)
	Constructor. More...

virtual	~MCMCBase ()=default
	Destructor. More...

void	RunMCMC () final
	Actual implementation of MCMC fitting algorithm. More...

void	StartFromPreviousFit (const std::string &FitName) final
	Allow to start from previous fit/chain. More...

void	setChainLength (unsigned int L)
	Set how long chain should be. More...

Public Member Functions inherited from FitterBase
	FitterBase (Manager *const fitMan)
	Constructor. More...

virtual	~FitterBase ()
	Destructor for the FitterBase class. More...

void	AddSampleHandler (SampleHandlerInterface *sample)
	This function adds a sample PDF object to the analysis framework. The sample PDF object will be utilized in fitting procedures or likelihood scans. More...

void	AddSystObj (ParameterHandlerBase *cov)
	This function adds a Covariance object to the analysis framework. The Covariance object will be utilized in fitting procedures or likelihood scans. More...

void	DragRace (const int NLaps=100)
	Calculates the required time for each sample or covariance object in a drag race simulation. Inspired by Dan's feature. More...

void	RunLLHScan ()
	Perform a 1D likelihood scan. More...

void	RunLLHMap ()
	Perform a general multi-dimensional likelihood scan. More...

void	GetStepScaleBasedOnLLHScan (const std::string &filename="")
	LLH scan is good first estimate of step scale. More...

void	Run2DLLHScan ()
	Perform a 2D likelihood scan. More...

void	RunSigmaVar ()
	Perform a 1D/2D sigma var for all samples. More...

std::string	GetName () const
	Get name of class. More...

Protected Member Functions
void	DoMCMCStep ()
	The full StartStep->DoStep->EndStep chain. More...

virtual void	ProposeStep ()=0
	Propose a step. More...

void	PreStepProcess ()
	Actions before step proposal [start stopwatch]. More...

void	PostStepProcess ()
	Actions after step proposal [end stopwatch, fill tree]. More...

virtual void	DoStep ()=0
	The MCMC step proposal and acceptance. More...

virtual double	AcceptanceProbability ()=0
	Step acceptance probability. More...

bool	IsStepAccepted (const double acc_prob)
	Is step accepted? More...

void	AcceptStep ()
	Accept a step. More...

void	AdaptiveStep ()
	Adaptive MCMC step. More...

void	PrintProgress (const bool StepsPrint=true)
	Print the progress. More...

void	CheckAcceptanceRates ()
	Check the acceptance rates. More...

Protected Member Functions inherited from FitterBase
void	ProcessMCMC ()
	Process MCMC output. More...

void	PrepareOutput ()
	Prepare the output file. More...

void	SaveOutput ()
	Save output and close files. More...

void	SanitiseInputs ()
	Remove obsolete memory and make other checks before fit starts. More...

void	SaveSettings ()
	Save the settings that the MCMC was run with. More...

bool	GetScanRange (std::map< std::string, std::vector< double >> &scanRanges) const
	YSP: Set up a mapping to store parameters with user-specified ranges, suggested by D. Barrow. More...

void	GetParameterScanRange (const ParameterHandlerBase *cov, const int i, double &CentralValue, double &lower, double &upper, const int n_points, const std::string &suffix="") const
	Helper function to get parameter scan range, central value. More...

bool	CheckSkipParameter (const std::vector< std::string > &SkipVector, const std::string &ParamName) const
	KS: Check whether we want to skip parameter using skip vector. More...

void	CustomRange (const std::string &ParName, const double sigma, double &ParamShiftValue) const
	For comparison with other fitting frameworks (like P-Theta) we usually have to apply different parameter values then usual 1, 3 sigma. More...

std::vector< double >	CalculateBinEdges (double lowerlimit, double upperlimit, int n_points) const

Protected Attributes
std::unique_ptr< MulticanonicalMCMCHandler >	multicanonicalHandler
	multicanonical handler for umbrella sampling More...

bool	out_of_bounds
	Do we reject based on hitting boundaries in systs. More...

bool	accept
	Accept. More...

unsigned int	chainLength
	number of steps in chain More...

bool	anneal
	simulated annealing More...

double	AnnealTemp
	simulated annealing temperature More...

bool	multicanonical
	multi-canonical method toggle on/off More...

Protected Attributes inherited from FitterBase
Manager *	fitMan
	The manager for configuration handling. More...

unsigned int	step
	current state More...

double	logLCurr
	current likelihood More...

double	logLProp
	proposed likelihood More...

double	accProb
	current acceptance prob More...

int	accCount
	counts accepted steps More...

unsigned int	stepStart
	step start, by default 0 if we start from previous chain then it will be different More...

std::vector< double >	sample_llh
	store the llh breakdowns More...

std::vector< double >	syst_llh
	systematic llh breakdowns More...

std::vector< SampleHandlerInterface * >	samples
	Sample holder. More...

unsigned int	TotalNSamples
	Total number of samples used, single SampleHandler can store more than one analysis sample! More...

std::vector< ParameterHandlerBase * >	systematics
	Systematic holder. More...

std::unique_ptr< TStopwatch >	clock
	tells global time how long fit took More...

std::unique_ptr< TStopwatch >	stepClock
	tells how long single step/fit iteration took More...

double	stepTime
	Time of single step. More...

std::unique_ptr< TRandom3 >	random
	Random number. More...

TFile *	outputFile
	Output. More...

TDirectory *	CovFolder
	Output cov folder. More...

TDirectory *	SampleFolder
	Output sample folder. More...

TTree *	outTree
	Output tree with posteriors. More...

int	auto_save
	auto save every N steps More...

bool	fTestLikelihood
	Necessary for some fitting algorithms like PSO. More...

bool	FileSaved
	Checks if file saved not repeat some operations. More...

bool	SettingsSaved
	Checks if setting saved not repeat some operations. More...

bool	OutputPrepared
	Checks if output prepared not repeat some operations. More...

std::string	AlgorithmName
	Name of fitting algorithm that is being used. More...

Detailed Description

Base class for MCMC fitting algorithms.

Inherits from FitterBase and defines the interface for MCMC-based fitting, including chain management and step handling.

Author: Asher Kaboth

Definition at line 8 of file MCMCBase.h.

Constructor & Destructor Documentation

◆ MCMCBase()

MCMCBase::MCMCBase ( Manager *const fitMan )

Constructor.

Parameters

fitMan A pointer to a manager object, which will handle all settings.

Definition at line 6 of file MCMCBase.cpp.

                                : FitterBase(man) {
 // *************************
     // Beginning step number
     stepStart = 0;
  
     // Starting parameters should be thrown
     out_of_bounds = false;
     chainLength = Get<unsigned>(fitMan->raw()["General"]["MCMC"]["NSteps"], __FILE__, __LINE__);
     if (chainLength < 10){
         MACH3LOG_ERROR("MCMC chain length must be at least 10 steps, otherwise this will result in a floating point exception.");
         throw MaCh3Exception(__FILE__, __LINE__);
     }
  
     AnnealTemp = GetFromManager<double>(fitMan->raw()["General"]["MCMC"]["AnnealTemp"], -999, __FILE__ , __LINE__);
     if (AnnealTemp < 0)
         anneal = false;
     else
     {
         MACH3LOG_INFO("Enabling simulated annealing with T = {}", AnnealTemp);
         anneal = true;
     }
 }

◆ ~MCMCBase()

virtual MCMCBase::~MCMCBase ( )

virtualdefault

Destructor.

Member Function Documentation

◆ AcceptanceProbability()

virtual double MCMCBase::AcceptanceProbability ( )

protectedpure virtual

Step acceptance probability.

Implemented in MR2T2, and DelayedMR2T2.

◆ AcceptStep()

void MCMCBase::AcceptStep ( )

protected

Accept a step.

Definition at line 251 of file MCMCBase.cpp.

                           {
 // *************************
     ++accCount;
     logLCurr = logLProp;
  
     // Loop over systematics and accept
     for (size_t s = 0; s < systematics.size(); ++s)
     {
         systematics[s]->AcceptStep();
     }
 }

◆ AdaptiveStep()

void MCMCBase::AdaptiveStep ( )

protected

Adaptive MCMC step.

Definition at line 221 of file MCMCBase.cpp.

                             {
 // *************************
     // Save the Adaptive output
     for (const auto &syst : systematics)
     {
         if (syst->GetDoAdaption()){
             syst->UpdateAdaptiveCovariance();
         }
     }
 }

◆ CheckAcceptanceRates()

void MCMCBase::CheckAcceptanceRates ( )

protected

Check the acceptance rates.

Definition at line 173 of file MCMCBase.cpp.

                                     {
 // *******************
     const auto StepEnd = stepStart + chainLength;
  
     // Skip check for reallllly early steps
     if(step - stepStart < std::min(10.0, static_cast<double>(StepEnd)/10)) {
         return;
     }
     // KS: Do not add "throw", this can crash CI or some short dummy chains for testing
     if(accCount==0 && step - stepStart > StepEnd/10) {
         MACH3LOG_CRITICAL("No steps were accepted in the MCMC chain after {} steps. Please check your  step sizes.", step - stepStart);
     }
  
     double acc_rate = static_cast<double>(accCount) / static_cast<double>(step - stepStart);
  
     if (acc_rate < 0.01) {
         MACH3LOG_WARN("Acceptance rate is very low: {:.4f}. This may indicate that the proposal distribution is not well-tuned. Consider reducing the step size.", acc_rate);
     } else if (acc_rate > 0.5) {
         MACH3LOG_WARN("Acceptance rate is very high: {:.4f}. This may indicate that the proposal distribution is too narrow. Consider increasing the step sizes.", acc_rate);
     }
 }

◆ DoMCMCStep()

void MCMCBase::DoMCMCStep ( )

protected

The full StartStep->DoStep->EndStep chain.

Starts step timer, prints progress

Step proposal, acceptance etc

Tree filling etc.

Definition at line 102 of file MCMCBase.cpp.

                           {
 // *******************
     PreStepProcess();
     DoStep();
     PostStepProcess();
 }

◆ DoStep()

virtual void MCMCBase::DoStep ( )

protectedpure virtual

The MCMC step proposal and acceptance.

Implemented in MR2T2, and DelayedMR2T2.

◆ IsStepAccepted()

bool MCMCBase::IsStepAccepted ( const double acc_prob )

protected

Is step accepted?

Parameters

acc_prob used for telling whether step is accepted or not

Definition at line 233 of file MCMCBase.cpp.

                                                    {
 // *************************
     // Get the random number
     const double fRandom = random->Rndm();
     // Do the accept/reject
     #ifdef MACH3_DEBUG
     debugFile << " logLProp: " << logLProp << " logLCurr: " << logLCurr << " acc_prob: " << acc_prob << " fRandom: " << fRandom << std::endl;
     #endif
  
     if (fRandom > acc_prob)
     {
         // Reject
         return false;
     }
     return true;
 }

◆ PostStepProcess()

void MCMCBase::PostStepProcess ( )

protected

Actions after step proposal [end stopwatch, fill tree].

Definition at line 127 of file MCMCBase.cpp.

                                {
 // *************************
     //KS: Some version of ROOT keep spamming about accessing already deleted object which is wrong and not helpful...
     int originalErrorLevel = gErrorIgnoreLevel;
     gErrorIgnoreLevel = kFatal;
  
     stepClock->Stop();
     stepTime = stepClock->RealTime();
  
     // Write step to output tree
     outTree->Fill();
  
     // Do Adaptive MCMC
     AdaptiveStep();
  
     if (step % auto_save == 0){
         outTree->AutoSave();
     }
     gErrorIgnoreLevel = originalErrorLevel;
 }

◆ PreStepProcess()

void MCMCBase::PreStepProcess ( )

protected

Actions before step proposal [start stopwatch].

Definition at line 113 of file MCMCBase.cpp.

                               {
 // *******************
     stepClock->Start();
     out_of_bounds = false;
  
     // Print 10 steps in total
     if ((step - stepStart) % (chainLength / 10) == 0)
     {
         CheckAcceptanceRates();
         PrintProgress();
     }
 }

◆ PrintProgress()

void MCMCBase::PrintProgress ( const bool StepsPrint = true )

protected

Print the progress.

Parameters

StepsPrint whether to print info about accepted steps and -LogL

Definition at line 150 of file MCMCBase.cpp.

                                             {
 // *******************
     if(StepsPrint) MACH3LOG_INFO("Step:\t{}/{}, current: {:.2f}, proposed: {:.2f}", step - stepStart, chainLength, logLCurr, logLProp);
     if(StepsPrint) MACH3LOG_INFO("Accepted/Total steps: {}/{} = {:.2f}", accCount, step - stepStart, static_cast<double>(accCount) / static_cast<double>(step - stepStart));
  
     for (size_t i = 0; i < samples.size(); ++i) {
         samples[i]->PrintRates();
     }
  
     for (ParameterHandlerBase *cov : systematics) {
         cov->PrintPreFitCurrPropValues();
     }
 #ifdef MACH3_DEBUG
     if (debug)
     {
         debugFile << "\n-------------------------------------------------------" << std::endl;
         debugFile << "Step:\t" << step + 1 << "/" << chainLength << "  |  current: " << logLCurr << " proposed: " << logLProp << std::endl;
     }
 #endif
 }

◆ ProposeStep()

virtual void MCMCBase::ProposeStep ( )

protectedpure virtual

Propose a step.

Implemented in MR2T2.

◆ RunMCMC()

void MCMCBase::RunMCMC ( )

finalvirtual

Actual implementation of MCMC fitting algorithm.

Initialise the multicanonical handler

Implements FitterBase.

Definition at line 32 of file MCMCBase.cpp.

                        {
 // *******************
     // Multicanonical method toggle from yaml config
     multicanonical = GetFromManager<bool>(fitMan->raw()["General"]["MCMC"]["Multicanonical"]["Enabled"], false, __FILE__, __LINE__);
     MACH3LOG_INFO("Multicanonical Method: {}", multicanonical);
  
     if (multicanonical) {
         multicanonicalHandler = std::make_unique<MulticanonicalMCMCHandler>();
     
         // Initialize the multicanonical handler with the systematics
         multicanonicalHandler->InitializeMulticanonicalHandlerConfig(fitMan, systematics);
         AlgorithmName += "_UmbrellaSampling"; // Append to the algorithm name
 #ifdef MACH3_DEBUG
     // Enable debug output stream for multicanonical handler if debug is enabled
     multicanonicalHandler->setDebugStream(&debugFile, debug);
 #endif
     }
  
     // Save the settings into the output file
     SaveSettings();
  
     // Prepare the output branches
     PrepareOutput();
  
     // Remove obsolete memory and make other checks before fit starts
     SanitiseInputs();
  
     // Print Progress before Propose Step
     PrintProgress(false);
  
     // Only propose step if we are running fresh chain. If we are running from previous chain then it's not needed and we can use usual pipeline
     if(stepStart == 0) {
         // Reconfigure the samples, systematics and oscillation for first weight
         // ProposeStep sets logLProp
         ProposeStep();
  
         // Initialise the value of the multicanonical parameter to the centres of the umbrellas
         if (multicanonical){
             multicanonicalHandler->InitializeMulticanonicalParams(systematics);
         }
         // Set the current logL to the proposed logL for the 0th step
         // Accept the first step to set logLCurr: this shouldn't affect the MCMC because we ignore the first N steps in burn-in
         logLCurr = logLProp;
     }
  
     // Begin MCMC
     const auto StepEnd = stepStart + chainLength;
     for (step = stepStart; step < StepEnd; ++step)
     {
         DoMCMCStep();
     }
     // Save all the MCMC output
     SaveOutput();
  
     // Process MCMC
     ProcessMCMC();
  
     // Save the adaptive MCMC
     for (const auto &syst : systematics)
     {
         if (syst->GetDoAdaption())
         {
             auto adaptive_handler = syst->GetAdaptiveHandler();
             adaptive_handler->SaveAdaptiveToFile(adaptive_handler->GetOutFileName(), syst->GetName(), true);
         }
     }
 }

◆ setChainLength()

void MCMCBase::setChainLength ( unsigned int L )

inline

Set how long chain should be.

Parameters

L	new chain length

Definition at line 26 of file MCMCBase.h.

26 { chainLength = L; };

◆ StartFromPreviousFit()

void MCMCBase::StartFromPreviousFit ( const std::string & FitName )

finalvirtual

Allow to start from previous fit/chain.

Parameters

FitName Name of previous chain

Reimplemented from FitterBase.

Definition at line 196 of file MCMCBase.cpp.

                                                             {
 // *******************
     // Use base class
     FitterBase::StartFromPreviousFit(FitName);
  
     // For MCMC we also need to set stepStart
     TFile *infile = M3::Open(FitName, "READ", __FILE__, __LINE__);
     TTree *posts = infile->Get<TTree>("posteriors");
     unsigned int step_val = 0;
  
     posts->SetBranchAddress("step", &step_val);
     posts->GetEntry(posts->GetEntries() - 1);
  
     stepStart = step_val;
     // KS: Also update number of steps if using adaption
     for (unsigned int i = 0; i < systematics.size(); ++i) {
         if (systematics[i]->GetDoAdaption()) {
             systematics[i]->SetNumberOfSteps(stepStart);
         }
     }
     infile->Close();
     delete infile;
 }

Member Data Documentation

◆ accept

bool MCMCBase::accept

protected

Accept.

Definition at line 70 of file MCMCBase.h.

◆ anneal

bool MCMCBase::anneal

protected

simulated annealing

Definition at line 76 of file MCMCBase.h.

◆ AnnealTemp

double MCMCBase::AnnealTemp

protected

simulated annealing temperature

Definition at line 78 of file MCMCBase.h.

◆ chainLength

unsigned int MCMCBase::chainLength

protected

number of steps in chain

Definition at line 73 of file MCMCBase.h.

◆ multicanonical

bool MCMCBase::multicanonical

protected

multi-canonical method toggle on/off

Definition at line 81 of file MCMCBase.h.

◆ multicanonicalHandler

std::unique_ptr<MulticanonicalMCMCHandler> MCMCBase::multicanonicalHandler

protected

multicanonical handler for umbrella sampling

Definition at line 64 of file MCMCBase.h.

◆ out_of_bounds

bool MCMCBase::out_of_bounds

protected

Do we reject based on hitting boundaries in systs.

Definition at line 67 of file MCMCBase.h.

The documentation for this class was generated from the following files:

Fitters/MCMCBase.h
Fitters/MCMCBase.cpp

Public Member Functions

Protected Member Functions

Protected Attributes

Detailed Description

Constructor & Destructor Documentation

◆ MCMCBase()

◆ ~MCMCBase()

Member Function Documentation

◆ AcceptanceProbability()

◆ AcceptStep()

◆ AdaptiveStep()

◆ CheckAcceptanceRates()

◆ DoMCMCStep()

◆ DoStep()

◆ IsStepAccepted()

◆ PostStepProcess()

◆ PreStepProcess()

◆ PrintProgress()

◆ ProposeStep()

◆ RunMCMC()

◆ setChainLength()

◆ StartFromPreviousFit()

Member Data Documentation

◆ accept

◆ anneal

◆ AnnealTemp

◆ chainLength

◆ multicanonical

◆ multicanonicalHandler

◆ out_of_bounds