![]() |
Stan
2.10.0
probability, sampling & optimization
|
Automatic Differentiation Variational Inference. More...
#include <advi.hpp>
Public Member Functions | |
advi (Model &m, Eigen::VectorXd &cont_params, BaseRNG &rng, int n_monte_carlo_grad, int n_monte_carlo_elbo, int eval_elbo, int n_posterior_samples) | |
Constructor. More... | |
double | calc_ELBO (const Q &variational, interface_callbacks::writer::base_writer &message_writer) const |
Calculates the Evidence Lower BOund (ELBO) by sampling from the variational distribution and then evaluating the log joint, adjusted by the entropy term of the variational distribution. More... | |
void | calc_ELBO_grad (const Q &variational, Q &elbo_grad, interface_callbacks::writer::base_writer &message_writer) const |
Calculates the "black box" gradient of the ELBO. More... | |
double | adapt_eta (Q &variational, int adapt_iterations, interface_callbacks::writer::base_writer &message_writer) const |
Heuristic grid search to adapt eta to the scale of the problem. More... | |
void | stochastic_gradient_ascent (Q &variational, double eta, double tol_rel_obj, int max_iterations, interface_callbacks::writer::base_writer &message_writer, interface_callbacks::writer::base_writer &diagnostic_writer) const |
Runs stochastic gradient ascent with an adaptive stepsize sequence. More... | |
int | run (double eta, bool adapt_engaged, int adapt_iterations, double tol_rel_obj, int max_iterations, interface_callbacks::writer::base_writer &message_writer, interface_callbacks::writer::base_writer ¶meter_writer, interface_callbacks::writer::base_writer &diagnostic_writer) const |
Runs ADVI and writes to output. More... | |
double | circ_buff_median (const boost::circular_buffer< double > &cb) const |
Compute the median of a circular buffer. More... | |
double | rel_difference (double prev, double curr) const |
Compute the relative difference between two double values. More... | |
Protected Attributes | |
Model & | model_ |
Eigen::VectorXd & | cont_params_ |
BaseRNG & | rng_ |
int | n_monte_carlo_grad_ |
int | n_monte_carlo_elbo_ |
int | eval_elbo_ |
int | n_posterior_samples_ |
Automatic Differentiation Variational Inference.
Implements "black box" variational inference using stochastic gradient ascent to maximize the Evidence Lower Bound for a given model and variational family.
Model | class of model |
Q | class of variational distribution |
BaseRNG | class of random number generator |
|
inline |
Constructor.
m | stan model |
cont_params | initialization of continuous parameters |
rng | random number generator |
n_monte_carlo_grad | number of samples for gradient computation |
n_monte_carlo_elbo | number of samples for ELBO computation |
eval_elbo | evaluate ELBO at every "eval_elbo" iters |
n_posterior_samples | number of samples to draw from posterior |
std::runtime_error | if n_monte_carlo_grad is not positive |
std::runtime_error | if n_monte_carlo_elbo is not positive |
std::runtime_error | if eval_elbo is not positive |
std::runtime_error | if n_posterior_samples is not positive |
|
inline |
Heuristic grid search to adapt eta to the scale of the problem.
[in] | variational | initial variational distribution. |
adapt_iterations | number of iterations to spend doing stochastic gradient ascent at each proposed eta value. | |
message_writer | writer for messages |
std::domain_error | If either (a) the initial ELBO cannot be computed at the initial variational distribution, (b) all step-size proposals in eta_sequence fail. |
|
inline |
Calculates the Evidence Lower BOund (ELBO) by sampling from the variational distribution and then evaluating the log joint, adjusted by the entropy term of the variational distribution.
[in] | variational | variational approximation at which to evaluate the ELBO. |
message_writer | writer for messages |
std::domain_error | If, after n_monte_carlo_elbo_ number of draws from the variational distribution all give non-finite log joint evaluations. This means that the model is severly ill conditioned or that the variational distribution has somehow collapsed. |
|
inline |
|
inline |
|
inline |
|
inline |
Runs ADVI and writes to output.
eta | eta parameter of stepsize sequence |
adapt_engaged | boolean flag for eta adaptation |
adapt_iterations | number of iterations for eta adaptation |
tol_rel_obj | relative tolerance parameter for convergence |
max_iterations | max number of iterations to run algorithm |
message_writer | writer for messages |
parameter_writer | writer for parameters (typically to file) |
diagnostic_writer | writer for diagnostic information |
|
inline |
Runs stochastic gradient ascent with an adaptive stepsize sequence.
[in,out] | variational | initia variational distribution |
eta | stepsize scaling parameter | |
tol_rel_obj | relative tolerance parameter for convergence | |
max_iterations | max number of iterations to run algorithm | |
message_writer | writer for mesasges | |
diagnostic_writer | writer for diagnostic information |
std::domain_error | If the ELBO or its gradient is ever non-finite, at any iteration |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |