Classes | Functions
stan::services::util Namespace Reference

Classes

class  gq_writer
 
class  mcmc_writer
 

Functions

boost::ecuyer1988 create_rng (unsigned int seed, unsigned int chain)
 
stan::io::dump create_unit_e_dense_inv_metric (size_t num_params)
 
stan::io::dump create_unit_e_diag_inv_metric (size_t num_params)
 
void experimental_message (stan::callbacks::logger &logger)
 
template<class Model , class RNG >
void generate_transitions (stan::mcmc::base_mcmc &sampler, int num_iterations, int start, int finish, int num_thin, int refresh, bool save, bool warmup, util::mcmc_writer &mcmc_writer, stan::mcmc::sample &init_s, Model &model, RNG &base_rng, callbacks::interrupt &callback, callbacks::logger &logger)
 
template<bool Jacobian = true, class Model , class RNG >
std::vector< double > initialize (Model &model, stan::io::var_context &init, RNG &rng, double init_radius, bool print_timing, stan::callbacks::logger &logger, stan::callbacks::writer &init_writer)
 
Eigen::MatrixXd read_dense_inv_metric (stan::io::var_context &init_context, size_t num_params, callbacks::logger &logger)
 
Eigen::VectorXd read_diag_inv_metric (stan::io::var_context &init_context, size_t num_params, callbacks::logger &logger)
 
template<class Sampler , class Model , class RNG >
void run_adaptive_sampler (Sampler &sampler, Model &model, std::vector< double > &cont_vector, int num_warmup, int num_samples, int num_thin, int refresh, bool save_warmup, RNG &rng, callbacks::interrupt &interrupt, callbacks::logger &logger, callbacks::writer &sample_writer, callbacks::writer &diagnostic_writer)
 
template<class Model , class RNG >
void run_sampler (stan::mcmc::base_mcmc &sampler, Model &model, std::vector< double > &cont_vector, int num_warmup, int num_samples, int num_thin, int refresh, bool save_warmup, RNG &rng, callbacks::interrupt &interrupt, callbacks::logger &logger, callbacks::writer &sample_writer, callbacks::writer &diagnostic_writer)
 
void validate_dense_inv_metric (const Eigen::MatrixXd &inv_metric, callbacks::logger &logger)
 
void validate_diag_inv_metric (const Eigen::VectorXd &inv_metric, callbacks::logger &logger)
 

Function Documentation

boost::ecuyer1988 stan::services::util::create_rng ( unsigned int  seed,
unsigned int  chain 
)
inline

Creates a pseudo random number generator from a random seed and a chain id by initializing the PRNG with the seed and then advancing past pow(2, 50) times the chain ID draws to ensure different chains sample from different segments of the pseudo random number sequence.

Chain IDs should be kept to larger values than one to ensure that the draws used to initialized transformed data are not duplicated.

Parameters
[in]seedthe random seed
[in]chainthe chain id
Returns
a boost::ecuyer1988 instance

Definition at line 25 of file create_rng.hpp.

Referenced by stan::services::optimize::bfgs(), stan::services::diagnose::diagnose(), stan::services::sample::fixed_param(), stan::services::experimental::advi::fullrank(), stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_nuts_unit_e(), stan::services::sample::hmc_nuts_unit_e_adapt(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_dense_e_adapt(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_diag_e_adapt(), stan::services::sample::hmc_static_unit_e(), stan::services::sample::hmc_static_unit_e_adapt(), stan::services::optimize::lbfgs(), stan::services::experimental::advi::meanfield(), stan::services::optimize::newton(), stan::services::standalone_generate(), TEST(), and TEST_F().

26  {
27  using boost::uintmax_t;
28  static uintmax_t DISCARD_STRIDE = static_cast<uintmax_t>(1) << 50;
29  boost::ecuyer1988 rng(seed);
30  rng.discard(DISCARD_STRIDE * chain);
31  return rng;
32  }
unsigned int seed
Definition: runWimpSim.h:102
chain
Check that an output directory exists.
stan::io::dump stan::services::util::create_unit_e_dense_inv_metric ( size_t  num_params)
inline

Create a stan::dump object which contains vector "metric" of specified size where all elements are ones.

Parameters
[in]num_paramsexpected number of denseonal elements
Returns
var_context

Definition at line 21 of file create_unit_e_dense_inv_metric.hpp.

References downloadDefinitionData::dump, MECModelEnuComparisons::i, and stan::math::num_elements().

Referenced by stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_dense_e_adapt(), TEST(), and TEST_F().

21  {
22  Eigen::MatrixXd inv_metric(num_params, num_params);
23  inv_metric.setIdentity();
24  size_t num_elements = num_params * num_params;
25  std::stringstream txt;
26  txt << "inv_metric <- structure(c(";
27  for (size_t i = 0; i < num_elements; i++) {
28  txt << inv_metric(i);
29  if (i < num_elements - 1)
30  txt << ", ";
31  }
32  txt << "),.Dim=c("
33  << num_params
34  << ", "
35  << num_params
36  << "))";
37  return stan::io::dump(txt);
38  }
int num_elements(const T &x)
stan::io::dump stan::services::util::create_unit_e_diag_inv_metric ( size_t  num_params)
inline

Create a stan::dump object which contains vector "metric" of specified size where all elements are ones.

Parameters
[in]num_paramsexpected number of diagonal elements
Returns
var_context

Definition at line 20 of file create_unit_e_diag_inv_metric.hpp.

References downloadDefinitionData::dump, and MECModelEnuComparisons::i.

Referenced by stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_diag_e_adapt(), TEST(), and TEST_F().

20  {
21  std::stringstream txt;
22  txt << "inv_metric <- structure(c(";
23  for (size_t i = 0; i < num_params; ++i) {
24  txt << "1.0";
25  if (i < num_params - 1)
26  txt << ", ";
27  }
28  txt << "),.Dim=c(" << num_params << "))";
29  return stan::io::dump(txt);
30  }
void stan::services::util::experimental_message ( stan::callbacks::logger logger)
inline

Writes an experimental message to the writer. All experimental algorithms should call this function.

Parameters
[in,out]loggerlogger for experimental algorithm message

Definition at line 16 of file experimental_message.hpp.

References stan::callbacks::logger::info().

Referenced by stan::services::experimental::advi::fullrank(), stan::services::experimental::advi::meanfield(), and TEST().

16  {
17  logger.info("------------------------------"
18  "------------------------------");
19  logger.info("EXPERIMENTAL ALGORITHM:");
20  logger.info(" This procedure has not been thoroughly tested"
21  " and may be unstable");
22  logger.info(" or buggy. The interface is subject to change.");
23  logger.info("------------------------------"
24  "------------------------------");
25  logger.info("");
26  logger.info("");
27  }
virtual void info(const std::string &message)
Definition: logger.hpp:47
template<class Model , class RNG >
void stan::services::util::generate_transitions ( stan::mcmc::base_mcmc sampler,
int  num_iterations,
int  start,
int  finish,
int  num_thin,
int  refresh,
bool  save,
bool  warmup,
util::mcmc_writer mcmc_writer,
stan::mcmc::sample init_s,
Model model,
RNG &  base_rng,
callbacks::interrupt callback,
callbacks::logger logger 
)

Generates MCMC transitions.

Template Parameters
Modelmodel class
RNGrandom number generator class
Parameters
[in,out]samplerMCMC sampler used to generate transitions
[in]num_iterationsnumber of MCMC transitions
[in]startstarting iteration number used for printing messages
[in]finishend iteration number used for printing messages
[in]num_thinwhen save is true, a draw will be written to the mcmc_writer every num_thin iterations
[in]refreshnumber of iterations to print a message. If refresh is zero, iteration number messages will not be printed
[in]saveif save is true, the transitions will be written to the mcmc_writer. If false, transitions will not be written
[in]warmupindicates whether these transitions are warmup. Used for printing iteration number messages
[in,out]mcmc_writerwriter to handle mcmc otuput
[in,out]init_sstarts as the initial unconstrained parameter values. When the function completes, this will have the final iteration's unconstrained parameter values
[in]modelmodel
[in,out]base_rngrandom number generator
[in,out]callbackinterrupt callback called once an iteration
[in,out]loggerlogger for messages

Definition at line 41 of file generate_transitions.hpp.

References stan::math::ceil(), stan::callbacks::logger::info(), std::log10(), m, datagram_client::message, stan::mcmc::base_mcmc::transition(), stan::services::util::mcmc_writer::write_diagnostic_params(), and stan::services::util::mcmc_writer::write_sample_params().

Referenced by stan::services::sample::fixed_param(), run_adaptive_sampler(), run_sampler(), and TEST_F().

49  {
50  for (int m = 0; m < num_iterations; ++m) {
51  callback();
52 
53  if (refresh > 0
54  && (start + m + 1 == finish
55  || m == 0
56  || (m + 1) % refresh == 0)) {
57  int it_print_width
58  = std::ceil(std::log10(static_cast<double>(finish)));
59  std::stringstream message;
60  message << "Iteration: ";
61  message << std::setw(it_print_width) << m + 1 + start
62  << " / " << finish;
63  message << " [" << std::setw(3)
64  << static_cast<int>( (100.0 * (start + m + 1)) / finish )
65  << "%] ";
66  message << (warmup ? " (Warmup)" : " (Sampling)");
67 
68  logger.info(message);
69  }
70 
71  init_s = sampler.transition(init_s, logger);
72 
73  if (save && ((m % num_thin) == 0)) {
74  mcmc_writer.write_sample_params(base_rng, init_s, sampler, model);
75  mcmc_writer.write_diagnostic_params(init_s, sampler);
76  }
77  }
78  }
virtual sample transition(sample &init_sample, callbacks::logger &logger)=0
void save(bool IsData, bool Bad)
Definition: cellShifts.C:453
T log10(T number)
Definition: d0nt_math.hpp:120
const XML_Char XML_Content * model
Definition: expat.h:151
void refresh()
Definition: show_event.C:21
fvar< T > ceil(const fvar< T > &x)
Definition: ceil.hpp:11
template<bool Jacobian = true, class Model , class RNG >
std::vector<double> stan::services::util::initialize ( Model model,
stan::io::var_context init,
RNG &  rng,
double  init_radius,
bool  print_timing,
stan::callbacks::logger logger,
stan::callbacks::writer init_writer 
)

Returns a valid initial value of the parameters of the model on the unconstrained scale.

For identical inputs (model, init, rng, init_radius), this function will produce the same initialization.

Initialization first tries to use the provided stan::io::var_context, then it will generate random uniform values from -init_radius to +init_radius for missing parameters.

When the var_context provides all variables or the init_radius is 0, this function will only evaluate the log probability of the model with the unconstrained parameters once to see if it's valid.

When at least some of the initialization is random, it will randomly initialize until it finds a set of unconstrained parameters that are valid or it hits MAX_INIT_TRIES = 100 (hard-coded).

Valid initialization is defined as a finite, non-NaN value for the evaluation of the log probability density function and all its gradients.

Template Parameters
Jacobianindicates whether to include the Jacobian term when evaluating the log density function
Modelthe type of the model class
RNGthe type of the random number generator
Parameters
[in]modelthe model
[in]inita var_context with initial values
[in,out]rngrandom number generator
[in]init_radiusthe radius for generating random values. A value of 0 indicates that the unconstrained parameters (not provided by init) should be initialized with 0.
[in]print_timingindicates whether a timing message should be printed to the logger
[in,out]loggerlogger for messages
[in,out]init_writerinit writer (on the unconstrained scale)
Exceptions
exceptionpassed through from the model if the model has a fatal error (not a std::domain_error)
std::domain_errorif the model can not be initialized and the model does not have a fatal error (only allows for std::domain_error)
Returns
valid unconstrained parameters for the model

Definition at line 68 of file initialize.hpp.

References stan::io::var_context::contains_r(), stan::math::domain_error(), e, stan::io::random_var_context::get_unconstrained(), stan::model::gradient(), stan::callbacks::logger::info(), boost::math::isfinite(), lem_server::msg, getGoodRuns4SAM::n, and stan::math::sum().

Referenced by stan::optimization::BFGSLineSearch< M, QNUpdateType, Scalar, DimAtCompile >::BFGSLineSearch(), stan::services::diagnose::diagnose(), stan::services::sample::fixed_param(), stan::services::experimental::advi::fullrank(), stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_nuts_unit_e(), stan::services::sample::hmc_nuts_unit_e_adapt(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_dense_e_adapt(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_diag_e_adapt(), stan::services::sample::hmc_static_unit_e(), stan::services::sample::hmc_static_unit_e_adapt(), stan::optimization::BFGSLineSearch< M, QNUpdateType, Scalar, DimAtCompile >::initialize(), stan::services::experimental::advi::meanfield(), stan::optimization::BFGSMinimizer< ModelAdaptor< M >, QNUpdateType, Scalar, DimAtCompile >::minimize(), and TEST_F().

75  {
76  std::vector<double> unconstrained;
77  std::vector<int> disc_vector;
78 
79  bool is_fully_initialized = true;
80  bool any_initialized = false;
81  std::vector<std::string> param_names;
82  model.get_param_names(param_names);
83  for (size_t n = 0; n < param_names.size(); n++) {
84  is_fully_initialized &= init.contains_r(param_names[n]);
85  any_initialized |= init.contains_r(param_names[n]);
86  }
87 
88  bool is_initialized_with_zero = init_radius == 0.0;
89 
90  int MAX_INIT_TRIES = is_fully_initialized || is_initialized_with_zero
91  ? 1 : 100;
92  int num_init_tries = 0;
93  for (; num_init_tries < MAX_INIT_TRIES; num_init_tries++) {
94  std::stringstream msg;
95  try {
97  random_context(model, rng, init_radius, is_initialized_with_zero);
98 
99  if (!any_initialized) {
100  unconstrained = random_context.get_unconstrained();
101  } else {
102  stan::io::chained_var_context context(init, random_context);
103 
104  model.transform_inits(context,
105  disc_vector,
106  unconstrained,
107  &msg);
108  }
109  } catch (std::domain_error& e) {
110  if (msg.str().length() > 0)
111  logger.info(msg);
112  logger.info("Rejecting initial value:");
113  logger.info(" Error evaluating the log probability"
114  " at the initial value.");
115  logger.info(e.what());
116  continue;
117  } catch (std::exception& e) {
118  if (msg.str().length() > 0)
119  logger.info(msg);
120  logger.info("Unrecoverable error evaluating the log probability"
121  " at the initial value.");
122  logger.info(e.what());
123  throw;
124  }
125 
126  msg.str("");
127  double log_prob(0);
128  try {
129  // we evaluate the log_prob function with propto=false
130  // because we're evaluating with `double` as the type of
131  // the parameters.
132  log_prob = model.template log_prob<false, Jacobian>
133  (unconstrained, disc_vector, &msg);
134  if (msg.str().length() > 0)
135  logger.info(msg);
136  } catch (std::domain_error& e) {
137  if (msg.str().length() > 0)
138  logger.info(msg);
139  logger.info("Rejecting initial value:");
140  logger.info(" Error evaluating the log probability"
141  " at the initial value.");
142  logger.info(e.what());
143  continue;
144  } catch (std::exception& e) {
145  if (msg.str().length() > 0)
146  logger.info(msg);
147  logger.info("Unrecoverable error evaluating the log probability"
148  " at the initial value.");
149  logger.info(e.what());
150  throw;
151  }
152  if (!boost::math::isfinite(log_prob)) {
153  logger.info("Rejecting initial value:");
154  logger.info(" Log probability evaluates to log(0),"
155  " i.e. negative infinity.");
156  logger.info(" Stan can't start sampling from this"
157  " initial value.");
158  continue;
159  }
160  std::stringstream log_prob_msg;
161  std::vector<double> gradient;
162  clock_t start_check = clock();
163  try {
164  // we evaluate this with propto=true since we're
165  // evaluating with autodiff variables
166  log_prob = stan::model::log_prob_grad<true, Jacobian>
167  (model, unconstrained, disc_vector,
168  gradient, &log_prob_msg);
169  } catch (const std::exception& e) {
170  if (log_prob_msg.str().length() > 0)
171  logger.info(log_prob_msg);
172  logger.info(e.what());
173  throw;
174  }
175  clock_t end_check = clock();
176  double deltaT = static_cast<double>(end_check - start_check)
177  / CLOCKS_PER_SEC;
178  if (log_prob_msg.str().length() > 0)
179  logger.info(log_prob_msg);
180 
181  bool gradient_ok = boost::math::isfinite(stan::math::sum(gradient));
182 
183  if (!gradient_ok) {
184  logger.info("Rejecting initial value:");
185  logger.info(" Gradient evaluated at the initial value"
186  " is not finite.");
187  logger.info(" Stan can't start sampling from this"
188  " initial value.");
189  }
190  if (gradient_ok && print_timing) {
191  logger.info("");
192  std::stringstream msg1;
193  msg1 << "Gradient evaluation took " << deltaT << " seconds";
194  logger.info(msg1);
195 
196  std::stringstream msg2;
197  msg2 << "1000 transitions using 10 leapfrog steps"
198  << " per transition would take"
199  << " " << 1e4 * deltaT << " seconds.";
200  logger.info(msg2);
201 
202  logger.info("Adjust your expectations accordingly!");
203  logger.info("");
204  logger.info("");
205  }
206  if (gradient_ok) {
207  init_writer(unconstrained);
208  return unconstrained;
209  }
210  }
211 
212  if (!is_initialized_with_zero) {
213  logger.info("");
214  std::stringstream msg;
215  msg << "Initialization between (-" << init_radius
216  << ", " << init_radius << ") failed after"
217  << " " << MAX_INIT_TRIES << " attempts. ";
218  logger.info(msg);
219  logger.info(" Try specifying initial values,"
220  " reducing ranges of constrained values,"
221  " or reparameterizing the model.");
222  }
223  throw std::domain_error("Initialization failed.");
224 }
fvar< T > sum(const std::vector< fvar< T > > &m)
Definition: sum.hpp:20
bool isfinite(const stan::math::var &v)
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
void gradient(const M &model, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &x, double &f, Eigen::Matrix< double, Eigen::Dynamic, 1 > &grad_f, std::ostream *msgs=0)
Definition: gradient.hpp:15
virtual bool contains_r(const std::string &name) const =0
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
const XML_Char * context
Definition: expat.h:434
virtual void info(const std::string &message)
Definition: logger.hpp:47
Float_t e
Definition: plot.C:35
const XML_Char XML_Content * model
Definition: expat.h:151
Eigen::MatrixXd stan::services::util::read_dense_inv_metric ( stan::io::var_context init_context,
size_t  num_params,
callbacks::logger logger 
)
inline

Extract dense inverse Euclidean metric from a var_context object.

Parameters
[in]init_contexta var_context with array of initial values
[in]num_paramsexpected number of row, column elements
[in,out]loggerLogger for messages
Exceptions
std::domain_errorif cannot read the Euclidean metric
Returns
inv_metric

Definition at line 27 of file read_dense_inv_metric.hpp.

References stan::math::domain_error(), e, stan::callbacks::logger::error(), stan::math::to_matrix(), stan::io::var_context::to_vec(), stan::io::var_context::validate_dims(), and stan::io::var_context::vals_r().

Referenced by stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_dense_e_adapt(), and TEST().

29  {
30  Eigen::MatrixXd inv_metric;
31  try {
32  init_context.validate_dims(
33  "read dense inv metric", "inv_metric", "matrix",
34  init_context.to_vec(num_params, num_params));
35  std::vector<double> dense_vals =
36  init_context.vals_r("inv_metric");
37  inv_metric =
38  stan::math::to_matrix(dense_vals, num_params, num_params);
39  } catch (const std::exception& e) {
40  logger.error("Cannot get inverse metric from input file.");
41  logger.error("Caught exception: ");
42  logger.error(e.what());
43  throw std::domain_error("Initialization failure");
44  }
45 
46  return inv_metric;
47  }
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
static std::vector< size_t > to_vec()
virtual std::vector< double > vals_r(const std::string &name) const =0
void validate_dims(const std::string &stage, const std::string &name, const std::string &base_type, const std::vector< size_t > &dims_declared) const
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > to_matrix(const Eigen::Matrix< T, R, C > &x)
Definition: to_matrix.hpp:25
Float_t e
Definition: plot.C:35
Eigen::VectorXd stan::services::util::read_diag_inv_metric ( stan::io::var_context init_context,
size_t  num_params,
callbacks::logger logger 
)
inline

Extract diagonal values for an inverse Euclidean metric from a var_context object.

Parameters
[in]init_contexta var_context with initial values
[in]num_paramsexpected number of diagonal elements
[in,out]loggerLogger for messages
Exceptions
std::domain_errorif the Euclidean metric is invalid
Returns
inv_metric vector of diagonal values

Definition at line 28 of file read_diag_inv_metric.hpp.

References stan::math::domain_error(), e, stan::callbacks::logger::error(), MECModelEnuComparisons::i, stan::io::var_context::to_vec(), stan::io::var_context::validate_dims(), and stan::io::var_context::vals_r().

Referenced by stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_diag_e_adapt(), and TEST().

30  {
31  Eigen::VectorXd inv_metric(num_params);
32  try {
33  init_context.validate_dims(
34  "read diag inv metric", "inv_metric", "vector_d",
35  init_context.to_vec(num_params));
36  std::vector<double> diag_vals =
37  init_context.vals_r("inv_metric");
38  for (size_t i=0; i < num_params; i++) {
39  inv_metric(i) = diag_vals[i];
40  }
41  } catch (const std::exception& e) {
42  logger.error("Cannot get inverse Euclidean metric from input file.");
43  logger.error("Caught exception: ");
44  logger.error(e.what());
45  throw std::domain_error("Initialization failure");
46  }
47  return inv_metric;
48  }
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
static std::vector< size_t > to_vec()
virtual std::vector< double > vals_r(const std::string &name) const =0
void validate_dims(const std::string &stage, const std::string &name, const std::string &base_type, const std::vector< size_t > &dims_declared) const
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
Float_t e
Definition: plot.C:35
template<class Sampler , class Model , class RNG >
void stan::services::util::run_adaptive_sampler ( Sampler &  sampler,
Model model,
std::vector< double > &  cont_vector,
int  num_warmup,
int  num_samples,
int  num_thin,
int  refresh,
bool  save_warmup,
RNG &  rng,
callbacks::interrupt interrupt,
callbacks::logger logger,
callbacks::writer sample_writer,
callbacks::writer diagnostic_writer 
)

Runs the sampler with adaptation.

Template Parameters
SamplerType of adaptive sampler.
ModelType of model
RNGType of random number generator
Parameters
[in,out]samplerthe mcmc sampler to use on the model
[in]modelthe model concept to use for computing log probability
[in]cont_vectorinitial parameter values
[in]num_warmupnumber of warmup draws
[in]num_samplesnumber of post warmup draws
[in]num_thinnumber to thin the draws. Must be greater than or equal to 1.
[in]refreshcontrols output to the logger
[in]save_warmupindicates whether the warmup draws should be sent to the sample writer
[in,out]rngrandom number generator
[in,out]interruptinterrupt callback
[in,out]loggerlogger for messages
[in,out]sample_writerwriter for draws
[in,out]diagnostic_writerwriter for diagnostic information

Definition at line 38 of file run_adaptive_sampler.hpp.

References e, febshutoff_auto::end, generate_transitions(), stan::callbacks::logger::info(), febshutoff_auto::start, stan::services::util::mcmc_writer::write_adapt_finish(), stan::services::util::mcmc_writer::write_diagnostic_names(), stan::services::util::mcmc_writer::write_sample_names(), and stan::services::util::mcmc_writer::write_timing().

Referenced by stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_nuts_unit_e_adapt(), stan::services::sample::hmc_static_dense_e_adapt(), stan::services::sample::hmc_static_diag_e_adapt(), stan::services::sample::hmc_static_unit_e_adapt(), and TEST_F().

46  {
47  Eigen::Map<Eigen::VectorXd> cont_params(cont_vector.data(),
48  cont_vector.size());
49 
50  sampler.engage_adaptation();
51  try {
52  sampler.z().q = cont_params;
53  sampler.init_stepsize(logger);
54  } catch (const std::exception& e) {
55  logger.info("Exception initializing step size.");
56  logger.info(e.what());
57  return;
58  }
59 
60  services::util::mcmc_writer
61  writer(sample_writer, diagnostic_writer, logger);
62  stan::mcmc::sample s(cont_params, 0, 0);
63 
64  // Headers
65  writer.write_sample_names(s, sampler, model);
66  writer.write_diagnostic_names(s, sampler, model);
67 
68  clock_t start = clock();
69  util::generate_transitions(sampler, num_warmup, 0,
70  num_warmup + num_samples, num_thin,
71  refresh, save_warmup, true,
72  writer,
73  s, model, rng,
74  interrupt, logger);
75  clock_t end = clock();
76  double warm_delta_t = static_cast<double>(end - start) / CLOCKS_PER_SEC;
77 
78  sampler.disengage_adaptation();
79  writer.write_adapt_finish(sampler);
80  sampler.write_sampler_state(sample_writer);
81 
82  start = clock();
83  util::generate_transitions(sampler, num_samples, num_warmup,
84  num_warmup + num_samples, num_thin,
85  refresh, true, false,
86  writer,
87  s, model, rng,
88  interrupt, logger);
89  end = clock();
90  double sample_delta_t
91  = static_cast<double>(end - start) / CLOCKS_PER_SEC;
92 
93  writer.write_timing(warm_delta_t, sample_delta_t);
94  }
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
const XML_Char * s
Definition: expat.h:262
Float_t e
Definition: plot.C:35
const XML_Char XML_Content * model
Definition: expat.h:151
void generate_transitions(stan::mcmc::base_mcmc &sampler, int num_iterations, int start, int finish, int num_thin, int refresh, bool save, bool warmup, util::mcmc_writer &mcmc_writer, stan::mcmc::sample &init_s, Model &model, RNG &base_rng, callbacks::interrupt &callback, callbacks::logger &logger)
void refresh()
Definition: show_event.C:21
template<class Model , class RNG >
void stan::services::util::run_sampler ( stan::mcmc::base_mcmc sampler,
Model model,
std::vector< double > &  cont_vector,
int  num_warmup,
int  num_samples,
int  num_thin,
int  refresh,
bool  save_warmup,
RNG &  rng,
callbacks::interrupt interrupt,
callbacks::logger logger,
callbacks::writer sample_writer,
callbacks::writer diagnostic_writer 
)

Runs the sampler without adaptation.

Template Parameters
ModelType of model
RNGType of random number generator
Parameters
[in,out]samplerthe mcmc sampler to use on the model
[in]modelthe model concept to use for computing log probability
[in]cont_vectorinitial parameter values
[in]num_warmupnumber of warmup draws
[in]num_samplesnumber of post warmup draws
[in]num_thinnumber to thin the draws. Must be greater than or equal to 1.
[in]refreshcontrols output to the logger
[in]save_warmupindicates whether the warmup draws should be sent to the sample writer
[in,out]rngrandom number generator
[in,out]interruptinterrupt callback
[in,out]loggerlogger for messages
[in,out]sample_writerwriter for draws
[in,out]diagnostic_writerwriter for diagnostic information

Definition at line 36 of file run_sampler.hpp.

References febshutoff_auto::end, generate_transitions(), febshutoff_auto::start, and stan::mcmc::base_mcmc::write_sampler_state().

Referenced by stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_unit_e(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_unit_e(), and TEST_F().

43  {
44  Eigen::Map<Eigen::VectorXd> cont_params(cont_vector.data(),
45  cont_vector.size());
46  services::util::mcmc_writer
47  writer(sample_writer, diagnostic_writer, logger);
48  stan::mcmc::sample s(cont_params, 0, 0);
49 
50  // Headers
51  writer.write_sample_names(s, sampler, model);
52  writer.write_diagnostic_names(s, sampler, model);
53 
54  clock_t start = clock();
55  util::generate_transitions(sampler, num_warmup, 0,
56  num_warmup + num_samples, num_thin,
57  refresh, save_warmup, true,
58  writer,
59  s, model, rng,
60  interrupt, logger);
61  clock_t end = clock();
62  double warm_delta_t = static_cast<double>(end - start) / CLOCKS_PER_SEC;
63 
64  writer.write_adapt_finish(sampler);
65  sampler.write_sampler_state(sample_writer);
66 
67  start = clock();
68  util::generate_transitions(sampler, num_samples, num_warmup,
69  num_warmup + num_samples, num_thin,
70  refresh, true, false,
71  writer,
72  s, model, rng,
73  interrupt, logger);
74  end = clock();
75  double sample_delta_t
76  = static_cast<double>(end - start) / CLOCKS_PER_SEC;
77 
78  writer.write_timing(warm_delta_t, sample_delta_t);
79  }
const XML_Char * s
Definition: expat.h:262
virtual void write_sampler_state(callbacks::writer &writer)
Definition: base_mcmc.hpp:28
const XML_Char XML_Content * model
Definition: expat.h:151
void generate_transitions(stan::mcmc::base_mcmc &sampler, int num_iterations, int start, int finish, int num_thin, int refresh, bool save, bool warmup, util::mcmc_writer &mcmc_writer, stan::mcmc::sample &init_s, Model &model, RNG &base_rng, callbacks::interrupt &callback, callbacks::logger &logger)
void refresh()
Definition: show_event.C:21
void stan::services::util::validate_dense_inv_metric ( const Eigen::MatrixXd &  inv_metric,
callbacks::logger logger 
)
inline

Validate that dense inverse Euclidean metric is positive definite

Parameters
[in]inv_metricinverse Euclidean metric
[in,out]loggerLogger for messages
Exceptions
std::domain_errorif matrix is not positive definite

Definition at line 20 of file validate_dense_inv_metric.hpp.

References stan::math::check_pos_definite(), stan::math::domain_error(), e, and stan::callbacks::logger::error().

Referenced by stan::services::sample::hmc_nuts_dense_e(), stan::services::sample::hmc_nuts_dense_e_adapt(), stan::services::sample::hmc_static_dense_e(), stan::services::sample::hmc_static_dense_e_adapt(), and TEST().

21  {
22  try {
23  stan::math::check_pos_definite("check_pos_definite",
24  "inv_metric", inv_metric);
25  } catch (const std::domain_error& e) {
26  logger.error("Inverse Euclidean metric not positive definite.");
27  throw std::domain_error("Initialization failure");
28  }
29  }
void check_pos_definite(const char *function, const char *name, const Eigen::Matrix< T_y,-1,-1 > &y)
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
Float_t e
Definition: plot.C:35
void stan::services::util::validate_diag_inv_metric ( const Eigen::VectorXd &  inv_metric,
callbacks::logger logger 
)
inline

Validate that diag inverse Euclidean metric is positive definite

Parameters
[in]inv_metricinverse Euclidean metric
[in,out]loggerLogger for messages
Exceptions
std::domain_errorif matrix is not positive definite

Definition at line 20 of file validate_diag_inv_metric.hpp.

References stan::math::check_finite(), stan::math::check_positive(), stan::math::domain_error(), e, and stan::callbacks::logger::error().

Referenced by stan::services::sample::hmc_nuts_diag_e(), stan::services::sample::hmc_nuts_diag_e_adapt(), stan::services::sample::hmc_static_diag_e(), stan::services::sample::hmc_static_diag_e_adapt(), and TEST().

21  {
22  try {
23  stan::math::check_finite("check_finite",
24  "inv_metric", inv_metric);
25  stan::math::check_positive("check_positive",
26  "inv_metric", inv_metric);
27  } catch (const std::domain_error& e) {
28  logger.error("Inverse Euclidean metric not positive definite.");
29  throw std::domain_error("Initialization failure");
30  }
31  }
void check_finite(const char *function, const char *name, const T_y &y)
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
void check_positive(const char *function, const char *name, const T_y &y)
Float_t e
Definition: plot.C:35