Public Member Functions | Public Attributes | Protected Attributes | List of all members
stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG > Class Template Reference

#include "stan/mcmc/hmc/nuts/base_nuts.hpp"

Inheritance diagram for stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >:
stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG > stan::mcmc::base_mcmc

Public Member Functions

 base_nuts (const Model &model, BaseRNG &rng)
 
 base_nuts (const Model &model, BaseRNG &rng, Eigen::VectorXd &inv_e_metric)
 
 base_nuts (const Model &model, BaseRNG &rng, Eigen::MatrixXd &inv_e_metric)
 
 ~base_nuts ()
 
void set_metric (const Eigen::MatrixXd &inv_e_metric)
 
void set_metric (const Eigen::VectorXd &inv_e_metric)
 
void set_max_depth (int d)
 
void set_max_delta (double d)
 
int get_max_depth ()
 
double get_max_delta ()
 
sample transition (sample &init_sample, callbacks::logger &logger)
 
void get_sampler_param_names (std::vector< std::string > &names)
 
void get_sampler_params (std::vector< double > &values)
 
virtual bool compute_criterion (Eigen::VectorXd &p_sharp_minus, Eigen::VectorXd &p_sharp_plus, Eigen::VectorXd &rho)
 
bool build_tree (int depth, ps_point &z_propose, Eigen::VectorXd &p_sharp_left, Eigen::VectorXd &p_sharp_right, Eigen::VectorXd &rho, double H0, double sign, int &n_leapfrog, double &log_sum_weight, double &sum_metro_prob, callbacks::logger &logger)
 
void write_sampler_stepsize (callbacks::writer &writer)
 
void write_sampler_metric (callbacks::writer &writer)
 
void write_sampler_state (callbacks::writer &writer)
 
void get_sampler_diagnostic_names (std::vector< std::string > &model_names, std::vector< std::string > &names)
 
void get_sampler_diagnostics (std::vector< double > &values)
 
void seed (const Eigen::VectorXd &q)
 
void init_hamiltonian (callbacks::logger &logger)
 
void init_stepsize (callbacks::logger &logger)
 
Hamiltonian< Model, BaseRNG >::PointType & z ()
 
virtual void set_nominal_stepsize (double e)
 
double get_nominal_stepsize ()
 
double get_current_stepsize ()
 
virtual void set_stepsize_jitter (double j)
 
double get_stepsize_jitter ()
 
void sample_stepsize ()
 

Public Attributes

int depth_
 
int max_depth_
 
double max_deltaH_
 
int n_leapfrog_
 
bool divergent_
 
double energy_
 

Protected Attributes

Hamiltonian< Model, BaseRNG >::PointType z_
 
Integrator< Hamiltonian< Model, BaseRNG > > integrator_
 
Hamiltonian< Model, BaseRNG > hamiltonian_
 
BaseRNG & rand_int_
 
boost::uniform_01< BaseRNG & > rand_uniform_
 
double nom_epsilon_
 
double epsilon_
 
double epsilon_jitter_
 

Detailed Description

template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
class stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >

The No-U-Turn sampler (NUTS) with multinomial sampling

Definition at line 22 of file base_nuts.hpp.

Constructor & Destructor Documentation

template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::base_nuts ( const Model model,
BaseRNG &  rng 
)
inline

Definition at line 24 of file base_nuts.hpp.

25  : base_hmc<Model, Hamiltonian, Integrator, BaseRNG>(model, rng),
26  depth_(0), max_depth_(5), max_deltaH_(1000),
27  n_leapfrog_(0), divergent_(false), energy_(0) {
28  }
const XML_Char XML_Content * model
Definition: expat.h:151
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::base_nuts ( const Model model,
BaseRNG &  rng,
Eigen::VectorXd &  inv_e_metric 
)
inline

specialized constructor for specified diag mass matrix

Definition at line 33 of file base_nuts.hpp.

35  : base_hmc<Model, Hamiltonian, Integrator, BaseRNG>(model, rng,
36  inv_e_metric),
37  depth_(0), max_depth_(5), max_deltaH_(1000),
38  n_leapfrog_(0), divergent_(false), energy_(0) {
39  }
const XML_Char XML_Content * model
Definition: expat.h:151
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::base_nuts ( const Model model,
BaseRNG &  rng,
Eigen::MatrixXd &  inv_e_metric 
)
inline

specialized constructor for specified dense mass matrix

Definition at line 44 of file base_nuts.hpp.

46  : base_hmc<Model, Hamiltonian, Integrator, BaseRNG>(model, rng,
47  inv_e_metric),
48  depth_(0), max_depth_(5), max_deltaH_(1000),
49  n_leapfrog_(0), divergent_(false), energy_(0) {
50  }
const XML_Char XML_Content * model
Definition: expat.h:151
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::~base_nuts ( )
inline

Definition at line 52 of file base_nuts.hpp.

52 {}

Member Function Documentation

template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
bool stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::build_tree ( int  depth,
ps_point z_propose,
Eigen::VectorXd &  p_sharp_left,
Eigen::VectorXd &  p_sharp_right,
Eigen::VectorXd &  rho,
double  H0,
double  sign,
int n_leapfrog,
double &  log_sum_weight,
double &  sum_metro_prob,
callbacks::logger logger 
)
inline

Recursively build a new subtree to completion or until the subtree becomes invalid. Returns validity of the resulting subtree.

Parameters
depthDepth of the desired subtree
z_proposeState proposed from subtree
p_sharp_leftp_sharp from left boundary of returned tree
p_sharp_rightp_sharp from the right boundary of returned tree
rhoSummed momentum across trajectory
H0Hamiltonian of initial state
signDirection in time to built subtree
n_leapfrogSummed number of leapfrog evaluations
log_sum_weightLog of summed weights across trajectory
sum_metro_probSummed Metropolis probabilities across trajectory
loggerLogger for messages

Definition at line 206 of file base_nuts.hpp.

Referenced by stan::mcmc::base_nuts< Model, diag_e_metric, expl_leapfrog, BaseRNG >::build_tree(), TEST(), and stan::mcmc::base_nuts< Model, diag_e_metric, expl_leapfrog, BaseRNG >::transition().

212  {
213  // Base case
214  if (depth == 0) {
215  this->integrator_.evolve(this->z_, this->hamiltonian_,
216  sign * this->epsilon_,
217  logger);
218  ++n_leapfrog;
219 
220  double h = this->hamiltonian_.H(this->z_);
221  if (boost::math::isnan(h))
222  h = std::numeric_limits<double>::infinity();
223 
224  if ((h - H0) > this->max_deltaH_) this->divergent_ = true;
225 
226  log_sum_weight = math::log_sum_exp(log_sum_weight, H0 - h);
227 
228  if (H0 - h > 0)
229  sum_metro_prob += 1;
230  else
231  sum_metro_prob += std::exp(H0 - h);
232 
233  z_propose = this->z_;
234  rho += this->z_.p;
235 
236  p_sharp_left = this->hamiltonian_.dtau_dp(this->z_);
237  p_sharp_right = p_sharp_left;
238 
239  return !this->divergent_;
240  }
241  // General recursion
242  Eigen::VectorXd p_sharp_dummy(this->z_.p.size());
243 
244  // Build the left subtree
245  double log_sum_weight_left = -std::numeric_limits<double>::infinity();
246  Eigen::VectorXd rho_left = Eigen::VectorXd::Zero(rho.size());
247 
248  bool valid_left
249  = build_tree(depth - 1, z_propose,
250  p_sharp_left, p_sharp_dummy, rho_left,
251  H0, sign, n_leapfrog,
252  log_sum_weight_left, sum_metro_prob,
253  logger);
254 
255  if (!valid_left) return false;
256 
257  // Build the right subtree
258  ps_point z_propose_right(this->z_);
259 
260  double log_sum_weight_right = -std::numeric_limits<double>::infinity();
261  Eigen::VectorXd rho_right = Eigen::VectorXd::Zero(rho.size());
262 
263  bool valid_right
264  = build_tree(depth - 1, z_propose_right,
265  p_sharp_dummy, p_sharp_right, rho_right,
266  H0, sign, n_leapfrog,
267  log_sum_weight_right, sum_metro_prob,
268  logger);
269 
270  if (!valid_right) return false;
271 
272  // Multinomial sample from right subtree
273  double log_sum_weight_subtree
274  = math::log_sum_exp(log_sum_weight_left, log_sum_weight_right);
275  log_sum_weight
276  = math::log_sum_exp(log_sum_weight, log_sum_weight_subtree);
277 
278  if (log_sum_weight_right > log_sum_weight_subtree) {
279  z_propose = z_propose_right;
280  } else {
281  double accept_prob
282  = std::exp(log_sum_weight_right - log_sum_weight_subtree);
283  if (this->rand_uniform_() < accept_prob)
284  z_propose = z_propose_right;
285  }
286 
287  Eigen::VectorXd rho_subtree = rho_left + rho_right;
288  rho += rho_subtree;
289 
290  return compute_criterion(p_sharp_left, p_sharp_right, rho_subtree);
291  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
fvar< T > log_sum_exp(const std::vector< fvar< T > > &v)
Definition: log_sum_exp.hpp:12
bool isnan(const stan::math::var &v)
Definition: boost_isnan.hpp:20
bool build_tree(int depth, ps_point &z_propose, Eigen::VectorXd &p_sharp_left, Eigen::VectorXd &p_sharp_right, Eigen::VectorXd &rho, double H0, double sign, int &n_leapfrog, double &log_sum_weight, double &sum_metro_prob, callbacks::logger &logger)
Definition: base_nuts.hpp:206
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
virtual bool compute_criterion(Eigen::VectorXd &p_sharp_minus, Eigen::VectorXd &p_sharp_plus, Eigen::VectorXd &rho)
Definition: base_nuts.hpp:182
boost::uniform_01< BaseRNG & > rand_uniform_
Definition: base_hmc.hpp:189
void Zero()
Hamiltonian< Model, BaseRNG > hamiltonian_
Definition: base_hmc.hpp:184
Integrator< Hamiltonian< Model, BaseRNG > > integrator_
Definition: base_hmc.hpp:183
def sign(x)
Definition: canMan.py:197
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
virtual bool stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::compute_criterion ( Eigen::VectorXd &  p_sharp_minus,
Eigen::VectorXd &  p_sharp_plus,
Eigen::VectorXd &  rho 
)
inlinevirtual

Reimplemented in stan::mcmc::rho_inspector_mock_nuts.

Definition at line 182 of file base_nuts.hpp.

Referenced by stan::mcmc::base_nuts< Model, diag_e_metric, expl_leapfrog, BaseRNG >::build_tree(), and stan::mcmc::base_nuts< Model, diag_e_metric, expl_leapfrog, BaseRNG >::transition().

184  {
185  return p_sharp_plus.dot(rho) > 0
186  && p_sharp_minus.dot(rho) > 0;
187  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::get_current_stepsize ( )
inlineinherited

Definition at line 161 of file base_hmc.hpp.

161  {
162  return this->epsilon_;
163  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::get_max_delta ( )
inline

Definition at line 72 of file base_nuts.hpp.

72 { return this->max_deltaH_; }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
int stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::get_max_depth ( )
inline

Definition at line 71 of file base_nuts.hpp.

Referenced by TEST().

71 { return this->max_depth_; }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::get_nominal_stepsize ( )
inlineinherited

Definition at line 157 of file base_hmc.hpp.

Referenced by stan::mcmc::base_hmc< Model, diag_e_metric, expl_leapfrog, BaseRNG >::write_sampler_stepsize().

157  {
158  return this->nom_epsilon_;
159  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::get_sampler_diagnostic_names ( std::vector< std::string > &  model_names,
std::vector< std::string > &  names 
)
inlinevirtualinherited

Reimplemented from stan::mcmc::base_mcmc.

Definition at line 64 of file base_hmc.hpp.

65  {
66  z_.get_param_names(model_names, names);
67  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::get_sampler_diagnostics ( std::vector< double > &  values)
inlinevirtualinherited

Reimplemented from stan::mcmc::base_mcmc.

Definition at line 69 of file base_hmc.hpp.

69  {
70  z_.get_params(values);
71  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::get_sampler_param_names ( std::vector< std::string > &  names)
inlinevirtual

Reimplemented from stan::mcmc::base_mcmc.

Definition at line 166 of file base_nuts.hpp.

166  {
167  names.push_back("stepsize__");
168  names.push_back("treedepth__");
169  names.push_back("n_leapfrog__");
170  names.push_back("divergent__");
171  names.push_back("energy__");
172  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::get_sampler_params ( std::vector< double > &  values)
inlinevirtual

Reimplemented from stan::mcmc::base_mcmc.

Definition at line 174 of file base_nuts.hpp.

174  {
175  values.push_back(this->epsilon_);
176  values.push_back(this->depth_);
177  values.push_back(this->n_leapfrog_);
178  values.push_back(this->divergent_);
179  values.push_back(this->energy_);
180  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::get_stepsize_jitter ( )
inlineinherited

Definition at line 170 of file base_hmc.hpp.

170  {
171  return this->epsilon_jitter_;
172  }
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::init_hamiltonian ( callbacks::logger logger)
inlineinherited

Definition at line 78 of file base_hmc.hpp.

78  {
79  this->hamiltonian_.init(this->z_, logger);
80  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
Hamiltonian< Model, BaseRNG > hamiltonian_
Definition: base_hmc.hpp:184
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::init_stepsize ( callbacks::logger logger)
inlineinherited

Definition at line 83 of file base_hmc.hpp.

83  {
84  ps_point z_init(this->z_);
85 
86  // Skip initialization for extreme step sizes
87  if (this->nom_epsilon_ == 0 || this->nom_epsilon_ > 1e7)
88  return;
89 
90  this->hamiltonian_.sample_p(this->z_, this->rand_int_);
91  this->hamiltonian_.init(this->z_, logger);
92 
93  // Guaranteed to be finite if randomly initialized
94  double H0 = this->hamiltonian_.H(this->z_);
95 
96  this->integrator_.evolve(this->z_, this->hamiltonian_,
97  this->nom_epsilon_,
98  logger);
99 
100  double h = this->hamiltonian_.H(this->z_);
101  if (boost::math::isnan(h))
102  h = std::numeric_limits<double>::infinity();
103 
104  double delta_H = H0 - h;
105 
106  int direction = delta_H > std::log(0.8) ? 1 : -1;
107 
108  while (1) {
109  this->z_.ps_point::operator=(z_init);
110 
111  this->hamiltonian_.sample_p(this->z_, this->rand_int_);
112  this->hamiltonian_.init(this->z_, logger);
113 
114  double H0 = this->hamiltonian_.H(this->z_);
115 
116  this->integrator_.evolve(this->z_, this->hamiltonian_,
117  this->nom_epsilon_,
118  logger);
119 
120  double h = this->hamiltonian_.H(this->z_);
121  if (boost::math::isnan(h))
122  h = std::numeric_limits<double>::infinity();
123 
124  double delta_H = H0 - h;
125 
126  if ((direction == 1) && !(delta_H > std::log(0.8)))
127  break;
128  else if ((direction == -1) && !(delta_H < std::log(0.8)))
129  break;
130  else
131  this->nom_epsilon_
132  = direction == 1
133  ? 2.0 * this->nom_epsilon_
134  : 0.5 * this->nom_epsilon_;
135 
136  if (this->nom_epsilon_ > 1e7)
137  throw std::runtime_error("Posterior is improper. "
138  "Please check your model.");
139  if (this->nom_epsilon_ == 0)
140  throw std::runtime_error("No acceptably small step size could "
141  "be found. Perhaps the posterior is "
142  "not continuous?");
143  }
144 
145  this->z_.ps_point::operator=(z_init);
146  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
bool isnan(const stan::math::var &v)
Definition: boost_isnan.hpp:20
Hamiltonian< Model, BaseRNG > hamiltonian_
Definition: base_hmc.hpp:184
Integrator< Hamiltonian< Model, BaseRNG > > integrator_
Definition: base_hmc.hpp:183
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::sample_stepsize ( )
inlineinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::seed ( const Eigen::VectorXd &  q)
inlineinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::set_max_delta ( double  d)
inline

Definition at line 67 of file base_nuts.hpp.

67  {
68  max_deltaH_ = d;
69  }
Float_t d
Definition: plot.C:236
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::set_max_depth ( int  d)
inline

Definition at line 62 of file base_nuts.hpp.

62  {
63  if (d > 0)
64  max_depth_ = d;
65  }
Float_t d
Definition: plot.C:236
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::set_metric ( const Eigen::MatrixXd &  inv_e_metric)
inline

Definition at line 54 of file base_nuts.hpp.

54  {
55  this->z_.set_metric(inv_e_metric);
56  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::set_metric ( const Eigen::VectorXd &  inv_e_metric)
inline

Definition at line 58 of file base_nuts.hpp.

58  {
59  this->z_.set_metric(inv_e_metric);
60  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
virtual void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::set_nominal_stepsize ( double  e)
inlinevirtualinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
virtual void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::set_stepsize_jitter ( double  j)
inlinevirtualinherited

Definition at line 165 of file base_hmc.hpp.

Referenced by TEST().

165  {
166  if (j > 0 && j < 1)
167  epsilon_jitter_ = j;
168  }
const double j
Definition: BetheBloch.cxx:29
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
sample stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::transition ( sample init_sample,
callbacks::logger logger 
)
inlinevirtual

Implements stan::mcmc::base_mcmc.

Definition at line 75 of file base_nuts.hpp.

Referenced by TEST().

75  {
76  // Initialize the algorithm
77  this->sample_stepsize();
78 
79  this->seed(init_sample.cont_params());
80 
81  this->hamiltonian_.sample_p(this->z_, this->rand_int_);
82  this->hamiltonian_.init(this->z_, logger);
83 
84  ps_point z_plus(this->z_);
85  ps_point z_minus(z_plus);
86 
87  ps_point z_sample(z_plus);
88  ps_point z_propose(z_plus);
89 
90  Eigen::VectorXd p_sharp_plus = this->hamiltonian_.dtau_dp(this->z_);
91  Eigen::VectorXd p_sharp_dummy = p_sharp_plus;
92  Eigen::VectorXd p_sharp_minus = p_sharp_plus;
93  Eigen::VectorXd rho = this->z_.p;
94 
95  double log_sum_weight = 0; // log(exp(H0 - H0))
96  double H0 = this->hamiltonian_.H(this->z_);
97  int n_leapfrog = 0;
98  double sum_metro_prob = 0;
99 
100  // Build a trajectory until the NUTS criterion is no longer satisfied
101  this->depth_ = 0;
102  this->divergent_ = false;
103 
104  while (this->depth_ < this->max_depth_) {
105  // Build a new subtree in a random direction
106  Eigen::VectorXd rho_subtree = Eigen::VectorXd::Zero(rho.size());
107  bool valid_subtree = false;
108  double log_sum_weight_subtree
109  = -std::numeric_limits<double>::infinity();
110 
111  if (this->rand_uniform_() > 0.5) {
112  this->z_.ps_point::operator=(z_plus);
113  valid_subtree
114  = build_tree(this->depth_, z_propose,
115  p_sharp_dummy, p_sharp_plus, rho_subtree,
116  H0, 1, n_leapfrog,
117  log_sum_weight_subtree, sum_metro_prob,
118  logger);
119  z_plus.ps_point::operator=(this->z_);
120  } else {
121  this->z_.ps_point::operator=(z_minus);
122  valid_subtree
123  = build_tree(this->depth_, z_propose,
124  p_sharp_dummy, p_sharp_minus, rho_subtree,
125  H0, -1, n_leapfrog,
126  log_sum_weight_subtree, sum_metro_prob,
127  logger);
128  z_minus.ps_point::operator=(this->z_);
129  }
130 
131  if (!valid_subtree) break;
132 
133  // Sample from an accepted subtree
134  ++(this->depth_);
135 
136  if (log_sum_weight_subtree > log_sum_weight) {
137  z_sample = z_propose;
138  } else {
139  double accept_prob
140  = std::exp(log_sum_weight_subtree - log_sum_weight);
141  if (this->rand_uniform_() < accept_prob)
142  z_sample = z_propose;
143  }
144 
145  log_sum_weight
146  = math::log_sum_exp(log_sum_weight, log_sum_weight_subtree);
147 
148  // Break when NUTS criterion is no longer satisfied
149  rho += rho_subtree;
150  if (!compute_criterion(p_sharp_minus, p_sharp_plus, rho))
151  break;
152  }
153 
154  this->n_leapfrog_ = n_leapfrog;
155 
156  // Compute average acceptance probabilty across entire trajectory,
157  // even over subtrees that may have been rejected
158  double accept_prob
159  = sum_metro_prob / static_cast<double>(n_leapfrog);
160 
161  this->z_.ps_point::operator=(z_sample);
162  this->energy_ = this->hamiltonian_.H(this->z_);
163  return sample(this->z_.q, -this->z_.V, accept_prob);
164  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
fvar< T > log_sum_exp(const std::vector< fvar< T > > &v)
Definition: log_sum_exp.hpp:12
bool build_tree(int depth, ps_point &z_propose, Eigen::VectorXd &p_sharp_left, Eigen::VectorXd &p_sharp_right, Eigen::VectorXd &rho, double H0, double sign, int &n_leapfrog, double &log_sum_weight, double &sum_metro_prob, callbacks::logger &logger)
Definition: base_nuts.hpp:206
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
void seed(const Eigen::VectorXd &q)
Definition: base_hmc.hpp:73
virtual bool compute_criterion(Eigen::VectorXd &p_sharp_minus, Eigen::VectorXd &p_sharp_plus, Eigen::VectorXd &rho)
Definition: base_nuts.hpp:182
boost::uniform_01< BaseRNG & > rand_uniform_
Definition: base_hmc.hpp:189
void Zero()
Hamiltonian< Model, BaseRNG > hamiltonian_
Definition: base_hmc.hpp:184
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::write_sampler_metric ( callbacks::writer writer)
inlineinherited

write elements of mass matrix

Definition at line 51 of file base_hmc.hpp.

Referenced by stan::mcmc::base_hmc< Model, diag_e_metric, expl_leapfrog, BaseRNG >::write_sampler_state().

51  {
52  z_.write_metric(writer);
53  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::write_sampler_state ( callbacks::writer writer)
inlinevirtualinherited

write stepsize and elements of mass matrix

Reimplemented from stan::mcmc::base_mcmc.

Definition at line 59 of file base_hmc.hpp.

Referenced by TEST().

59  {
60  write_sampler_stepsize(writer);
61  write_sampler_metric(writer);
62  }
void write_sampler_metric(callbacks::writer &writer)
Definition: base_hmc.hpp:51
void write_sampler_stepsize(callbacks::writer &writer)
Definition: base_hmc.hpp:41
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
void stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::write_sampler_stepsize ( callbacks::writer writer)
inlineinherited

format and write stepsize

Definition at line 41 of file base_hmc.hpp.

Referenced by stan::mcmc::base_hmc< Model, diag_e_metric, expl_leapfrog, BaseRNG >::write_sampler_state().

41  {
42  std::stringstream nominal_stepsize;
43  nominal_stepsize << "Step size = " << get_nominal_stepsize();
44  writer(nominal_stepsize.str());
45  }
double get_nominal_stepsize()
Definition: base_hmc.hpp:157
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
Hamiltonian<Model, BaseRNG>::PointType& stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::z ( )
inlineinherited

Definition at line 148 of file base_hmc.hpp.

Referenced by TEST().

148  {
149  return z_;
150  }
Hamiltonian< Model, BaseRNG >::PointType z_
Definition: base_hmc.hpp:182

Member Data Documentation

template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
int stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::depth_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
bool stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::divergent_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::energy_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::epsilon_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::epsilon_jitter_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
Hamiltonian<Model, BaseRNG> stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::hamiltonian_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
Integrator<Hamiltonian<Model, BaseRNG> > stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::integrator_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::max_deltaH_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
int stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::max_depth_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
int stan::mcmc::base_nuts< Model, Hamiltonian, Integrator, BaseRNG >::n_leapfrog_
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
double stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::nom_epsilon_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
BaseRNG& stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::rand_int_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
boost::uniform_01<BaseRNG&> stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::rand_uniform_
protectedinherited
template<class Model, template< class, class > class Hamiltonian, template< class > class Integrator, class BaseRNG>
Hamiltonian<Model, BaseRNG>::PointType stan::mcmc::base_hmc< Model, Hamiltonian, Integrator, BaseRNG >::z_
protectedinherited

The documentation for this class was generated from the following file: