Metropolis Hastings Sampling
Metropolis Hastings sampling is a widely used Markov Chain Monte Carlo (MCMC) algorithm for generating a sequence of samples from a target probability distribution where direct sampling is difficult. The samples can be used to estimate properties of the target distribution, like its mean, variance, and higher moments. The samples can also be used to approximate the target distribution, for example via a histogram. The algorithm works by constructing a Markov chain with a stationary distribution equal to the target distribution.
Steps of the algorithm are as follows:
Choose an initial state for the Markov chain, usually from the target distribution.
At each iteration, propose a new state by sampling from a proposal distribution.
Calculate the acceptance probability.
Accept or reject the proposed state based on the acceptance probability.
Repeat steps 2-4 for a large number of iterations to obtain a sequence of samples.
Apply “burn-in” and “thining”.
Final samples from the Markov chain can then be used to estimate the target distribution or its properties.
MetropolisHastings Class
The MetropolisHastings
class is imported by:
from psimpy.sampler.metropolis_hastings import MetropolisHastings
Methods
- class MetropolisHastings(ndim, init_state, f_sample, target=None, ln_target=None, bounds=None, f_density=None, symmetric=True, nburn=0, nthin=1, seed=None, args_target=None, kwgs_target=None, args_f_sample=None, kwgs_f_sample=None, args_f_density=None, kwgs_f_density=None)[source]
Metropolis Hastings sampling.
- Parameters
ndim (int) – Dimension of parameter
x
.init_state (numpy array) – Contains
ndim
numeric values representing the starting point of the Metropolis hastings chain. Shape(ndim,)
.f_sample (Callable) – A function which proposes a new state
x'
given a current statex
. Call withf_sample(x, *args_f_sample, **kwgs_f_sample)
. Return the proposed statex'
as anumpy.ndarray
of shape(ndim,)
. If \(ndim=1\),f_sample
may also return a scalar value.target (Callable, optional) – Target probability density function, up to a normalizing constant. Call with
target(x, *args_target, **kwgs_target)
. Return the probability density value atx
. Iftarget
is not given,ln_target
must be provided.ln_target (Callable, optional) – Natural logarithm of target probability density function. Call with
ln_target(x, *args_target, **kwgs_target)
. Return the value of natural logarithm of target probability density value atx
. Ifln_target
is not given,target
must be provided.bounds (numpy array, optional) – Upper and lower boundaries of each parameter. Shape
(ndim, 2)
. bounds[:, 0] corresponds to lower boundaries of each parameter and bounds[:, 1] to upper boundaries of each parameter.f_density (Callable, optional) – Call with
f_density(x1, x2, *args_f_density, **kwgs_f_density)
. Return the probability density at statex1
given statex2
. It must be provided ifsymmetric
is set to False.symmetric (bool, optional) – Whether
f_density
is symmetric or asymmetric.nburn (int, optional) – Number of burnin. Default no burnin.
nthin (int, optional) – Number of thining to reduce dependency. Default no thining.
seed (int, optional) – Seed to initialize the pseudo-random number generator. Default None.
args_target (list, optional) – Positional arguments for
target
orln_target
.kwgs_target (dict, optional) – Keyword arguments for
target
orln_target
.args_f_sample (list, optional) – Positional arguments for
f_sample
.kwgs_f_sample (dict, optional) – Keyword arguments for
f_sample
.args_f_density (list, optional) – Positional arguments for
f_density
.kwgs_f_density (dict, optional) – Keyword arguments for
f_density
.
- sample(nsamples)[source]
Draw samples.
- Parameters
nsamples (int) – Number of samples to be drawn.
- Return type
tuple
[ndarray
,ndarray
]- Returns
mh_samples (numpy array) – Samples drawn from
target
. Shape(nsamples, ndim)
.mh_accept (numpy array) – Shape of
(nsamples,)
. Each element indicates whether the corresponding sample is the proposed new state (value 1) or the old state (value 0).np.mean(mh_accept)
thus gives the overall acceptance ratio.