Bayesian Updating
Methods for Bayesian updating.
Index
UncertaintyQuantification.AbstractBayesianMethodUncertaintyQuantification.AbstractBayesianPointEstimateUncertaintyQuantification.LaplaceEstimateBayesianUncertaintyQuantification.MaximumAPosterioriBayesianUncertaintyQuantification.MaximumLikelihoodBayesianUncertaintyQuantification.SingleComponentMetropolisHastingsUncertaintyQuantification.TransitionalMarkovChainMonteCarloUncertaintyQuantification.TransportMapBayesianUncertaintyQuantification.bayesianupdating
Types
UncertaintyQuantification.AbstractBayesianMethod Type
AbstractBayesianMethodSubtypes are used to dispatch to the different MCMC methods in bayesianupdating.
Subtypes are:
sourceUncertaintyQuantification.SingleComponentMetropolisHastings Type
SingleComponentMetropolisHastings(proposal, x0, n, burnin, islog)Passed to bayesianupdating to run the single-component Metropolis-Hastings algorithm starting from x0 with univariate proposal distibution proposal (or vector of proposal distributions per dimension). Will generate n samples after performing burnin steps of the Markov chain and discarding the samples. The flag islog specifies whether the prior and likelihood functions passed to the bayesianupdating method are already given as logarithms.
Alternative constructor
SingleComponentMetropolisHastings(proposal, x0, n, burnin) # `islog` = trueUncertaintyQuantification.TransitionalMarkovChainMonteCarlo Type
TransitionalMarkovChainMonteCarlo(prior, n, burnin, β, islog)
Passed to [`bayesianupdating`](@ref) to run thetransitional Markov chain Monte Carlo algorithm with [`RandomVariable'](@ref) vector `prior`. At each transitional level, one sample will be generated from `n` independent Markov chains after `burnin` steps have been discarded. The flag `islog` specifies whether the prior and likelihood functions passed to the [`bayesianupdating`](@ref) method are already given as logarithms.Alternative constructors
TransitionalMarkovChainMonteCarlo(prior, n, burnin, β) # `islog` = true
TransitionalMarkovChainMonteCarlo(prior, n, burnin) # `β` = 0.2, `islog` = trueReferences
[36]
sourceUncertaintyQuantification.AbstractBayesianPointEstimate Type
AbstractBayesianPointEstimateSubtypes are used to dispatch to the different point estimation methods in bayesianupdating.
Subtypes are:
sourceUncertaintyQuantification.MaximumAPosterioriBayesian Type
MaximumAPosterioriBayesian(prior, optimmethod, x0; islog, lowerbounds, upperbounds)Passed to bayesianupdating to estimate one or more maxima of the posterior distribution starting from x0. The optimization uses the method specified in optimmethod. Will calculate one estimation per point in x0. The flag islog specifies whether the prior and likelihood functions passed to the bayesianupdating method are already given as logarithms. lowerbounds and upperbounds specify optimization intervals.
Alternative constructors
MaximumAPosterioriBayesian(prior, optimmethod, x0; islog) # `lowerbounds` = [-Inf], # `upperbounds` = [Inf]
MaximumAPosterioriBayesian(prior, optimmethod, x0) # `islog` = trueSee also MaximumLikelihoodBayesian, bayesianupdating, TransitionalMarkovChainMonteCarlo.
UncertaintyQuantification.MaximumLikelihoodBayesian Type
MaximumLikelihoodBayesian(prior, x0; islog, lowerbounds, upperbounds)Passed to bayesianupdating to estimate one or more maxima of the likelihood starting from x0. The optimization uses the method specified in optimmethod. Will calculate one estimation per point in x0. The flag islog specifies whether the prior and likelihood functions passed to the bayesianupdating method are already given as logarithms. lowerbounds and upperbounds specify optimization intervals.
Alternative constructors
MaximumLikelihoodBayesian(prior, x0; islog) # `lowerbounds` = [-Inf], # `upperbounds` = [Inf]
MaximumLikelihoodBayesian(prior, x0) # `islog` = trueNotes
The method uses prior only as information on which parameters are supposed to be optimized. The prior itself does not influence the result of the maximum likelihood estimate and can be given as a dummy distribution. For example, if two parameters a and b are supposed to be optimized, the prior could look like this
prior = RandomVariable.(Uniform(0,1), [:a, :b])See also MaximumAPosterioriBayesian, bayesianupdating, TransitionalMarkovChainMonteCarlo.
UncertaintyQuantification.LaplaceEstimateBayesian Type
LaplaceEstimateBayesian(prior, optimmethod, x0; islog, lowerbounds, upperbounds)Estimates means and covariances of a mixture of Gaussians to approximate the posterior density. Passed to bayesianupdating to estimate one or more maxima of the posterior starting from x0. The optimization uses the method specified in optimmethod. Will calculate one estimation per point in x0, these are then filtered, s.t. multiple estimates of the same point are discarded. The flag islog specifies whether the prior and likelihood functions passed to the bayesianupdating method are already given as logarithms. Also specifies whether the posterior is given as log-function. lowerbounds and upperbounds specify optimization intervals.
Alternative constructors
LaplaceEstimateBayesian(prior, optimmethod, x0; islog) # `lowerbounds` = [-Inf], # `upperbounds` = [Inf]
LaplaceEstimateBayesian(prior, optimmethod, x0) # `islog` = trueNotes
The method makes use of the MaximumAPosterioriBayesian method to estimate the maximum a posteriori (MAP) estimate, and then calculates the Hessian of the posterior at the MAP estimate to construct a Gaussian approximation of the posterior distribution.
See also MaximumAPosterioriBayesian, bayesianupdating, TransitionalMarkovChainMonteCarlo.
UncertaintyQuantification.TransportMapBayesian Type
TransportMapBayesian(prior, transportmap, quadrature, gradient, optimizer, optim_options, islog, transformprior)Passed to bayesianupdating to perform variational inference with a transport map that maps between the posterior and standard normal space, fitted with mapfromdensity. The prior is a vector of prior distributions, transportmap is the transport map to be optimized, and quadrature specifies the quadrature points in the standard normal space. The gradient can be specified as either an AbstractADType or a Function (default: AutoFiniteDiff()). The optimizer specifies the optimization method from Optim.jl (default: LBFGS()), and optim_options allows passing options to the optimizer (default: Optim.Options()). The flag islog specifies whether densities are already log-densities. The flag transformprior specifies whether the prior should be transformed to standard normal space.
Alternative constructor
TransportMapBayesian(prior, transportmap, quadrature) # gradient = AutoMooncake(), optimizer = LBFGS(), optim_options = Optim.Options(), `islog` = true, `transformprior` = true,References
See also MaximumAPosterioriBayesian, MaximumLikelihoodBayesian, TransitionalMarkovChainMonteCarlo.
Methods
UncertaintyQuantification.bayesianupdating Function
bayesianupdating(likelihood, models, pointestimate; prior, filtertolerance, optimizer, optimoptions)Perform bayesian updating using the given likelihood, models and any point estimation method AbstractBayesianPointEstimate.
Notes
Method can be called with an empty Vector of models, i.e.
bayesianupdating(likelihood, [], pointestimate)If prior is not given, the method will construct a prior distribution from the prior specified in AbstractBayesianPointEstimate.prior.
likelihood is a Julia function which must be defined in terms of a DataFrame of samples, and must evaluate the likelihood for each row of the DataFrame
For example, a loglikelihood based on normal distribution using 'Data':
likelihood(df) = [sum(logpdf.(Normal.(df_i.x, 1), Data)) for df_i in eachrow(df)]If a model evaluation is required to evaluate the likelihood, a vector of UQModels must be passed to bayesianupdating. For example if the variable x above is the output of a numerical model.
filtertolerance is a tolerance value to filter out multiple estimates of the same point. If the distance between two points is smaller than filtertolerance, one of them will be discarded. This is useful if the optimization method finds multiple local maxima that are very close to each other. If all points should be kept, filtertolerance can be set to 0.
For a general overview of the function, see bayesianupdating.
bayesianupdating(likelihood, models, lpestimate; prior, filtertolerance, optimizer, optimoptions, adbackend)Perform bayesian updating with Laplace estimation using the given likelihood, models and the MAP estimation MaximumAPosterioriBayesian. Laplace estimation is basically an extension of the MAP estimation, where the Hessian of the posterior is calculated at the MAP estimate and used to construct a Gaussian approximation of the posterior distribution. Returns a JointDistribution built from the estimated mean values and covariances. The Hessian is estimated using a backend defined from DiffereniationInteface.jl and can be changed using ADTypes.
Notes
Method can be called with an empty Vector of models, i.e.
bayesianupdating(likelihood, [], pointestimate)If prior is not given, the method will construct a prior distribution from the prior specified in AbstractBayesianPointEstimate.prior.
likelihood is a Julia function which must be defined in terms of a DataFrame of samples, and must evaluate the likelihood for each row of the DataFrame
For example, a loglikelihood based on normal distribution using 'Data':
likelihood(df) = [sum(logpdf.(Normal.(df_i.x, 1), Data)) for df_i in eachrow(df)]If a model evaluation is required to evaluate the likelihood, a vector of UQModels must be passed to bayesianupdating. For example if the variable x above is the output of a numerical model.
filtertolerance is a tolerance value to filter out multiple estimates of the same point. If the distance between two points is smaller than filtertolerance, one of them will be discarded. This is useful if the optimization method finds multiple local maxima that are very close to each other.
For a general overview of the function, see bayesianupdating.
bayesianupdating(prior, likelihood, models, tm)Perform bayesian updating using the given prior, likelihood, models and tm TransportMapBayesian. Returns a JointDistribution with the optimized TransportMap.
Alternatively, the prior function can be omitted and the prior is constructed from the given vector of random variables in tm.
bayesianupdating(likelihood, models, tm)See also TransportMapBayesian, MaximumAPosterioriBayesian, `TransportMap.
bayesianupdating(prior, likelihood, models, mcmc)Perform bayesian updating using the given prior, likelihood, models and any MCMC sampler AbstractBayesianMethod.
Alternatively the method can be called without models.
bayesianupdating(prior, likelihood, mcmc)When using TransitionalMarkovChainMonteCarlo the prior can automatically be constructed.
bayesinupdating(likelihood, models, tmcmc)
bayesianupdating(likelihood, tmcmc)Notes
likelihood is a Julia function which must be defined in terms of a DataFrame of samples, and must evaluate the likelihood for each row of the DataFrame
For example, a loglikelihood based on normal distribution using 'Data':
likelihood(df) = [sum(logpdf.(Normal.(df_i.x, 1), Data)) for df_i in eachrow(df)]If a model evaluation is required to evaluate the likelihood, a vector of UQModels must be passed to bayesianupdating. For example if the variable x above is the output of a numerical model.