ThermodynamicIntegration
Thermodynamic Integration
Thermodynamic integration is a technique from physics to get an accurate estimate of the log evidence. By creating a schedule going from the prior to the posterior and estimating the log likelihood at each step one gets a stable ad robust estimate of the log evidence.
A simple example
A simple package to compute Thermodynamic Integration for computing the evidence in a Bayesian setting.
You need to provide the logprior
and the loglikelihood
as well as an initial sample:
using Distributions, ThermodynamicIntegration
D = 5
prior = MvNormal(0.5 * ones(D)) # The prior distribution
likelihood = MvNormal(2.0 * ones(D))
logprior(x) = logpdf(prior, x) # The log-prior function
loglikelihood(x) = logpdf(likelihood, x) # The log-likelihood function
alg = ThermInt(n_samples=5000) # We are going to sample 5000 samples at every step
logZ = alg(logprior, loglikelihood, rand(prior)) # Compute the log evidence
# -8.244829688529377
true_logZ = -0.5 * (logdet(cov(prior) + cov(likelihood)) + D * log(2π)) # we compare twith the true value
# -8.211990123364176
You can also simply pass a Turing model :
using Turing
@model function gauss(y)
x ~ prior
y ~ MvNormal(x, cov(likelihood))
end
alg = ThermInt(n_samples=5000)
model = gauss(zeros(D))
turing_logZ = alg(model)
# # -8.211990123364176
Parallel sampling
The algorithm also works on multiple threads by calling :
alg = ThermInt(n_samples=5000)
logZ = alg(logprior, loglikelihood, rand(prior), TIParallelThreads())
Sampling methods
Right now sampling is based on AdvancedHMC.jl
, with the ForwardDiff
AD backend.
To change the backend to Zygote
or ReverseDiff
(recommended for variables with large dimensions you can do:
using Zygote # (or ReverseDiff)
ThermoDynamicIntegration.set_adbackend(:Zygote) # (or :ReverseDiff)
More samplers will be available in the future.
Further options
You can disactivate the progress by calling progress=false
alg()