DreamMCMC

class pints.DreamMCMC(chains, x0, sigma0=None)

Uses differential evolution adaptive Metropolis (DREAM) MCMC as described in [1] to perform posterior sampling from the posterior.

In each step of the algorithm N chains are evolved using the following steps:

1. Select proposal:

   x_proposed = x[i,r] + (1 + e) * gamma(delta, d, p_g) *
                sum_j=1^delta (X[i,r1[j]] - x[i,r2[j]])
                + epsilon
   

where [r1[j], r2[j]] are random chain indices chosen (without replacement) from the N available chains, which must not equal each other or i, where i indicates the current time step; delta ~ uniform_discrete(1,D) determines the number of terms to include in the summation:

e ~ U(-b*, b*) in d dimensions;
gamma(delta, d, p_g) =
  if p_g < u1 ~ U(0,1):
    2.38 / sqrt(2 * delta * d)
  else:
    1

epsilon ~ N(0,b) in d dimensions (where d is the dimensionality of the parameter vector).

2. Modify random subsets of the proposal according to a crossover probability CR:

for j in 1:N:
  if 1 - CR > u2 ~ U(0,1):
    x_proposed[j] = x[j]
  else:
    x_proposed[j] = x_proposed[j] from 1

If x_proposed / x[i,r] > u ~ U(0,1), then x[i+1,r] = x_proposed; otherwise, x[i+1,r] = x[i].

Here b > 0, b* > 0, 1 >= p_g >= 0, 1 >= CR >= 0.

Extends MultiChainMCMC, does not use sensitivities.

References

[1]

“Accelerating Markov Chain Monte Carlo Simulation by Differential Evolution with Self-Adaptive Randomized Subspace Sampling”, 2009, Vrugt et al., International Journal of Nonlinear Sciences and Numerical Simulation. https://doi.org/10.1515/IJNSNS.2009.10.3.273

CR()

Returns the probability of crossover occurring if constant crossover mode is enabled (see set_CR()).

ask()

See pints.MultiChainMCMC.ask().

b()

Returns the Gaussian scale coefficient used in updating the position of each chain.

b_star()

Returns b*, which determines the weight given to other chains’ positions in determining new positions (see set_b_star()).

constant_crossover()

Returns True if constant crossover mode is enabled.

current_log_pdfs()

Returns the log pdf values of the current points (i.e. of the most recent points returned by tell()).

delta_max()

Returns the maximum number of other chains’ positions to use to determine the next sampler position (see set_delta_max()).

in_initial_phase()

See pints.MCMCSampler.in_initial_phase().

nCR()

Returns the size of the discrete crossover probability distribution (only used if constant crossover mode is disabled), see set_nCR().

n_hyper_parameters()

See TunableMethod.n_hyper_parameters().

classmethod name()

See pints.MCMCSampler.name().

needs_initial_phase()

See pints.MCMCSampler.needs_initial_phase().

needs_sensitivities()

Returns True if this methods needs sensitivities to be passed in to tell along with the evaluated logpdf.

p_g()

Returns p_g. See set_p_g().

set_CR(CR)

Sets the probability of crossover occurring if constant crossover mode is enabled. CR is a probability and so must be in the range [0, 1].

set_b(b)

Sets the Gaussian scale coefficient used in updating the position of each chain (must be non-negative).

set_b_star(b_star)

Sets b*, which determines the weight given to other chains’ positions in determining new positions (must be non-negative).

set_constant_crossover(enabled)

Enables/disables constant-crossover mode (must be bool).

set_delta_max(delta_max)

Sets the maximum number of other chains’ positions to use to determine the next sampler position. delta_max must be in the range [1, nchains - 2].

set_hyper_parameters(x)

The hyper-parameter vector is [b, b_star, p_g, delta_max, initial_phase, constant_crossover, CR, nCR].

See TunableMethod.set_hyper_parameters().

set_initial_phase(initial_phase)

See pints.MCMCSampler.needs_initial_phase().

set_nCR(nCR)

Sets the size of the discrete crossover probability distribution (only used if constant crossover mode is disabled). nCR must be greater than or equal to 2.

set_p_g(p_g)

Sets p_g which is the probability of choosing a higher gamma versus regular (a higher gamma means that other chains are given more weight). p_g must be in the range [0, 1].

tell(proposed_log_pdfs)

See pints.MultiChainMCMC.tell().