lbfgs.hpp
Go to the documentation of this file.
1 #ifndef STAN_SERVICES_OPTIMIZE_LBFGS_HPP
2 #define STAN_SERVICES_OPTIMIZE_LBFGS_HPP
3 
14 #include <fstream>
15 #include <iostream>
16 #include <iomanip>
17 #include <string>
18 #include <vector>
19 
20 namespace stan {
21  namespace services {
22  namespace optimize {
23 
24  /**
25  * Runs the L-BFGS algorithm for a model.
26  *
27  * @tparam Model A model implementation
28  * @param[in] model Input model to test (with data already instantiated)
29  * @param[in] init var context for initialization
30  * @param[in] random_seed random seed for the random number generator
31  * @param[in] chain chain id to advance the pseudo random number generator
32  * @param[in] init_radius radius to initialize
33  * @param[in] history_size amount of history to keep for L-BFGS
34  * @param[in] init_alpha line search step size for first iteration
35  * @param[in] tol_obj convergence tolerance on absolute changes in
36  * objective function value
37  * @param[in] tol_rel_obj convergence tolerance on relative changes
38  * in objective function value
39  * @param[in] tol_grad convergence tolerance on the norm of the gradient
40  * @param[in] tol_rel_grad convergence tolerance on the relative norm of
41  * the gradient
42  * @param[in] tol_param convergence tolerance on changes in parameter
43  * value
44  * @param[in] num_iterations maximum number of iterations
45  * @param[in] save_iterations indicates whether all the interations should
46  * be saved to the parameter_writer
47  * @param[in] refresh how often to write output to logger
48  * @param[in,out] interrupt callback to be called every iteration
49  * @param[in,out] logger Logger for messages
50  * @param[in,out] init_writer Writer callback for unconstrained inits
51  * @param[in,out] parameter_writer output for parameter values
52  * @return error_codes::OK if successful
53  */
54  template <class Model>
56  unsigned int random_seed, unsigned int chain,
57  double init_radius, int history_size, double init_alpha,
58  double tol_obj, double tol_rel_obj, double tol_grad,
59  double tol_rel_grad, double tol_param, int num_iterations,
60  bool save_iterations, int refresh,
61  callbacks::interrupt& interrupt,
62  callbacks::logger& logger,
63  callbacks::writer& init_writer,
64  callbacks::writer& parameter_writer) {
65  boost::ecuyer1988 rng = util::create_rng(random_seed, chain);
66 
67  std::vector<int> disc_vector;
68  std::vector<double> cont_vector
69  = util::initialize<false>(model, init, rng, init_radius, false,
70  logger, init_writer);
71 
72  std::stringstream lbfgs_ss;
75  Optimizer lbfgs(model, cont_vector, disc_vector, &lbfgs_ss);
76  lbfgs.get_qnupdate().set_history_size(history_size);
77  lbfgs._ls_opts.alpha0 = init_alpha;
78  lbfgs._conv_opts.tolAbsF = tol_obj;
79  lbfgs._conv_opts.tolRelF = tol_rel_obj;
80  lbfgs._conv_opts.tolAbsGrad = tol_grad;
81  lbfgs._conv_opts.tolRelGrad = tol_rel_grad;
82  lbfgs._conv_opts.tolAbsX = tol_param;
83  lbfgs._conv_opts.maxIts = num_iterations;
84 
85  double lp = lbfgs.logp();
86 
87  std::stringstream initial_msg;
88  initial_msg << "Initial log joint probability = " << lp;
89  logger.info(initial_msg);
90 
91  std::vector<std::string> names;
92  names.push_back("lp__");
93  model.constrained_param_names(names, true, true);
94  parameter_writer(names);
95 
96  if (save_iterations) {
97  std::vector<double> values;
98  std::stringstream msg;
99  model.write_array(rng, cont_vector, disc_vector, values,
100  true, true, &msg);
101  if (msg.str().length() > 0)
102  logger.info(msg);
103 
104  values.insert(values.begin(), lp);
105  parameter_writer(values);
106  }
107  int ret = 0;
108 
109  while (ret == 0) {
110  interrupt();
111  if (refresh > 0
112  && (lbfgs.iter_num() == 0
113  || ((lbfgs.iter_num() + 1) % refresh == 0)))
114  logger.info(" Iter"
115  " log prob"
116  " ||dx||"
117  " ||grad||"
118  " alpha"
119  " alpha0"
120  " # evals"
121  " Notes ");
122 
123  ret = lbfgs.step();
124  lp = lbfgs.logp();
125  lbfgs.params_r(cont_vector);
126 
127  if (refresh > 0
128  && (ret != 0
129  || !lbfgs.note().empty()
130  || lbfgs.iter_num() == 0
131  || ((lbfgs.iter_num() + 1) % refresh == 0))) {
132  std::stringstream msg;
133  msg << " " << std::setw(7) << lbfgs.iter_num() << " ";
134  msg << " " << std::setw(12) << std::setprecision(6)
135  << lp << " ";
136  msg << " " << std::setw(12) << std::setprecision(6)
137  << lbfgs.prev_step_size() << " ";
138  msg << " " << std::setw(12) << std::setprecision(6)
139  << lbfgs.curr_g().norm() << " ";
140  msg << " " << std::setw(10) << std::setprecision(4)
141  << lbfgs.alpha() << " ";
142  msg << " " << std::setw(10) << std::setprecision(4)
143  << lbfgs.alpha0() << " ";
144  msg << " " << std::setw(7)
145  << lbfgs.grad_evals() << " ";
146  msg << " " << lbfgs.note() << " ";
147  logger.info(msg);
148  }
149 
150  if (lbfgs_ss.str().length() > 0) {
151  logger.info(lbfgs_ss);
152  lbfgs_ss.str("");
153  }
154 
155  if (save_iterations) {
156  std::vector<double> values;
157  std::stringstream msg;
158  model.write_array(rng, cont_vector, disc_vector, values,
159  true, true, &msg);
160  if (msg.str().length() > 0)
161  logger.info(msg);
162 
163  values.insert(values.begin(), lp);
164  parameter_writer(values);
165  }
166  }
167 
168  if (!save_iterations) {
169  std::vector<double> values;
170  std::stringstream msg;
171  model.write_array(rng, cont_vector, disc_vector, values,
172  true, true, &msg);
173  if (msg.str().length() > 0)
174  logger.info(msg);
175 
176  values.insert(values.begin(), lp);
177  parameter_writer(values);
178  }
179 
180  int return_code;
181  if (ret >= 0) {
182  logger.info("Optimization terminated normally: ");
183  return_code = error_codes::OK;
184  } else {
185  logger.info("Optimization terminated with error: ");
186  return_code = error_codes::SOFTWARE;
187  }
188  logger.info(" " + lbfgs.get_code_string(ret));
189 
190  return return_code;
191  }
192 
193  }
194  }
195 }
196 #endif
QNUpdateType & get_qnupdate()
Definition: bfgs.hpp:88
Scalar prev_step_size() const
Definition: bfgs.hpp:100
const VectorT & curr_g() const
Definition: bfgs.hpp:93
stan::optimization::BFGSMinimizer< stan::optimization::ModelAdaptor< Model >, stan::optimization::BFGSUpdate_HInv<> > Optimizer
const Scalar & alpha() const
Definition: bfgs.hpp:112
const Scalar & alpha0() const
Definition: bfgs.hpp:111
rosenbrock_model_namespace::rosenbrock_model Model
std::string get_code_string(int retCode)
Definition: bfgs.hpp:117
chain
Check that an output directory exists.
LSOptions< Scalar > _ls_opts
Definition: bfgs.hpp:85
ConvergenceOptions< Scalar > _conv_opts
Definition: bfgs.hpp:86
const std::string & note() const
Definition: bfgs.hpp:115
const size_t iter_num() const
Definition: bfgs.hpp:113
const XML_Char XML_Content * model
Definition: expat.h:151
int lbfgs(Model &model, stan::io::var_context &init, unsigned int random_seed, unsigned int chain, double init_radius, int history_size, double init_alpha, double tol_obj, double tol_rel_obj, double tol_grad, double tol_rel_grad, double tol_param, int num_iterations, bool save_iterations, int refresh, callbacks::interrupt &interrupt, callbacks::logger &logger, callbacks::writer &init_writer, callbacks::writer &parameter_writer)
Definition: lbfgs.hpp:55
boost::ecuyer1988 create_rng(unsigned int seed, unsigned int chain)
Definition: create_rng.hpp:25
void refresh()
Definition: show_event.C:21