Got a probability density function you want to draw samples from? Don't want to learn all the fancy stuff of the fancy sampler packages? The KissMCMC (Keep it simple, stupid, MCMC) package intends to provide a few simple MCMC samplers.
using KissMCMC
# the distribution to sample from,
logpdf(x::T) where {T} = x<0 ? -convert(T,Inf) : -x
# initial point of walker
theta0 = 0.5
# Metropolis MCMC sampler:
sample_prop_normal(theta) = 1.5*randn() + theta # samples the proposal (or jump) distribution
thetas, accept_ratio = metropolis(logpdf, sample_prop_normal, theta0, niter=10^5)
println("Accept ratio Metropolis: $accept_ratio")
# emcee MCMC sampler:
thetase, accept_ratioe = emcee(logpdf, make_theta0s(theta0, 0.1, logpdf, 100), niter=10^5)
# check convergence using integrated autocorrelation
thetase, accept_ratioe = squash_walkers(thetase, accept_ratioe) # puts all walkers into one
println("Accept ratio emcee: $accept_ratio")
using Plots
histogram(thetas, normalize=true, fillalpha=0.4)
histogram!(thetase, normalize=true, fillalpha=0.1)
plot!(0:0.01:5, map(x->exp(logpdf(x)[1]), 0:0.01:5), lw=3)
outputs:
MCMC samplers:
- Metropolis (serial)
metropolis
- Affine invariant MCMC, aka emcee
emcee
(threaded)
Other, probably better Julia MCMC packages:
- https://github.com/tpapp/DynamicHMC.jl
- https://github.com/madsjulia/AffineInvariantMCMC.jl
- https://github.com/brian-j-smith/Mamba.jl
- and many others
The (original) emcee python package: https://github.com/dfm/emcee