Learn R Programming

RBaM (version 1.1.1)

MCMC_AM: Adaptive Metropolis sampler

Description

An adaptive Metropolis sampler largely inspired by Haario et al. (2001, https://www.jstor.org/stable/3318737). The jump covariance is adapted using the empirical covariance of previously-sampled values, and the scaling factor is adapted in order to comply with a specified move rate interval.

Usage

MCMC_AM(
  logPdf,
  x0,
  C0 = diag((0.01 * (abs(x0) + 0.1))^2),
  scaleFactor = 2.4/sqrt(length(x0)),
  nAdapt = 50,
  nCycles = 20,
  minMoveRate = 0.2,
  maxMoveRate = 0.5,
  downMult = 0.9,
  upMult = 1.1,
  burnCov = 0.2,
  dofCovMin = 10,
  nCovMax = 1000,
  ...
)

Value

A list with the following components:

samples

data frame, MCMC simulations.

components

data frame, corresponding values of the log-posterior, the log-prior and the log-likelihood.

C

matrix, the adapted jump covariance matrix.

scaleFactor

numeric, the adapted scaling factor.

Arguments

logPdf

function, evaluating the log-density of the distribution to sample from (up to a proportionality constant). logPdf can return either a single numeric value, interpreted as the target log-pdf, or a list containing components named 'logPosterior', 'logLikelihood' and 'logPrior'.

x0

numeric vector, starting point

C0

numeric matrix, covariance matrix of the Gaussian jump distribution (up to a scale factor, see next).

scaleFactor

numeric >0, used to scale the jump covariance. The covariance of the jump distribution is equal to (scaleFactor^2)*C0

nAdapt

integer > 1, number of iterations before adapting covariance C and scaleFactor.

nCycles

integer > 1, number of adaption cycles. Total number of iterations is hence equal to nAdapt*nCycles. nCycles=1 leads to the standard non-adaptive Metropolis sampler.

minMoveRate

numeric in (0;1), lower bound for the desired move rate interval.

maxMoveRate

numeric in (0;1), upper bound for the desired move rate interval.

downMult

numeric in (0;1), multiplicative factor used to decrease scaleFactor when move rate is too low.

upMult

numeric (>1, avoid 1/downMult) multiplicative factor used to increase scaleFactor when move rate is too high.

burnCov

numeric in (0;1), fraction of initial values to be discarded before computing the empirical covariance of sampled vectors, which is used to adapt the jump covariance.

dofCovMin

integer, minimum number of degrees of freedom required to compute the empirical covariance of sampled vectors and hence to adapt the jump covariance. If D denotes the length of x0, at least dofCovMin*(D+0.5*(D-1)*(D-2)) iterations are required before adapting the jump covariance (i.e. dofCovMin times the number of unknown elements in the covariance matrix).

nCovMax

integer, maximum number of iterations used to compute the empirical covariance. If the number of available iterations is larger than nCovMax, iterations are 'slimmed' to reach nCovMax.

...

other arguments passed to function logPdf

Examples

Run this code
# Define a 2-dimensional target log-pdf
logPdf <- function(x){
    p1=log(0.6*dnorm(x[1],0,1)+0.4*dnorm(x[1],2,0.5))
    p2=log(dlnorm(x[2],0,1))
    return(p1+p2)
}
# Sample from it
mcmc=MCMC_AM(logPdf,c(1,1))
plot(mcmc$samples)

Run the code above in your browser using DataLab