API: Turing.Inference

CSMC(...)

Equivalent to PG.

ESS

Elliptical slice sampling algorithm.

Examples

julia> @model function gdemo(x)
           m ~ Normal()
           x ~ Normal(m, 0.5)
       end
gdemo (generic function with 2 methods)

julia> sample(gdemo(1.0), ESS(), 1_000) |> mean
Mean

│ Row │ parameters │ mean     │
│     │ Symbol     │ Float64  │
├─────┼────────────┼──────────┤
│ 1   │ m          │ 0.824853 │

Emcee(n_walkers::Int, stretch_length=2.0)

Affine-invariant ensemble sampling algorithm.

Reference

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. (2013). emcee: The MCMC Hammer. Publications of the Astronomical Society of the Pacific, 125 (925), 306. https://doi.org/10.1086/670067

ExternalSampler{S<:AbstractSampler,AD<:ADTypes.AbstractADType,Unconstrained}

Represents a sampler that does not have a custom implementation of AbstractMCMC.step(rng, ::DynamicPPL.Model, spl).

The Unconstrained type-parameter is to indicate whether the sampler requires unconstrained space.

Fields

• sampler::AbstractMCMC.AbstractSampler: the sampler to wrap

• adtype::ADTypes.AbstractADType: the automatic differentiation (AD) backend to use

Turing.jl's interface for external samplers

If you implement a new MySampler <: AbstractSampler and want it to work with Turing.jl models, there are two options:

1. Directly implement the AbstractMCMC.step methods for DynamicPPL.Model. This is the most powerful option and is what Turing.jl's in-house samplers do. Implementing this means that you can directly call sample(model, MySampler(), N).

2. Implement a generic AbstractMCMC.step method for AbstractMCMC.LogDensityModel. This struct wraps an object that obeys the LogDensityProblems.jl interface, so your step implementation does not need to know anything about Turing.jl or DynamicPPL.jl. To use this with Turing.jl, you will need to wrap your sampler: sample(model, externalsampler(MySampler()), N).

This section describes the latter.

MySampler must implement the following methods:

• AbstractMCMC.step (the main function for taking a step in MCMC sampling; this is documented in AbstractMCMC.jl)
• Turing.Inference.getparams(::DynamicPPL.Model, external_transition): How to extract the parameters from the transition returned by your sampler (i.e., the first return value of step). There is a default implementation for this method, which is to return external_transition.θ.

There are a few more optional functions which you can implement to improve the integration with Turing.jl:

• Turing.Inference.isgibbscomponent(::MySampler): If you want your sampler to function as a component in Turing's Gibbs sampler, you should make this evaluate to true.

• Turing.Inference.requires_unconstrained_space(::MySampler): If your sampler requires unconstrained space, you should return true. This tells Turing to perform linking on the VarInfo before evaluation, and ensures that the parameter values passed to your sampler will always be in unconstrained (Euclidean) space.

Gibbs

A type representing a Gibbs sampler.

Constructors

Gibbs needs to be given a set of pairs of variable names and samplers. Instead of a single variable name per sampler, one can also give an iterable of variables, all of which are sampled by the same component sampler.

Each variable name can be given as either a Symbol or a VarName.

Some examples of valid constructors are:

Gibbs(:x => NUTS(), :y => MH())
Gibbs(@varname(x) => NUTS(), @varname(y) => MH())
Gibbs((@varname(x), :y) => NUTS(), :z => MH())

Fields

• varnames::NTuple{N, AbstractVector{<:AbstractPPL.VarName}} where N: varnames representing variables for each sampler

• samplers::NTuple{N, Any} where N: samplers for each entry in varnames

GibbsContext(target_varnames, global_varinfo, context)

A context used in the implementation of the Turing.jl Gibbs sampler.

There will be one GibbsContext for each iteration of a component sampler.

target_varnames is a a tuple of VarNames that the current component sampler is sampling. For those VarNames, GibbsContext will just pass tilde_assume!! calls to its child context. For other variables, their values will be fixed to the values they have in global_varinfo.

Fields

• target_varnames: the VarNames being sampled

• global_varinfo: a Ref to the global AbstractVarInfo object that holds values for all variables, both those fixed and those being sampled. We use a Ref because this field may need to be updated if new variables are introduced.

• context: the child context that tilde calls will eventually be passed onto.

HMC(ϵ::Float64, n_leapfrog::Int; adtype::ADTypes.AbstractADType = AutoForwardDiff())

Hamiltonian Monte Carlo sampler with static trajectory.

Arguments

• ϵ: The leapfrog step size to use.
• n_leapfrog: The number of leapfrog steps to use.
• adtype: The automatic differentiation (AD) backend. If not specified, ForwardDiff is used, with its chunksize automatically determined.

Usage

HMC(0.05, 10)

Tips

If you are receiving gradient errors when using HMC, try reducing the leapfrog step size ϵ, e.g.

# Original step size
sample(gdemo([1.5, 2]), HMC(0.1, 10), 1000)

# Reduced step size
sample(gdemo([1.5, 2]), HMC(0.01, 10), 1000)

HMCDA(
    n_adapts::Int, δ::Float64, λ::Float64; ϵ::Float64 = 0.0;
    adtype::ADTypes.AbstractADType = AutoForwardDiff(),
)

Hamiltonian Monte Carlo sampler with Dual Averaging algorithm.

Usage

HMCDA(200, 0.65, 0.3)

Arguments

• n_adapts: Numbers of samples to use for adaptation.
• δ: Target acceptance rate. 65% is often recommended.
• λ: Target leapfrog length.
• ϵ: Initial step size; 0 means automatically search by Turing.
• adtype: The automatic differentiation (AD) backend. If not specified, ForwardDiff is used, with its chunksize automatically determined.

Reference

For more information, please view the following paper (arXiv link):

Hoffman, Matthew D., and Andrew Gelman. "The No-U-turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo." Journal of Machine Learning Research 15, no. 1 (2014): 1593-1623.

IS()

Importance sampling algorithm.

Usage:

IS()

Example:

# Define a simple Normal model with unknown mean and variance.
@model function gdemo(x)
    s² ~ InverseGamma(2,3)
    m ~ Normal(0,sqrt.(s))
    x[1] ~ Normal(m, sqrt.(s))
    x[2] ~ Normal(m, sqrt.(s))
    return s², m
end

sample(gdemo([1.5, 2]), IS(), 1000)

MH(proposals...)

Construct a Metropolis-Hastings algorithm.

The arguments proposals can be

• Blank (i.e. MH()), in which case MH defaults to using the prior for each parameter as the proposal distribution.
• An iterable of pairs or tuples mapping a Symbol to a AdvancedMH.Proposal, Distribution, or Function that returns a conditional proposal distribution.
• A covariance matrix to use as for mean-zero multivariate normal proposals.

Examples

The default MH will draw proposal samples from the prior distribution using AdvancedMH.StaticProposal.

@model function gdemo(x, y)
    s² ~ InverseGamma(2,3)
    m ~ Normal(0, sqrt(s²))
    x ~ Normal(m, sqrt(s²))
    y ~ Normal(m, sqrt(s²))
end

chain = sample(gdemo(1.5, 2.0), MH(), 1_000)
mean(chain)

Specifying a single distribution implies the use of static MH:

# Use a static proposal for s² (which happens to be the same
# as the prior) and a static proposal for m (note that this
# isn't a random walk proposal).
chain = sample(
    gdemo(1.5, 2.0),
    MH(
        :s² => InverseGamma(2, 3),
        :m => Normal(0, 1)
    ),
    1_000
)
mean(chain)

Specifying explicit proposals using the AdvancedMH interface:

# Use a static proposal for s² and random walk with proposal
# standard deviation of 0.25 for m.
chain = sample(
    gdemo(1.5, 2.0),
    MH(
        :s² => AdvancedMH.StaticProposal(InverseGamma(2,3)),
        :m => AdvancedMH.RandomWalkProposal(Normal(0, 0.25))
    ),
    1_000
)
mean(chain)

Using a custom function to specify a conditional distribution:

# Use a static proposal for s and and a conditional proposal for m,
# where the proposal is centered around the current sample.
chain = sample(
    gdemo(1.5, 2.0),
    MH(
        :s² => InverseGamma(2, 3),
        :m => x -> Normal(x, 1)
    ),
    1_000
)
mean(chain)

Providing a covariance matrix will cause MH to perform random-walk sampling in the transformed space with proposals drawn from a multivariate normal distribution. The provided matrix must be positive semi-definite and square:

# Providing a custom variance-covariance matrix
chain = sample(
    gdemo(1.5, 2.0),
    MH(
        [0.25 0.05;
         0.05 0.50]
    ),
    1_000
)
mean(chain)

MHState(varinfo::AbstractVarInfo, logjoint_internal::Real)

State for Metropolis-Hastings sampling.

varinfo must have the correct parameters set inside it, but its other fields (e.g. accumulators, which track logp) can in general be missing or incorrect.

logjoint_internal is the log joint probability of the model, evaluated using the parameters and linking status of varinfo. It should be equal to DynamicPPL.getlogjoint_internal(varinfo). This information is returned by the MH sampler so we store this here to avoid re-evaluating the model unnecessarily.

NUTS(n_adapts::Int, δ::Float64; max_depth::Int=10, Δ_max::Float64=1000.0, init_ϵ::Float64=0.0; adtype::ADTypes.AbstractADType=AutoForwardDiff()

No-U-Turn Sampler (NUTS) sampler.

Usage:

NUTS()            # Use default NUTS configuration.
NUTS(1000, 0.65)  # Use 1000 adaption steps, and target accept ratio 0.65.

Arguments:

• n_adapts::Int : The number of samples to use with adaptation.
• δ::Float64 : Target acceptance rate for dual averaging.
• max_depth::Int : Maximum doubling tree depth.
• Δ_max::Float64 : Maximum divergence during doubling tree.
• init_ϵ::Float64 : Initial step size; 0 means automatically searching using a heuristic procedure.
• adtype::ADTypes.AbstractADType : The automatic differentiation (AD) backend. If not specified, ForwardDiff is used, with its chunksize automatically determined.

struct PG{R} <: Turing.Inference.ParticleInference

Particle Gibbs sampler.

Fields

• nparticles::Int64: Number of particles.

• resampler::Any: Resampling algorithm.

PG(n, [resampler = AdvancedPS.ResampleWithESSThreshold()]) PG(n, [resampler = AdvancedPS.resample_systematic, ]threshold)

Create a Particle Gibbs sampler of type PG with n particles.

If the algorithm for the resampling step is not specified explicitly, systematic resampling is performed if the estimated effective sample size per particle drops below 0.5.

PolynomialStepsize(a[, b=0, γ=0.55])

Create a polynomially decaying stepsize function.

At iteration t, the step size is

\[a (b + t)^{-γ}.\]

Prior()

Algorithm for sampling from the prior.

ProduceLogLikelihoodAccumulator{T<:Real} <: AbstractAccumulator

Exactly like LogLikelihoodAccumulator, but calls Libtask.produce on change of value.

Fields

• logp::Real: the scalar log likelihood value

RepeatSampler <: AbstractMCMC.AbstractSampler

A RepeatSampler is a container for a sampler and a number of times to repeat it.

Fields

• sampler: The sampler to repeat

• num_repeat: The number of times to repeat the sampler

Examples

repeated_sampler = RepeatSampler(sampler, 10)
AbstractMCMC.step(rng, model, repeated_sampler) # take 10 steps of `sampler`

SGHMC{AD}

Stochastic Gradient Hamiltonian Monte Carlo (SGHMC) sampler.

Fields

• learning_rate::Real

• momentum_decay::Real

• adtype::Any

Reference

Tianqi Chen, Emily Fox, & Carlos Guestrin (2014). Stochastic Gradient Hamiltonian Monte Carlo. In: Proceedings of the 31st International Conference on Machine Learning (pp. 1683–1691).

SGHMC(;
    learning_rate::Real,
    momentum_decay::Real,
    adtype::ADTypes.AbstractADType = AutoForwardDiff(),
)

Create a Stochastic Gradient Hamiltonian Monte Carlo (SGHMC) sampler.

If the automatic differentiation (AD) backend adtype is not provided, ForwardDiff with automatically determined chunksize is used.

Reference

Tianqi Chen, Emily Fox, & Carlos Guestrin (2014). Stochastic Gradient Hamiltonian Monte Carlo. In: Proceedings of the 31st International Conference on Machine Learning (pp. 1683–1691).

SGLD

Stochastic gradient Langevin dynamics (SGLD) sampler.

Fields

• stepsize::Any: Step size function.

• adtype::Any

Reference

Max Welling & Yee Whye Teh (2011). Bayesian Learning via Stochastic Gradient Langevin Dynamics. In: Proceedings of the 28th International Conference on Machine Learning (pp. 681–688).

SGLD(;
    stepsize = PolynomialStepsize(0.01),
    adtype::ADTypes.AbstractADType = AutoForwardDiff(),
)

Stochastic gradient Langevin dynamics (SGLD) sampler.

By default, a polynomially decaying stepsize is used.

If the automatic differentiation (AD) backend adtype is not provided, ForwardDiff with automatically determined chunksize is used.

Reference

Max Welling & Yee Whye Teh (2011). Bayesian Learning via Stochastic Gradient Langevin Dynamics. In: Proceedings of the 28th International Conference on Machine Learning (pp. 681–688).

See also: PolynomialStepsize

struct SMC{R} <: Turing.Inference.ParticleInference

Sequential Monte Carlo sampler.

Fields

• resampler::Any

SMC([resampler = AdvancedPS.ResampleWithESSThreshold()]) SMC([resampler = AdvancedPS.resample_systematic, ]threshold)

Create a sequential Monte Carlo sampler of type SMC.

If the algorithm for the resampling step is not specified explicitly, systematic resampling is performed if the estimated effective sample size per particle drops below 0.5.

_convert_initial_params(initial_params)

Convert initial_params to a DynamicPPl.AbstractInitStrategy if it is not already one, or throw a useful error message.

default_varinfo(rng, model, sampler)

Return a default varinfo object for the given model and sampler. The default method for this returns a NTVarInfo (i.e. 'typed varinfo').

dist_val_tuple(spl::MH, vi::VarInfo)

Return two NamedTuples.

The first NamedTuple has symbols as keys and distributions as values. The second NamedTuple has model symbols as keys and their stored values as values.

externalsampler(sampler::AbstractSampler; adtype=AutoForwardDiff(), unconstrained=true)

Wrap a sampler so it can be used as an inference algorithm.

Arguments

• sampler::AbstractSampler: The sampler to wrap.

Keyword Arguments

• adtype::ADTypes.AbstractADType=ADTypes.AutoForwardDiff(): The automatic differentiation (AD) backend to use.
• unconstrained::Bool=true: Whether the sampler requires unconstrained space.

gettracelocalresampledmaybe(fallback_resampled::Bool)

Get the Trace local resampled if one exists.

If executed within a TapedTask, return the resampled stored in the "taped globals" of the task, otherwise return fallback_resampled.

gettracelocalrngmaybe(rng::Random.AbstractRNG)

Get the Trace local rng if one exists.

If executed within a TapedTask, return the rng stored in the "taped globals" of the task, otherwise return vi.

gettracelocalvarinfomaybe(vi::AbstractVarInfo)

Get the Trace local varinfo if one exists.

If executed within a TapedTask, return the varinfo stored in the "taped globals" of the task, otherwise return vi.

Turing.Inference.getparams(model::DynamicPPL.Model, t::Any)

Return a vector of parameter values from the given sampler transition t (i.e., the first return value of AbstractMCMC.step). By default, returns the t.θ field.

removed in future releases and replaced with AbstractMCMC.getparams.

Turing.Inference.getparams(model::DynamicPPL.Model, t::AbstractVarInfo)

Return a key-value map of parameters from the varinfo.

Take the first step of MCMC for the first component sampler, and call the same function recursively on the remaining samplers, until no samplers remain. Return the global VarInfo and a tuple of initial states for all component samplers.

The step_function argument should always be either AbstractMCMC.step or AbstractMCMC.step_warmup.

Run a Gibbs step for the first varname/sampler/state tuple, and recursively call the same function on the tail, until there are no more samplers left.

The step_function argument should always be either AbstractMCMC.step or AbstractMCMC.step_warmup.

group_varnames_by_symbol(vns)

Group the varnames by their symbol.

Arguments

• vns: Iterable of VarName.

Returns

• OrderedDict{Symbol, Vector{VarName}}: A dictionary mapping symbol to a vector of varnames.

Turing.Inference.init_strategy(spl::AbstractSampler)

Get the default initialization strategy for a given sampler spl, i.e. how initial parameters for sampling are chosen if not specified by the user. By default, this is InitFromPrior(), which samples initial parameters from the prior distribution.

Initialise a VarInfo for the Gibbs sampler.

This is straight up copypasta from DynamicPPL's src/sampler.jl. It is repeated here to support calling both step and stepwarmup as the initial step. DynamicPPL initialstep is incompatible with stepwarmup.

isgibbscomponent(spl::AbstractSampler)

Return a boolean indicating whether spl is a valid component for a Gibbs sampler.

Defaults to false if no method has been defined for a particular algorithm type.

loadstate(chain::MCMCChains.Chains)

Load the final state of the sampler from a MCMCChains.Chains object.

To save the final state of the sampler, you must use sample(...; save_state=true). If this argument was not used during sampling, calling loadstate will throw an error.

make_conditional(model, target_variables, varinfo)

Return a new, conditioned model for a component of a Gibbs sampler.

Arguments

• model::DynamicPPL.Model: The model to condition.
• target_variables::AbstractVector{<:VarName}: The target variables of the component

sampler. These will not be conditioned.

• varinfo::DynamicPPL.AbstractVarInfo: Values for all variables in the model. All the

values in varinfo but not in target_variables will be conditioned to the values they have in varinfo.

Returns

• A new model with the variables not in target_variables conditioned.
• The GibbsContext object that will be used to condition the variables. This is necessary

because evaluation can mutate its global_varinfo field, which we need to access later.

match_linking!!(varinfo_local, prev_state_local, model)

Make sure the linked/invlinked status of varinfo_local matches that of the previous state for this sampler. This is relevant when multiple samplers are sampling the same variables, and one might need it to be linked while the other doesn't.

mh_accept(logp_current::Real, logp_proposal::Real, log_proposal_ratio::Real)

Decide if a proposal $x'$ with log probability $\log p(x') = logp_proposal$ and log proposal ratio $\log k(x', x) - \log k(x, x') = log_proposal_ratio$ in a Metropolis-Hastings algorithm with Markov kernel $k(x_t, x_{t+1})$ and current state $x$ with log probability $\log p(x) = logp_current$ is accepted by evaluating the Metropolis-Hastings acceptance criterion

\[\log U \leq \log p(x') - \log p(x) + \log k(x', x) - \log k(x, x')\]

for a uniform random number $U \in [0, 1)$.

requires_unconstrained_space(sampler::ExternalSampler)

Return true if the sampler requires unconstrained space, and false otherwise.

set_namedtuple!(vi::VarInfo, nt::NamedTuple)

Places the values of a NamedTuple into the relevant places of a VarInfo.

settracelocalvarinfomaybe(vi::AbstractVarInfo)

Set the Trace local varinfo if executing within a Trace. Return nothing.

If executed within a TapedTask, set the varinfo stored in the "taped globals" of the task. Otherwise do nothing.

setparams_varinfo!!(model, sampler::AbstractSampler, state, params::AbstractVarInfo)

A lot like AbstractMCMC.setparams!!, but instead of taking a vector of parameters, takes an AbstractVarInfo object. Also takes the sampler as an argument. By default, falls back to AbstractMCMC.setparams!!(model, state, params[:]).

model is typically a DynamicPPL.Model, but can also be e.g. an AbstractMCMC.LogDensityModel.