Adaptive Mixture Metropolis (AMM)¶
Implementation of the Roberts and Rosenthal [75] adaptive (multivariate) mixture Metropolis [42][43][60] sampler for simulating autocorrelated draws from a distribution that can be specified up to a constant of proportionality.
Model-Based Constructor¶
-
AMM
(params::ElementOrVector{Symbol}, Sigma::Matrix{T<:Real}; adapt::Symbol=:all)¶ Construct a
Sampler
object for adaptive mixture Metropolis sampling. Parameters are assumed to be continuous, but may be constrained or unconstrained.Arguments
params
: stochastic node(s) to be updated with the sampler. Constrained parameters are mapped to unconstrained space according to transformations defined by the Stochasticunlist()
function.Sigma
: covariance matrix for the non-adaptive multivariate normal proposal distribution. The covariance matrix is relative to the unconstrained parameter space, where candidate draws are generated.adapt
: type of adaptation phase. Options are:all
: adapt proposal during all iterations.:burnin
: adapt proposal during burn-in iterations.:none
: no adaptation (multivariate Metropolis sampling with fixed proposal).
Value
Returns aSampler{AMMTune}
type object.Example
Stand-Alone Function¶
-
amm!
(v::AMMVariate, SigmaF::Cholesky{Float64}, logf::Function; adapt::Bool=true)¶ Simulate one draw from a target distribution using an adaptive mixture Metropolis sampler. Parameters are assumed to be continuous and unconstrained.
Arguments
v
: current state of parameters to be simulated. When running the sampler in adaptive mode, thev
argument in a successive call to the function should contain thetune
field returned by the previous call.SigmaF
: Cholesky factorization of the covariance matrix for the non-adaptive multivariate normal proposal distribution.logf
: function that takes a singleDenseVector
argument of parameter values at which to compute the log-transformed density (up to a normalizing constant).adapt
: whether to adaptively update the proposal distribution.
Value
Returnsv
updated with simulated values and associated tuning parameters.Example
The following example samples parameters in a simple linear regression model. Details of the model specification and posterior distribution can be found in the Supplement.
################################################################################ ## Linear Regression ## y ~ N(b0 + b1 * x, s2) ## b0, b1 ~ N(0, 1000) ## s2 ~ invgamma(0.001, 0.001) ################################################################################ using Mamba ## Data data = Dict( :x => [1, 2, 3, 4, 5], :y => [1, 3, 3, 3, 5] ) ## Log-transformed Posterior(b0, b1, log(s2)) + Constant logf = function(x::DenseVector) b0 = x[1] b1 = x[2] logs2 = x[3] r = data[:y] - b0 - b1 * data[:x] (-0.5 * length(data[:y]) - 0.001) * logs2 - (0.5 * dot(r, r) + 0.001) / exp(logs2) - 0.5 * b0^2 / 1000 - 0.5 * b1^2 / 1000 end ## MCMC Simulation with Adaptive Multivariate Metopolis Sampling n = 5000 burnin = 1000 sim = Chains(n, 3, names = ["b0", "b1", "s2"]) theta = AMMVariate([0.0, 0.0, 0.0]) SigmaF = cholfact(eye(3)) for i in 1:n amm!(theta, SigmaF, logf, adapt = (i <= burnin)) sim[i, :, 1] = [theta[1:2]; exp(theta[3])] end describe(sim)
AMMVariate Type¶
Declaration¶
typealias AMMVariate SamplerVariate{AMMTune}
Fields¶
value::Vector{Float64}
: simulated values.tune::AMMTune
: tuning parameters for the sampling algorithm.
Constructors¶
-
AMMVariate
(x::AbstractVector{T<:Real})¶ -
AMMVariate
(x::AbstractVector{T<:Real}, tune::AMMTune) Construct a
AMMVariate
object that stores simulated values and tuning parameters for adaptive mixture Metropolis sampling.Arguments
x
: simulated values.tune
: tuning parameters for the sampling algorithm. If not supplied, parameters are set to their defaults.
Value
Returns aAMMVariate
type object with fields set to the values supplied to argumentsx
andtune
.
AMMTune Type¶
Declaration¶
type AMMTune <: SamplerTune
Fields¶
adapt::Bool
: whether the proposal distribution has been adaptively tuned.beta::Real
: proportion of weight given to draws from the non-adaptive proposal with covariance factorizationSigmaF
, relative to draws from the adaptively tuned proposal with covariance factorizationSigmaLm
, during adaptive updating. Fixed atbeta = 0.05
.m::Int
: number of adaptive update iterations that have been performed.Mv::Vector{Float64}
: running mean of drawsv
during adaptive updating. Used in the calculation ofSigmaLm
.Mvv::Vector{Float64}
: running mean ofv * v'
during adaptive updating. Used in the calculation ofSigmaLm
.scale::Real
: fixed value2.38^2
in the factor (scale / length(v)
) by which the adaptively updated covariance matrix is scaled—adopted from Gelman, Roberts, and Gilks [30].SigmaF::Cholesky{Float64}
: factorization of the non-adaptive covariance matrix.SigmaLm::Matrix{Float64}
: lower-triangular factorization of the adaptively tuned covariance matrix.