Posterior Sampling Simulations

Class providing simulators for Bayesian posterior sampling—algorithms that draw samples from the joint posterior distribution of model parameters given observational constraints. Implementations orchestrate the Markov Chain Monte Carlo or other stochastic exploration strategy, managing the collection of walkers or chains, proposing new parameter vectors, evaluating the likelihood and prior, and deciding whether to accept or reject proposed moves. The default implementation uses the Differential Evolution MCMC algorithm for efficient high-dimensional sampling.

Default implementation: posteriorSampleSimulationDifferentialEvolution

Methods

simulate

Execute the full posterior sampling run, iterating the proposal, acceptance, and logging cycle until either the stopping criterion is met or the maximum step count is reached.

posteriorSampleSimulationAnnealedDffrntlEvltn

This class extends the differentialEvolution class to include an annealing schedule—the simulation begins at high temperature, waits for convergence, lowers the temperature and repeats until convergence at \(T=1\) is reached. In addition to the options for the differentialEvolution algorithm, the details of the algorithm are controlled by the following sub-parameters:

[temperatureMaximum]

The maximum temperature to use when tempering.

[temperatureLevels]

The number of temperature levels to use.

The temperature at level \(i\) is given by:

\[\log T_i = {i-1 \over N-1} \log T_\mathrm{max},\]

where \(T_\mathrm{max}=\)[temperatureMaximum] and \(N=\)[temperatureLevels].

Methods

initialize

Initialize the object.

• integer temperatureLevelCount [in]

• double precision temperatureMaximum [in]

Parameters

• [temperatureLevelCount] (integer; default 10) — The number temperature levels to use.

• [temperatureMaximum] (real) — The maximum temperature to reach.

• [stepsMaximum] (integer; default huge(0)) — The maximum number of steps to take. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [acceptanceAverageCount] (integer; default 10) — The number of steps over which to average the acceptance rate. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [stateSwapCount] (integer; default 10) — The number of steps between state swap steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [slowStepCount] (integer; default 1) — The number of steps between slow parameter update steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [recomputeCount] (integer; default -1) — The number of steps between forced recomputations of the likelihood. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [reportCount] (integer; default 10) — The number of steps between issuing reports. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [sampleOutliers] (boolean; default .true.) — If true, sample from outlier states. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [interactionRoot] (string; default none) — Root file name for interaction files, or ``none`’ if interaction is not required. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFileRoot] (string) — Root file name for log files. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logProposals] (boolean; default .false.) — If true, write a log of every proposed state to <logFileRoot>Proposals_<rank>.log, whether or not that proposal was accepted. The chain log records only the state retained at each step, so a rejected proposal’s parameters are otherwise never written anywhere. This matters for any analysis which pairs recorded model outputs (for example those written by posteriorSampleLikelihoodHaloMassFunction via pathSamples) with the parameters that produced them: without this log only accepted steps can be used, discarding the majority of the evaluations and, with them, the wider coverage of parameter space that rejected proposals provide. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [loadBalance] (boolean; default .true.) — If true, attempt to balance the workload across different compute nodes. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [ignoreChainNumberAdvice] (boolean; default .false.) — If true, ignore warnings and errors about not being able to span the full parameter space with the number of chains used. (inherited from posteriorSampleSimulationDifferentialEvolution)

posteriorSampleSimulationDifferentialEvolution

This class uses the differential evolution algorithm of Terr Braak (2006). Multiple, parallel chains are run and proposals are constructed by selecting two chains at random, taking a fraction, \(\gamma\), of the vector connecting the two chain states and adding this to the state of the current chain. The details of the algorithm are controlled by the following parameters:

[stepsMaximum]

The maximum number of steps to take.

[acceptanceAverageCount]

The number of steps over which to average the acceptance rate.

[stateSwapCount]

The number of steps after which to set \(\gamma=1\) to allow chains to swap states.

[logFileRoot]

The full path and root name of a file to log results to. The actual file name will have the rank of the MPI process appended to it.

[sampleOutliers]

If set to false then proposals for non-outlier chains post-convergence are constructed only from other non-outlier chains. Otherwise, proposals for non-outlier chains post-convergence are constructed from all other chains.

(Default implementation)

Methods

logging → logical

Return true if the simulator is currently logging state.

posterior

Return the log of posterior probability for the given posteriorSampleState.

• class(posteriorSampleStateClass) posteriorSampleState_ [inout]

• double precision logLikelihoodCurrent [in]

• double precision logPriorCurrent [in]

• double precision logPosterior [out]

• double precision logLikelihood [out]

• double precision logLikelihoodVariance [out]

• real timeEvaluate [inout]

• real timeEvaluatePrevious [in]

• logical forceAcceptance [inout]

update

Update the simulator to the new stateVector after a step.

• double precision(self%parametercount) stateVector [in]

temperature → double precision

Return the current temperature.

acceptProposal → logical

Return true if the proposed state should be accepted.

• double precision logPosterior [in]

• double precision logPosteriorProposed [in]

• double precision logLikelihoodVariance [in]

• double precision logLikelihoodVarianceProposed [in]

stepSize → double precision(self%parametercount)

Return the step size parameter, \(\gamma\), for the differential evolution proposal vector.

• logical forceAcceptance [inout]

chainSelect → integer

Select a chain.

• integer(:) blockedChains (optional) [in]

descriptorSpecial

Handle adding special parameters to the descriptor.

• type(inputParameters) descriptor [inout]

Parameters

• [stepsMaximum] (integer; default huge(0)) — The maximum number of steps to take.

• [acceptanceAverageCount] (integer; default 10) — The number of steps over which to average the acceptance rate.

• [stateSwapCount] (integer; default 10) — The number of steps between state swap steps.

• [slowStepCount] (integer; default 1) — The number of steps between slow parameter update steps.

• [recomputeCount] (integer; default -1) — The number of steps between forced recomputations of the likelihood.

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file.

• [reportCount] (integer; default 10) — The number of steps between issuing reports.

• [sampleOutliers] (boolean; default .true.) — If true, sample from outlier states.

• [interactionRoot] (string; default none) — Root file name for interaction files, or ``none`’ if interaction is not required.

• [logFileRoot] (string) — Root file name for log files.

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them.

• [logProposals] (boolean; default .false.) — If true, write a log of every proposed state to <logFileRoot>Proposals_<rank>.log, whether or not that proposal was accepted. The chain log records only the state retained at each step, so a rejected proposal’s parameters are otherwise never written anywhere. This matters for any analysis which pairs recorded model outputs (for example those written by posteriorSampleLikelihoodHaloMassFunction via pathSamples) with the parameters that produced them: without this log only accepted steps can be used, discarding the majority of the evaluations and, with them, the wider coverage of parameter space that rejected proposals provide.

• [loadBalance] (boolean; default .true.) — If true, attempt to balance the workload across different compute nodes.

• [ignoreChainNumberAdvice] (boolean; default .false.) — If true, ignore warnings and errors about not being able to span the full parameter space with the number of chains used.

posteriorSampleSimulationGrid

A posterior sampling simulation class which implements a simple grid search.

Methods

posterior

Return the log of posterior probability for the given simulationState.

• class(posteriorSampleStateClass) posteriorSampleState_ [inout]

• double precision logPosterior [out]

• double precision logLikelihood [out]

• real timeEvaluate [inout]

descriptorSpecial

Handle adding special parameters to the descriptor.

• type(inputParameters) descriptor [inout]

Parameters

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file.

• [logFileRoot] (string) — Root file name for log files.

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them.

• [outputLikelihoods] (boolean; default .false.) — If true, write likelihoods (and corresponding state vectors) to the output file.

posteriorSampleSimulationInitialState

A posterior sampling simulation class which implements evaluation of the model likelihood in the initial state only.

Methods

posterior

Return the log of posterior probability for the given simulationState.

• class(posteriorSampleStateClass) posteriorSampleState_ [inout]

• double precision logPosterior [out]

• double precision logLikelihood [out]

• real timeEvaluate [inout]

descriptorSpecial

Handle adding special parameters to the descriptor.

• type(inputParameters) descriptor [inout]

Parameters

• [logFileRoot] (string) — Root file name for log files.

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them.

posteriorSampleSimulationParticleSwarm

A posterior sampling simulation class which implements the particle swarm algorithm.

Methods

posterior

Return the log of posterior probability for the given simulationState.

• class(posteriorSampleStateClass) posteriorSampleState_ [inout]

• double precision logPosterior [out]

• double precision logLikelihood [out]

• double precision logLikelihoodVariance [out]

• real timeEvaluate [inout]

• real timeEvaluatePrevious [in]

• logical forceAcceptance [inout]

descriptorSpecial

Handle adding special parameters to the descriptor.

• type(inputParameters) descriptor [inout]

Parameters

• [stepsMaximum] (integer; default huge(0)) — The maximum number of steps to take.

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file.

• [reportCount] (integer; default 10) — The number of steps between issuing reports.

• [interactionRoot] (string; default none) — Root file name for interaction files, or ``none`’ if interaction is not required.

• [logFileRoot] (string) — Root file name for log files.

• [logFilePreviousRoot] (string; default none) — Root file name for log files from which to resume.

• [resume] (boolean; default .false.) — If true, resume from a previous set of log files.

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them.

• [inertiaWeight] (real; default 0.72d0) — The inertia weight \(w\) controlling how much of the particle’s current velocity is retained at each step; values less than 1 provide damping that aids convergence.

• [accelerationCoefficientPersonal] (real; default 1.193d0) — The personal (cognitive) acceleration coefficient \(c_1\) scaling the attraction of each particle toward its own historical best-known position.

• [accelerationCoefficientGlobal] (real; default 1.193d0) — The global (social) acceleration coefficient \(c_2\) scaling the attraction of each particle toward the swarm’s global best-known position.

• [velocityCoefficient] (real; default 0.5d0) — The maximum velocity coefficient \(v_\mathrm{max}\) expressed as a fraction of the parameter space range, limiting the step size of particle velocities during ongoing swarm iterations.

• [velocityCoefficientInitial] (real; default 0.0d0) — The initial velocity coefficient used to scale particle velocities at the start of the swarm simulation, before the main velocity update rule applies.

posteriorSampleSimulationStochasticDffrntlEvltn

This option extends the differentialEvolution option to run chains at a temperature matched to the uncertainty in the log-likelihood. This is designed to work with stochastic likelihood functions where an estimate of the uncertainty in the log-likelihood is available, and prevents the chains from becoming trapped in local maxima arising purely from random fluctuations. In addition to the options for the differentialEvolution algorithm, the details of the algorithm are controlled by the following parameters:

temperatureScale

The temperature scaling factor, \(alpha\), described below.

In computing the acceptance probability for transitions between states, the chain temperature is set to

\[T = 1 + \alpha C \sqrt{\sigma^2_\mathrm{current}+\sigma^2_\mathrm{proposed}},\]

where \(C = \log(\hat{R}/\hat{R}_\mathrm{t})\), \(\hat{R}\) is the current maximum (across all parameters) Gelman-Rubin convergence statistic, \(\hat{R}_\mathrm{t}\) is the target Gelman-Rubin convergence statistic at which convergence will be declared, and \(\sigma^2_\mathrm{current}\) and \(\sigma^2_\mathrm{proposed}\) are the variances in the log-likelihood of the current and proposed states respectively. This form ensures that the temperature declines to unity once the chains are converged (such that they will sample from the true posterior distribution), while ensuring that \(T\) is of order the size of the expected random fluctuations in the difference in log-likelihoods between states prior to chain convergence.

Methods

initialize

Initialize the object.

• double precision temperatureScale [in]

Parameters

• [temperatureScale] (real) — The temperature scale.

• [stepsMaximum] (integer; default huge(0)) — The maximum number of steps to take. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [acceptanceAverageCount] (integer; default 10) — The number of steps over which to average the acceptance rate. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [stateSwapCount] (integer; default 10) — The number of steps between state swap steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [slowStepCount] (integer; default 1) — The number of steps between slow parameter update steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [recomputeCount] (integer; default -1) — The number of steps between forced recomputations of the likelihood. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [reportCount] (integer; default 10) — The number of steps between issuing reports. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [sampleOutliers] (boolean; default .true.) — If true, sample from outlier states. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [interactionRoot] (string; default none) — Root file name for interaction files, or ``none`’ if interaction is not required. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFileRoot] (string) — Root file name for log files. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logProposals] (boolean; default .false.) — If true, write a log of every proposed state to <logFileRoot>Proposals_<rank>.log, whether or not that proposal was accepted. The chain log records only the state retained at each step, so a rejected proposal’s parameters are otherwise never written anywhere. This matters for any analysis which pairs recorded model outputs (for example those written by posteriorSampleLikelihoodHaloMassFunction via pathSamples) with the parameters that produced them: without this log only accepted steps can be used, discarding the majority of the evaluations and, with them, the wider coverage of parameter space that rejected proposals provide. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [loadBalance] (boolean; default .true.) — If true, attempt to balance the workload across different compute nodes. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [ignoreChainNumberAdvice] (boolean; default .false.) — If true, ignore warnings and errors about not being able to span the full parameter space with the number of chains used. (inherited from posteriorSampleSimulationDifferentialEvolution)

posteriorSampleSimulationTemperedDffrntlEvltn

This class extends the differentialEvolution option to include tempering during which the likelihood function is heated up and cooled down to allow chains to more easily walk through the likelihood landscape. In addition to the options for the differentialEvolution algorithm, the details of the algorithm are controlled by the following sub-parameters:

[untemperedStepCount]

The number of untempered (i.e. \(T=1\)) steps to take between tempering cycles.

[temperatureMaximum]

The maximum temperature to use when tempering.

[temperedLevels]

The number of tempered levels to use.

[stepsPerLevel]

The number of differential evolution steps to take at each tempering level.

[logFlushCount]

The number of steps after which the log file will be flushed to disk.

In each tempering cycle, the temperature is raised through levels \(1\)…\(N\) (where \(N=\)temperedLevels), and then back down through levels \(N-1\)…\(1\). The temperature at level \(i\) is given by:

\[\log T_i = {i \over N} \log T_\mathrm{max},\]

where \(T_\mathrm{max}=\)temperatureMaximum. During tempered steps, the \(\gamma\) parameter of the differential evolution algorithm is increased by a factor \(T^\alpha\), where \(\alpha\) is provided by the proposalSizeTemperatureExponent class. A value of \(\alpha=1/2\) is optimal for a Gaussian likelihood.

Methods

initialize

Return the current tempering level.

• class(posteriorSampleDffrntlEvltnPrpslSzTmpExpClass) posteriorSampleDffrntlEvltnPrpslSzTmpExp_ [in]

• integer temperingLevelCount [in]

• integer untemperedStepCount [in]

• integer stepsPerLevel [in]

• double precision temperatureMaximum [in]

level → integer

Return the current tempering level.

Parameters

• [untemperedStepCount] (integer; default 10) — The number of untempered steps to take.

• [stepsPerLevel] (integer; default 10) — The number of steps to take at each tempering level.

• [temperingLevelCount] (integer; default 10) — The number tempering levels to use.

• [temperatureMaximum] (real) — The maximum temperature to reach.

• [stepsMaximum] (integer; default huge(0)) — The maximum number of steps to take. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [acceptanceAverageCount] (integer; default 10) — The number of steps over which to average the acceptance rate. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [stateSwapCount] (integer; default 10) — The number of steps between state swap steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [slowStepCount] (integer; default 1) — The number of steps between slow parameter update steps. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [recomputeCount] (integer; default -1) — The number of steps between forced recomputations of the likelihood. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFlushCount] (integer; default 10) — The number of steps between flushing the log file. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [reportCount] (integer; default 10) — The number of steps between issuing reports. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [sampleOutliers] (boolean; default .true.) — If true, sample from outlier states. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [interactionRoot] (string; default none) — Root file name for interaction files, or ``none`’ if interaction is not required. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logFileRoot] (string) — Root file name for log files. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [appendLogs] (boolean; default .false.) — If true, do not overwrite existing log files, but instead append to them. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [logProposals] (boolean; default .false.) — If true, write a log of every proposed state to <logFileRoot>Proposals_<rank>.log, whether or not that proposal was accepted. The chain log records only the state retained at each step, so a rejected proposal’s parameters are otherwise never written anywhere. This matters for any analysis which pairs recorded model outputs (for example those written by posteriorSampleLikelihoodHaloMassFunction via pathSamples) with the parameters that produced them: without this log only accepted steps can be used, discarding the majority of the evaluations and, with them, the wider coverage of parameter space that rejected proposals provide. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [loadBalance] (boolean; default .true.) — If true, attempt to balance the workload across different compute nodes. (inherited from posteriorSampleSimulationDifferentialEvolution)

• [ignoreChainNumberAdvice] (boolean; default .false.) — If true, ignore warnings and errors about not being able to span the full parameter space with the number of chains used. (inherited from posteriorSampleSimulationDifferentialEvolution)