API Documentation
This is the stable public API of BAT. Forward/backward compatibility follows Julia's semantic versioning rules.
Types
BAT.AbstractMCMCWeightingScheme
BAT.AbstractMedianEstimator
BAT.AbstractModeEstimator
BAT.AbstractPosteriorMeasure
BAT.AbstractSamplingAlgorithm
BAT.AbstractTransformTarget
BAT.AdaptiveMHTuning
BAT.AssumeConvergence
BAT.AutocorLenAlgorithm
BAT.BATContext
BAT.BATHDF5IO
BAT.BATIOAlgorithm
BAT.BinningAlgorithm
BAT.BrooksGelmanConvergence
BAT.ConvergenceTest
BAT.CuhreIntegration
BAT.DensitySample
BAT.DensitySampleVector
BAT.DivonneIntegration
BAT.DoNotTransform
BAT.EffSampleSizeAlgorithm
BAT.EffSampleSizeFromAC
BAT.EvaluatedMeasure
BAT.ExplicitInit
BAT.FixedNBins
BAT.FreedmanDiaconisBinning
BAT.GelmanRubinConvergence
BAT.GeyerAutocorLen
BAT.HamiltonianMC
BAT.IIDSampling
BAT.IdentityTransformAlgorithm
BAT.InitFromIID
BAT.InitFromSamples
BAT.InitFromTarget
BAT.InitvalAlgorithm
BAT.IntegrationAlgorithm
BAT.KishESS
BAT.MCMCAlgorithm
BAT.MCMCBurninAlgorithm
BAT.MCMCChainPoolInit
BAT.MCMCInitAlgorithm
BAT.MCMCMultiCycleBurnin
BAT.MCMCNoOpTuning
BAT.MCMCSampling
BAT.MCMCTuningAlgorithm
BAT.MHProposalDistTuning
BAT.MaxDensitySearch
BAT.MetropolisHastings
BAT.ModeAsDefined
BAT.OptimAlg
BAT.OptimizationAlg
BAT.OrderedResampling
BAT.PosteriorMeasure
BAT.PriorSubstitution
BAT.PriorToGaussian
BAT.PriorToUniform
BAT.RandResampling
BAT.RepetitionWeighting
BAT.RiceBinning
BAT.SampleMedianEstimator
BAT.ScottBinning
BAT.SokalAutocorLen
BAT.SquareRootBinning
BAT.SturgesBinning
BAT.SuaveIntegration
BAT.ToRealVector
BAT.TransformAlgorithm
BAT.VEGASIntegration
Functions and macros
BAT.bat_convergence
BAT.bat_default
BAT.bat_eff_sample_size
BAT.bat_findmedian
BAT.bat_findmode
BAT.bat_initval
BAT.bat_integrate
BAT.bat_read
BAT.bat_report
BAT.bat_sample
BAT.bat_transform
BAT.bat_write
BAT.distbind
BAT.distprod
BAT.get_batcontext
BAT.lbqintegral
BAT.set_batcontext
Documentation
BAT.bat_convergence
— Functionbat_convergence(
algoutput,
[algorithm::ConvergenceTest],
[context::BATContext]
)
Check if an algorithm has converged, based on it's output algoutput
Returns a NamedTuple of the shape
(result, ...)
result
indicates whether algoutput
and must either be a Bool
or support convert(Bool, result)
. It should typically contains measures of algorithm convergence, like a convergence value and it's threshold, etc.
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_convergence
, add methods to bat_convergence_impl
instead.
BAT.bat_default
— Functionbat_default(f::Base.Callable, argname::Symbol, objectives...)
bat_default(f::Base.Callable, argname::Val, objectives...)
Get the default value for argument argname
of function f
to use for objective
(s).
objective
(s) are mandatory arguments of function f
that semantically constitute it's main objective(s), and that that a good default choice of optional arguments (e.g. choice of algorithm(s), etc.) may depend on. Which arguments are considered to be objectives is function-specific.
For example:
bat_default(bat_sample, :algorithm, density::PosteriorMeasure) == MetropolisHastings()
bat_default(bat_sample, Val(:algorithm), samples::DensitySampleVector) == OrderedResampling()
BAT.bat_eff_sample_size
— Functionbat_eff_sample_size(
v::Union{AbstractVector{<:Real},AbstractVectorOfSimilarVectors{<:Real}},
[algorithm::EffSampleSizeAlgorithm],
[context::BATContext]
)
bat_eff_sample_size(
smpls::DensitySampleVector,
[algorithm::EffSampleSizeAlgorithm],
[context::BATContext]
)
Estimate effective sample size estimation for variate series v
, resp. density samples smpls
, separately for each degree of freedom.
Returns a NamedTuple of the shape
(result = eff_sample_size, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_eff_sample_size
, add methods to bat_eff_sample_size_impl
instead.
BAT.bat_findmedian
— Functionbat_findmedian(
samples::DensitySampleVector
)
The function computes the median of marginalized samples
.
Returns a NamedTuple of the shape
(result = v, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_findmedian
, add methods to bat_findmedian_impl
instead.
BAT.bat_findmode
— Functionbat_findmode(
target::BAT.AnySampleable,
[algorithm::BAT.AbstractModeEstimator],
[context::BATContext]
)::DensitySampleVector
Estimate the global mode of target
.
Returns a NamedTuple of the shape
(result = X::DensitySampleVector, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_findmode
, add methods to bat_findmode_impl
instead.
BAT.bat_initval
— Functionbat_initval(
target::BAT.MeasureLike,
[algorithm::BAT.InitvalAlgorithm],
[context::BATContext]
)::V
bat_initval(
target::BAT.MeasureLike,
n::Integer,
[algorithm::BAT.InitvalAlgorithm],
[context::BATContext]
)::AbstractVector{<:V}
Generate one or n
random initial/starting value(s) suitable for target
.
Assuming the variates of target
are of type T
, returns a NamedTuple of the shape
(result = X::AbstractVector{T}, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_initval
, add methods like
bat_initval_impl(target::MeasureLike, algorithm::InitvalAlgorithm, context::BATContext)
bat_initval_impl(target::MeasureLike, n::Integer, algorithm::InitvalAlgorithm, context::BATContext)
to bat_initval_impl
instead.
BAT.bat_integrate
— Functionbat_integrate(
target::AnySampleable,
[algorithm::IntegrationAlgorithm],
[context::BATContext]
)::DensitySampleVector
Calculate the integral (evidence) of target
.
Returns a NamedTuple of the shape
(result = X::Measurements.Measurement, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_integrate
, add methods to bat_integrate_impl
instead.
BAT.bat_read
— Functionbat_read(
filename::AbstractString,
[key,]
[algorithm::BATIOAlgorithm]
)
Read data (optionally selected by key
) from filename
using algorithm
.
Example:
smpls = bat_read("samples.hdf5", smpls).result
Returns (result = content, ...)
Result properties not listed here are specific to the output algorithm
and are not part of the stable public API.
See bat_write
.
Currently supported file formats are:
- HDF5 with file extension ".h5" or ".hdf5"
HDF5 I/O functionality is only available when the HDF5 package is loaded (e.g. via import HDF5
).
Do not add add algorithms to bat_read
, add algorithms to bat_read_impl
instead.
BAT.bat_report
— Functionbat_report(obj...)::Markdown.MD
Generate a report on the given object(s).
Specialize bat_report!
instead of bat_report
.
BAT.bat_sample
— Functionbat_sample(
target::BAT.AnySampleable,
[algorithm::BAT.AbstractSamplingAlgorithm],
[context::BATContext]
)::DensitySampleVector
Draw samples from target
using algorithm
.
Depending on sampling algorithm, the samples may be independent or correlated (e.g. when using MCMC).
Returns a NamedTuple of the shape
(result = X::DensitySampleVector, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
Do not add add methods to bat_sample
, add methods to bat_sample_impl
instead.
BAT.bat_write
— Functionbat_write(
filename::AbstractString,
content,
[algorithm::BATIOAlgorithm]
)
Write content
to file filename
using algorithm
.
Example:
smpls = bat_sample(posterior, ...).result
bat_write("samples.hdf5", smpls)
Returns (result = filename, ...)
Result properties not listed here are specific to the output algorithm
and are not part of the stable public API.
See bat_read
.
Currently supported file formats are:
- HDF5 with file extension ".h5" or ".hdf5"
HDF5 I/O functionality is only available when the HDF5 package is loaded (e.g. via import HDF5
).
Do not add add algorithms to bat_write
, add algorithms to bat_write_impl
instead.
BAT.bat_transform
— Functionbat_transform(
how::AbstractTransformTarget,
object,
[algorithm::TransformAlgorithm]
)
bat_transform(
f,
object,
[algorithm::TransformAlgorithm]
)
Transform object
to another variate space defined/implied by target
, res. using the transformation function f
.
Returns a NamedTuple of the shape
(result = newdensity, trafo = vartrafo::Function, ...)
Result properties not listed here are algorithm-specific and are not part of the stable public API.
As a convenience,
flat_smpls, f_flatten = bat_transform(Vector, measure)
flat_smpls, f_flatten = bat_transform(Vector, samples)
can be used to flatten a the variate type of a measure (res. samples of a measure) to something like Vector{<:Real}
.
BAT.get_batcontext
— Functionget_batcontext()::BATContext
get_batcontext(obj)::BATContext
Gets resp. sets the default computational context for BAT.
Will create and set a new default context if none exists.
Note: get_batcontext()
does not have a stable return type. Code that needs type stability should pass a context to algorithms explicitly. BAT algorithms that call other algorithms must forward their context automatically, so context is always type stable within nested BAT algorithms.
See BATContext
and set_batcontext
.
BAT.set_batcontext
— Functionset_batcontext(new_context::BATContext)
set_batcontext(;
precision = ...,
rng = ...,
ad = ...
)
Sets the default computational context for BAT.
See BATContext
and get_batcontext
.
BAT.distbind
— Functiondistbind(f_k, dist, ::typeof(merge))
Performs a generalized monadic bind, in the functional programming sense, with a transition kernel f_k
, a distribution dist
, using merge
to control the type of "flattening".
BAT.distprod
— Functiondistprod(;a = some_dist, b = some_other_dist, ...)
distprod(();a = some_dist, b = some_other_dist, ...))
distprod([dist1, dist2, dist2, ...])
Generate a product of distributions, returning either a distribution that has NamedTuples as variates, or arrays as variates.
BAT.lbqintegral
— Functionlbqintegral(integrand, measure)
lbqintegral(likelihood, prior)
Returns an object that represents the Lebesgue integral over a function in respect to s reference measure. It is also the non-normalized posterior measure that results from integrating the likelihood of a given observation in respect to a prior measure.
BAT.AbstractMCMCWeightingScheme
— Typeabstract type AbstractMCMCWeightingScheme{T<:Real}
Abstract class for weighting schemes for MCMC samples.
Weight values will have type T
.
BAT.AbstractPosteriorMeasure
— Typeabstract type AbstractPosteriorMeasure <: BATMeasure end
Abstract type for posterior probability densities.
BAT.AbstractTransformTarget
— Typeabstract type AbstractTransformTarget
Abstract type for probability density transformation targets.
BAT.AdaptiveMHTuning
— Typestruct AdaptiveMHTuning <: MHProposalDistTuning
Adaptive MCMC tuning strategy for Metropolis-Hastings samplers.
Adapts the proposal function based on the acceptance ratio and covariance of the previous samples.
Constructors:
AdaptiveMHTuning(; fields...)
Fields:
λ::Float64
: Controls the weight given to new covariance information in adapting the proposal distribution. Default: 0.5α::IntervalSets.ClosedInterval{Float64}
: Metropolis-Hastings acceptance ratio target, tuning will try to adapt the proposal distribution to bring the acceptance ratio inside this interval. Default: ClosedInterval(0.15, 0.35)β::Float64
: Controls how much the spread of the proposal distribution is widened/narrowed depending on the current MH acceptance ratio. Default: 1.5c::IntervalSets.ClosedInterval{Float64}
: Interval for allowed scale/spread of the proposal distribution. Default: ClosedInterval(0.0001, 100.0)r::Real
: Reweighting factor. Take accumulated sample statistics of previous tuning cycles into account with a relative weight ofr
. Set to0
to completely reset sample statistics between each tuning cycle. Default: 0.5
BAT.AssumeConvergence
— Typestruct AssumeConvergence <: ConvergenceTest
No-op convergence algorithm for bat_convergence
, will always declare convergence.
Constructors:
AssumeConvergence(converged::Bool = true)
Fields:
converged::Bool
: Default: true
BAT.AutocorLenAlgorithm
— Typeabstract type AutocorLenAlgorithm
Abstract type for integrated autocorrelation length estimation algorithms.
BAT.BATContext
— Typestruct BATContext{T}
Experimental feature, not yet part of stable public API.
Set the default computational context for BAT.
Constructors:
BATContext{T}(rng::AbstractRNG, cunit::AbstractComputeUnit, ADSelector::AD)
BATContext(;
precision::Type{<:AbstractFloat} = ...,
rng::AbstractRNG = ...,
cunit::HeterogeneousComputing.AbstractComputeUnit = ...,
ad::AutoDiffOperators.ADSelector = ...,
)
See get_batcontext
and set_batcontext
.
BAT.BATHDF5IO
— Typestruct BATHDF5IO <: BATIOAlgorithm
Selects the BAT HDF5 format as the output format.
Constructors:
BATHDF5IO()
BAT.BATIOAlgorithm
— Typeabstract type BATIOAlgorithm
Abstract type for density transformation algorithms.
BAT.BrooksGelmanConvergence
— Typestruct BrooksGelmanConvergence <: ConvergenceTest
Brooks-Gelman maximum R^2 convergence test.
Constructors:
BrooksGelmanConvergence(; fields...)
Fields:
threshold::Float64
: Default: 1.1corrected::Bool
: Default: false
BAT.CuhreIntegration
— Typestruct CuhreIntegration <: IntegrationAlgorithm
CuhreIntegration integration algorithm.
Constructors:
CuhreIntegration(; fields...)
Fields:
trafo::AbstractTransformTarget
: Default: PriorToUniform()rtol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:RTOL))atol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:ATOL))minevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MINEVALS))maxevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXEVALS))key::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:KEY))nthreads::Int64
: Default: Base.Threads.nthreads()strict::Bool
: Default: true
This functionality is only available when the Cuba package is loaded (e.g. via import CUBA
).
BAT.DensitySample
— Typestruct DensitySample
A weighted sample drawn according to an statistical density, e.g. a BAT.MeasureLike
.
Constructors:
DensitySampleVector(v::Any, logd::Real, weight::Real, info::Any, aux::Any)
Fields:
v::Any
: variate valuelogd::Real
: log(density) value atv
weight::Real
: Weight of the sampleinfo::Any
: Additional info on the provenance of the sample. Content depends on the sampling algorithm.aux::Any
: Custom user-defined information attached to the sample.
Use DensitySampleVector
to store vectors of multiple samples with an efficient column-based memory layout.
BAT.DensitySampleVector
— Typestruct DensitySampleVector <: AbstractVector{<:DensitySample}
A vector of DensitySample
elements.
DensitySampleVector
is currently a type alias for StructArrays.StructArray{<:DensitySample,...}
, though this is subject to change without deprecation.
Constructors:
DensitySampleVector(
(
v::AbstractVector{<:AbstractVector{<:Real}},
logd::AbstractVector{<:Real},
weight::AbstractVector{<:Real},
info::AbstractVector{<:Any},
aux::AbstractVector{<:Any}
)
)
DensitySampleVector(
v::AbstractVector,
logval::AbstractVector{<:Real};
weight::Union{AbstractVector{<:Real}, Symbol} = fill(1, length(eachindex(v))),
info::AbstractVector = fill(nothing, length(eachindex(v))),
aux::AbstractVector = fill(nothing, length(eachindex(v)))
)
With weight = :multiplicity
repeated samples will be replaced by a single sample, with a weight equal to the number of repetitions.
BAT.DivonneIntegration
— Typestruct DivonneIntegration <: IntegrationAlgorithm
DivonneIntegration integration algorithm.
Constructors:
DivonneIntegration(; fields...)
Fields:
trafo::AbstractTransformTarget
: Default: PriorToUniform()rtol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:RTOL))atol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:ATOL))minevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MINEVALS))maxevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXEVALS))key1::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:KEY1))key2::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:KEY2))key3::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:KEY3))maxpass::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXPASS))border::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:BORDER))maxchisq::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXCHISQ))mindeviation::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MINDEVIATION))ngiven::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NGIVEN))ldxgiven::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:LDXGIVEN))nextra::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NEXTRA))nthreads::Int64
: Default: Base.Threads.nthreads()strict::Bool
: Default: true
This functionality is only available when the Cuba package is loaded (e.g. via import CUBA
).
BAT.DoNotTransform
— Typestruct DoNotTransform <: AbstractTransformTarget
The identity density transformation target, specifies that densities should not be transformed.
Constructors:
DoNotTransform()
BAT.EffSampleSizeAlgorithm
— Typeabstract type EffSampleSizeAlgorithm
Abstract type for integrated autocorrelation length estimation algorithms.
BAT.EffSampleSizeFromAC
— Typestruct EffSampleSizeFromAC <: EffSampleSizeAlgorithm
Effective sample size estimation based on the integrated autocorrelation length of the samples.
Constructors:
EffSampleSizeFromAC(; fields...)
Fields:
acalg::AutocorLenAlgorithm
: Default: GeyerAutocorLen()
BAT.EvaluatedMeasure
— Typestruct EvaluatedMeasure <: BATMeasure
Combined a measure with samples, and other information on it.
Constructors:
EvaluatedMeasure(
measure;
samples = ..., approx = ..., mass = ..., mode = ...,
_generator = ...
)
Field _generator
does not form part of the stable public API and is subject to change without deprecation.
BAT.ExplicitInit
— Typestruct ExplicitInit <: InitvalAlgorithm
Uses initial values from a given vector of one or more values/variates. The values are used in the order they appear in the vector, not randomly.
Constructors:
ExplicitInit(; fields...)
Fields:
xs::AbstractVector
BAT.GelmanRubinConvergence
— Typestruct GelmanRubinConvergence <: ConvergenceTest
Gelman-Rubin maximum R^2 convergence test.
Constructors:
GelmanRubinConvergence(; fields...)
Fields:
threshold::Float64
: Default: 1.1
BAT.GeyerAutocorLen
— Typestruct GeyerAutocorLen <: AutocorLenAlgorithm
Integrated autocorrelation length estimation based on Geyer’s initial monotone sequence criterion
See C. J. Geyer, "Practical Markov Chain Monte Carlo" (1992) and C. J. Geyer, "Introduction to Markov Chain Monte Carlo" (2011).
Constructors:
GeyerAutocorLen()
The same algorithm is used by STAN (v2.21) and MCMCChains.jl (v3.0, function ess_rhat
).
BAT.HamiltonianMC
— Typestruct HamiltonianMC <: MCMCAlgorithm
The Hamiltonian Monte Carlo (HMC) sampling algorithm.
Uses the HMC implementation provided by the package AdvancedHMC.
HMC uses gradients of the target measure's density, so your BATContext
needs to include an ADSelector
to specify which automatic differentiation backend should be used.
- Note: The fields of
HamiltonianMC
are still subject to change, and not
yet part of stable public BAT API!*
Constructors:
HamiltonianMC(; fields...)
Fields:
metric::BAT.HMCMetric
: Default: DiagEuclideanMetric()integrator::Any
: Default: extdefault(pkgext(Val(:AdvancedHMC)), Val(:DEFAULTINTEGRATOR))termination::Any
: Default: extdefault(pkgext(Val(:AdvancedHMC)), Val(:DEFAULTTERMINATION_CRITERION))tuning::BAT.HMCTuningAlgorithm
: Default: StanHMCTuning()
HamiltonianMC
is only available if the AdvancedHMC package is loaded (e.g. via import AdvancedHMC
).
BAT.IdentityTransformAlgorithm
— Typestruct IdentityTransformAlgorithm <: TransformAlgorithm
A no-op density transform algorithm that leaves any density unchanged.
Constructors:
IdentityTransformAlgorithm()
BAT.IIDSampling
— Typestruct IIDSampling <: AbstractSamplingAlgorithm
Sample via Random.rand
.
Constructors:
IIDSampling(; fields...)
Fields:
nsamples::Int64
: Default: 10 ^ 5
BAT.InitFromIID
— Typestruct InitFromIID <: InitvalAlgorithm
Generates initial values for sampling, optimization, etc. by random resampling from a given set of samples.
Constructors:
InitFromIID()
BAT.InitFromSamples
— Typestruct InitFromSamples <: InitvalAlgorithm
Generates initial values for sampling, optimization, etc. by direct sampling from a given i.i.d. sampleable source.
Constructors:
InitFromSamples()
BAT.InitFromTarget
— Typestruct InitFromTarget <: InitvalAlgorithm
Generates initial values for sampling, optimization, etc. by direct i.i.d. sampling a suitable component of that target density (e.g. it's prior) that supports it.
If the target is supports direct i.i.d. sampling, e.g. because it is a distribution, initial values are sampled directly from the target.
If the target is a posterior density, initial values are sampled from the prior (or the prior's prior if the prior is a posterior itself, etc.).
If the target is a sampled density, initial values are (re-)sampled from the available samples.
Constructors:
InitFromTarget()
BAT.InitvalAlgorithm
— Typeabstract type BAT.InitvalAlgorithm
Abstract type for BAT initial/starting value generation algorithms.
Many algorithms in BAT, like MCMC and optimization, need initial/starting values.
BAT.IntegrationAlgorithm
— Typeabstract type IntegrationAlgorithm
Abstract type for integration algorithms.
BAT.KishESS
— Typestruct KishESS <: EffSampleSizeAlgorithm
Kish's effective sample size estimator, uses only the sample weights.
Constructors:
KishESS()
BAT.MaxDensitySearch
— TypeMaxDensitySearch <: AbstractModeEstimator
Constructors:
MaxDensitySearch()
Estimate the mode as the variate with the highest posterior density value within a given set of samples.
BAT.MCMCAlgorithm
— Typeabstract type MCMCAlgorithm
Abstract type for Markov chain Monte Carlo algorithms.
To implement a new MCMC algorithm, subtypes of both MCMCAlgorithm
and MCMCIterator
are required.
The details of the MCMCIterator
and MCMCAlgorithm
API required to implement a new MCMC algorithm currently do not (yet) form part of the stable API and are subject to change without deprecation.
BAT.MCMCBurninAlgorithm
— Typeabstract type MCMCBurninAlgorithm
Abstract type for MCMC burn-in algorithms.
BAT.MCMCChainPoolInit
— Typestruct MCMCChainPoolInit <: MCMCInitAlgorithm
MCMC chain pool initialization strategy.
Constructors:
MCMCChainPoolInit(; fields...)
Fields:
init_tries_per_chain::IntervalSets.ClosedInterval{Int64}
: Default: ClosedInterval(8, 128)nsteps_init::Int64
: Default: 1000initval_alg::InitvalAlgorithm
: Default: InitFromTarget()
BAT.MCMCInitAlgorithm
— Typeabstract type MCMCInitAlgorithm
Abstract type for MCMC initialization algorithms.
BAT.MCMCMultiCycleBurnin
— Typestruct MCMCMultiCycleBurnin <: MCMCBurninAlgorithm
A multi-cycle MCMC burn-in algorithm.
Constructors:
MCMCMultiCycleBurnin(; fields...)
Fields:
nsteps_per_cycle::Int64
: Default: 10000max_ncycles::Int64
: Default: 30nsteps_final::Int64
: Default: div(nstepspercycle, 10)
BAT.MCMCNoOpTuning
— TypeMCMCNoOpTuning <: MCMCTuningAlgorithm
No-op tuning, marks MCMC chains as tuned without performing any other changes on them. Useful if chains are pre-tuned or tuning is an internal part of the MCMC sampler implementation.
BAT.MCMCSampling
— Typestruct MCMCSampling <: AbstractSamplingAlgorithm
Samples a probability density using Markov chain Monte Carlo.
Constructors:
MCMCSampling(; fields...)
Fields:
mcalg::MCMCAlgorithm
: Default: MetropolisHastings()trafo::AbstractTransformTarget
: Default: bat_default(MCMCSampling, Val(:trafo), mcalg)nchains::Int64
: Default: 4nsteps::Int64
: Default: bat_default(MCMCSampling, Val(:nsteps), mcalg, trafo, nchains)init::MCMCInitAlgorithm
: Default: bat_default(MCMCSampling, Val(:init), mcalg, trafo, nchains, nsteps)burnin::MCMCBurninAlgorithm
: Default: bat_default(MCMCSampling, Val(:burnin), mcalg, trafo, nchains, nsteps)convergence::BAT.ConvergenceTest
: Default: BrooksGelmanConvergence()strict::Bool
: Default: truestore_burnin::Bool
: Default: falsenonzero_weights::Bool
: Default: truecallback::Function
: Default: nop_func
BAT.MCMCTuningAlgorithm
— Typeabstract type MCMCTuningAlgorithm
Abstract type for MCMC tuning algorithms.
BAT.MetropolisHastings
— Typestruct MetropolisHastings <: MCMCAlgorithm
Metropolis-Hastings MCMC sampling algorithm.
Constructors:
MetropolisHastings(; fields...)
Fields:
proposal::Distributions.ContinuousDistribution
: Default: TDist(1.0)weighting::AbstractMCMCWeightingScheme
: Default: RepetitionWeighting()tuning::MHProposalDistTuning
: Default: AdaptiveMHTuning()
BAT.MHProposalDistTuning
— Typeabstract type MHProposalDistTuning
Abstract type for Metropolis-Hastings tuning strategies for proposal distributions.
BAT.ModeAsDefined
— Typestruct ModeAsDefined <: AbstractModeEstimator
Get the mode as defined by the density, resp. the underlying distribution (if available), via StatsBase.mode
.
Constructors:
ModeAsDefined()
BAT.OptimAlg
— TypeOptimAlg
Selects an optimization algorithm from the Optim.jl package.
Note that when using first order algorithms like Optim.LBFGS
, your BATContext
needs to include an ADSelector
that specifies which automatic differentiation backend should be used.
Constructors:
OptimAlg(; fields...)
optimalg
must be an Optim.AbstractOptimizer
.
Fields:
optalg::Any
: Default: extdefault(pkgext(Val(:Optim)), Val(:DEFAULTOPTALG))trafo::AbstractTransformTarget
: Default: PriorToGaussian()init::InitvalAlgorithm
: Default: InitFromTarget()maxiters::Int64
: Default: 1000maxtime::Float64
: Default: NaNabstol::Float64
: Default: NaNreltol::Float64
: Default: 0.0kwargs::NamedTuple
: Default: (;)
This algorithm is only available if the Optim package is loaded (e.g. via import Optim
.
BAT.OptimizationAlg
— TypeOptimizationAlg
Selects an optimization algorithm from the Optimization.jl package. Note that when using first order algorithms like OptimizationOptimJL.LBFGS
, your BATContext
needs to include an ADSelector
that specifies which automatic differentiation backend should be used. Constructors:
OptimizationAlg(; fields...)
optalg
must be an Optimization.AbstractOptimizer
. The field kwargs
can be used to pass additional keywords to the optimizers See the Optimization.jl documentation for the available keyword arguments. Fields:
optalg::Any
: Default: extdefault(pkgext(Val(:Optimization)), Val(:DEFAULTOPTALG))trafo::AbstractTransformTarget
: Default: PriorToGaussian()init::InitvalAlgorithm
: Default: InitFromTarget()maxiters::Int64
: Default: 1000maxtime::Float64
: Default: NaNabstol::Float64
: Default: NaNreltol::Float64
: Default: 0.0kwargs::NamedTuple
: Default: (;)
This algorithm is only available if the Optimization
package or any of its submodules, like OptimizationOptimJL
, is loaded (e.g. via import Optimization
).
BAT.OrderedResampling
— Typestruct OrderedResampling <: AbstractSamplingAlgorithm
Efficiently resamples from a given series of samples, keeping the order of samples.
Can be used to efficiently convert weighted samples into samples with unity weights.
Constructors:
OrderedResampling(; fields...)
Fields:
nsamples::Int64
: Default: 10 ^ 5
BAT.PosteriorMeasure
— Typestruct PosteriorMeasure{Li,Pr<:AbstractMeasure} <: AbstractPosteriorMeasure
A representation of a PosteriorMeasure, based a likelihood and prior. Likelihood and prior be accessed via
getlikelihood(posterior::PosteriorMeasure)::Li
getprior(posterior::PosteriorMeasure)::Pr
Constructors:
PosteriorMeasure(likelihood, prior)
PosteriorMeasure{T<:Real}(likelihood, prior)
Fields:
likelihood::Any
prior::MeasureBase.AbstractMeasure
BAT.PriorSubstitution
— Typestruct PriorSubstitution <: TransformAlgorithm
Substitute the prior by a given distribution and transform the likelihood accordingly. The log(abs(jacobian))
of the transformation does not need to be auto-differentiable even for operations that use the gradient of the posterior.
Constructors:
PriorSubstitution()
BAT.PriorToGaussian
— Typestruct PriorToGaussian <: AbstractTransformTarget
Specifies that posterior densities should be transformed in a way that makes their pior equivalent to a standard multivariate normal distribution with an identity covariance matrix.
Constructors:
PriorToGaussian()
BAT.PriorToUniform
— Typestruct PriorToUniform <: AbstractTransformTarget
Specifies that posterior densities should be transformed in a way that makes their pior equivalent to a uniform distribution over the unit hypercube.
Constructors:
PriorToUniform()
BAT.RandResampling
— Typestruct RandResampling <: AbstractSamplingAlgorithm
Resamples from a given set of samples.
Constructors:
RandResampling(; fields...)
Fields:
nsamples::Int64
: Default: 10 ^ 5
BAT.RepetitionWeighting
— Typestruct RepetitionWeighting{T<:AbstractFloat} <: AbstractMCMCWeightingScheme{T}
Sample weighting scheme suitable for sampling algorithms which may repeated samples multiple times in direct succession (e.g. MetropolisHastings
). The repeated sample is stored only once, with a weight equal to the number of times it has been repeated (e.g. because a Markov chain has not moved during a sampling step).
Constructors:
RepetitionWeighting()
BAT.SampleMedianEstimator
— Typestruct SampleMedianEstimator <: AbstractMedianEstimator
Get median values from samples using standard Julia statistics functions.
Constructors:
SampleMedianEstimator()
BAT.SokalAutocorLen
— Typestruct SokalAutocorLen <: AutocorLenAlgorithm
Integrated autocorrelation length estimation based on the automated windowing procedure descibed in A. D. Sokal, "Monte Carlo Methods in Statistical Mechanics" (1996)
Same procedure is used by the emcee Python package (v3.0).
Constructors:
SokalAutocorLen(; fields...)
Fields:
c::Int64
: Step size for window search Default: 5
BAT.SuaveIntegration
— Typestruct SuaveIntegration <: IntegrationAlgorithm
SuaveIntegration integration algorithm.
Constructors:
SuaveIntegration(; fields...)
Fields:
trafo::AbstractTransformTarget
: Default: PriorToUniform()rtol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:RTOL))atol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:ATOL))minevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MINEVALS))maxevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXEVALS))nnew::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NNEW))nmin::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NMIN))flatness::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:FLATNESS))nthreads::Int64
: Default: Base.Threads.nthreads()strict::Bool
: Default: true
This functionality is only available when the Cuba package is loaded (e.g. via import CUBA
).
BAT.TransformAlgorithm
— Typeabstract type TransformAlgorithm
Abstract type for density transformation algorithms.
BAT.VEGASIntegration
— Typestruct VEGASIntegration <: IntegrationAlgorithm
VEGASIntegration integration algorithm.
Constructors:
VEGASIntegration(; fields...)
Fields:
trafo::AbstractTransformTarget
: Default: PriorToUniform()rtol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:RTOL))atol::Float64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:ATOL))minevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MINEVALS))maxevals::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:MAXEVALS))nstart::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NSTART))nincrease::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NINCREASE))nbatch::Int64
: Default: ext_default(pkgext(Val(:Cuba)), Val(:NBATCH))nthreads::Int64
: Default: Base.Threads.nthreads()strict::Bool
: Default: true
This functionality is only available when the Cuba package is loaded (e.g. via import CUBA
).
BAT.BinningAlgorithm
— Typeabstract type BinningAlgorithm
Abstract type for binning algorithms.
BAT.FixedNBins
— TypeFixedNBins(nbins::Int)
Selects a fixed number of bins.
Constructor: FixedNBins(; fields...)
Fields:
nbins::Int64
: Default: 200
BAT.FreedmanDiaconisBinning
— Typestruct FreedmanDiaconisBinning <: BinningAlgorithm
Selects automatic binning based on the Freedman–Diaconis rule.
Constructor: FreedmanDiaconisBinning()
BAT.RiceBinning
— Typestruct RiceBinning <: BinningAlgorithm
Selects automatic binning based on the Rice rule.
Constructor: RiceBinning()
BAT.ScottBinning
— Typestruct ScottBinning <: BinningAlgorithm
Selects automatic binning based on Scott's normal reference rule.
Constructor: ScottBinning()
BAT.SquareRootBinning
— Typestruct SquareRootBinning <: BinningAlgorithm
Selects automatic binning based on the Square-root choice.
Constructor: SquareRootBinning()
BAT.SturgesBinning
— Typestruct SturgesBinning <: BinningAlgorithm
Selects automatic binning based on Sturges' formula.
Constructor: SturgesBinning()
BAT.ToRealVector
— Typestruct ToRealVector <: AbstractTransformTarget
Specifies that the input should be transformed into a measure over the space of real-valued flat vectors.
Constructors:
ToRealVector()
BAT.AbstractMedianEstimator
— Typeabstract type BAT.AbstractMedianEstimator
Abstract type for BAT optimization algorithms.
A typical application for optimization in BAT is mode estimation (see bat_findmode
),
BAT.AbstractModeEstimator
— Typeabstract type BAT.AbstractModeEstimator
Abstract type for BAT optimization algorithms.
A typical application for optimization in BAT is mode estimation (see bat_findmode
),
BAT.AbstractSamplingAlgorithm
— Typeabstract type BAT.AbstractSamplingAlgorithm
Abstract type for BAT sampling algorithms. See bat_sample
.
BAT.ConvergenceTest
— Typeabstract type ConvergenceTest
Abstract type for integrated autocorrelation length estimation algorithms.