General Usage

This package implements the AbstractMCMC interface. AbstractMCMC provides a unifying interface for MCMC algorithms applied to LogDensityProblems.

Examples

Drawing Samples From a LogDensityProblems Through AbstractMCMC

SliceSampling.jl implements the AbstractMCMC interface through LogDensityProblems. That is, one simply needs to define a LogDensityProblems and pass it to AbstractMCMC:

using AbstractMCMC
using Distributions
using LinearAlgebra
using LogDensityProblems
using Plots

using SliceSampling

struct Target{D}
	dist::D
end

LogDensityProblems.logdensity(target::Target, x) = logpdf(target.dist, x)

LogDensityProblems.dimension(target::Target) = length(target.distx)

LogDensityProblems.capabilities(::Type{<:Target}) = LogDensityProblems.LogDensityOrder{0}()

sampler         = GibbsPolarSlice(2.0)
n_samples       = 10000
model           = Target(MvTDist(5, zeros(10), Matrix(I, 10, 10)))
logdensitymodel = AbstractMCMC.LogDensityModel(model)

chain   = sample(logdensitymodel, sampler, n_samples; initial_params=randn(10))
samples = hcat([transition.params for transition in chain]...)

plot(samples[1,:], xlabel="Iteration", ylabel="Trace")
savefig("abstractmcmc_demo.svg")

"/home/runner/work/SliceSampling.jl/SliceSampling.jl/docs/build/abstractmcmc_demo.svg"

Drawing Samples From Turing Models

SliceSampling.jl can also be used to sample from Turing models through Turing's externalsampler interface:

using Distributions
using SliceSampling
using StatsBase
using Turing

@model function demo()
    s ~ InverseGamma(3, 3)
    m ~ Normal(0, sqrt(s))
end

sampler   = RandPermGibbs(SliceSteppingOut(2.))
n_samples = 10000
model     = demo()
chain     = sample(model, externalsampler(sampler), n_samples; progress=false)
summarystats(chain)

╭─FlexiSummary (9 statistics) ─────────────────────────────────────────────────╮
│   iter    collapsed                                                          │
│   chain   collapsed                                                          │
│ ↓ stat  = [mean, std, mcse, ess_bulk, ess_tail, rhat, q5, q50, q95]          │
│                                                                              │
│ Parameters (2) ── AbstractPPL.VarName                                        │
│  Float64  s, m                                                               │
│                                                                              │
│ Extras (3)                                                                   │
│  Float64  logprior, loglikelihood, logjoint                                  │
│                                                                              │
│ Summary                                                                      │
│   param     mean     std    mcse   ess_bulk   ess_tail    rhat       q5  …   │
│       s   1.4835  1.3947  0.0201  5718.1330  5534.6658  0.9999   0.4741  …   │
│       m  -0.0042  1.2138  0.0126  9046.7750  5956.1879  1.0000  -1.9209  …   │
╰──────────────────────────────────────────────────────────────────────────────╯

Conditional sampling in a Turing.Gibbs sampler

SliceSampling.jl be used as a conditional sampler in Turing.Gibbs.

using Distributions
using SliceSampling
using StatsBase
using Turing

@model function simple_choice(xs)
    p ~ Beta(2, 2)
    z ~ Bernoulli(p)
    for i in 1:length(xs)
        if z == 1
            xs[i] ~ Normal(0, 1)
        else
            xs[i] ~ Normal(2, 1)
        end
    end
end

sampler = Turing.Gibbs(
    :p => externalsampler(SliceSteppingOut(2.0)),
    :z => PG(20),
)

n_samples = 1000
model     = simple_choice([1.5, 2.0, 0.3])
chain     = sample(model, sampler, n_samples; progress=false)
summarystats(chain)

╭─FlexiSummary (9 statistics) ─────────────────────────────────────────────────╮
│   iter    collapsed                                                          │
│   chain   collapsed                                                          │
│ ↓ stat  = [mean, std, mcse, ess_bulk, ess_tail, rhat, q5, q50, q95]          │
│                                                                              │
│ Parameters (2) ── AbstractPPL.VarName                                        │
│  Float64  p, z                                                               │
│                                                                              │
│ Extras (3)                                                                   │
│  Float64  logprior, loglikelihood, logjoint                                  │
│                                                                              │
│ Summary                                                                      │
│   param    mean     std    mcse  ess_bulk  ess_tail    rhat      q5  …       │
│       p  0.4399  0.2170  0.0076  793.2805  496.8385  1.0020  0.1002  …       │
│       z  0.1800  0.3844  0.0142  732.7185       NaN  0.9991  0.0000  …       │
╰──────────────────────────────────────────────────────────────────────────────╯

Drawing Samples

For drawing samples using the algorithms provided by SliceSampling, the user only needs to call:

sample([rng,] model, slice, N; initial_params)

• slice::AbstractSliceSampling: Any slice sampling algorithm provided by SliceSampling.
• model: A model implementing the LogDensityProblems interface.
• N: The number of samples

The output is a vector of SliceSampling.Transitions, which contains the following:

struct Transition

For the keyword arguments, SliceSampling allows:

• initial_params: The initial state of the Markov chain (default: nothing).

If initial_params is nothing, the following function can be implemented to provide an initialization:

initial_sample(rng, model)

Performing a Single Transition

For more fined-grained control, the user can call AbstractMCMC.step. That is, the chain can be initialized by calling:

transition, state = AbstractMCMC.steps([rng,] model, slice; initial_params)

and then each MCMC transition on state can be performed by calling:

transition, state = AbstractMCMC.steps([rng,] model, slice, state)

For more details, refer to the documentation of AbstractMCMC.