Stan  2.14.0
probability, sampling & optimization
chains.hpp
Go to the documentation of this file.
1 #ifndef STAN_MCMC_CHAINS_HPP
2 #define STAN_MCMC_CHAINS_HPP
3 
5 #include <stan/math/prim/mat.hpp>
6 #include <boost/accumulators/accumulators.hpp>
7 #include <boost/accumulators/statistics/stats.hpp>
8 #include <boost/accumulators/statistics/mean.hpp>
9 #include <boost/accumulators/statistics/tail_quantile.hpp>
10 #include <boost/accumulators/statistics/p_square_quantile.hpp>
11 #include <boost/accumulators/statistics/variance.hpp>
12 #include <boost/accumulators/statistics/covariance.hpp>
13 #include <boost/accumulators/statistics/variates/covariate.hpp>
14 #include <boost/random/uniform_int_distribution.hpp>
15 #include <boost/random/additive_combine.hpp>
16 #include <algorithm>
17 #include <cmath>
18 #include <iostream>
19 #include <map>
20 #include <stdexcept>
21 #include <string>
22 #include <sstream>
23 #include <utility>
24 #include <vector>
25 #include <cstdlib>
26 
27 namespace stan {
28  namespace mcmc {
29  using Eigen::Dynamic;
30 
45  template <class RNG = boost::random::ecuyer1988>
46  class chains {
47  private:
48  Eigen::Matrix<std::string, Dynamic, 1> param_names_;
49  Eigen::Matrix<Eigen::MatrixXd, Dynamic, 1> samples_;
50  Eigen::VectorXi warmup_;
51 
52  static double mean(const Eigen::VectorXd& x) {
53  return (x.array() / x.size()).sum();
54  }
55 
56  static double variance(const Eigen::VectorXd& x) {
57  double m = mean(x);
58  return ((x.array() - m) / std::sqrt((x.size() - 1.0))).square().sum();
59  }
60 
61  static double sd(const Eigen::VectorXd& x) {
62  return std::sqrt(variance(x));
63  }
64 
65 
66  static double covariance(const Eigen::VectorXd& x,
67  const Eigen::VectorXd& y,
68  std::ostream* err = 0) {
69  if (x.rows() != y.rows() && err)
70  *err << "warning: covariance of different length chains";
71  using boost::accumulators::accumulator_set;
72  using boost::accumulators::stats;
73  using boost::accumulators::tag::variance;
74  using boost::accumulators::tag::covariance;
75  using boost::accumulators::tag::covariate1;
76 
77  accumulator_set<double, stats<covariance<double, covariate1> > > acc;
78 
79  int M = std::min(x.size(), y.size());
80  for (int i = 0; i < M; i++)
81  acc(x(i), boost::accumulators::covariate1 = y(i));
82 
83  return boost::accumulators::covariance(acc) * M / (M-1);
84  }
85 
86  static double correlation(const Eigen::VectorXd& x,
87  const Eigen::VectorXd& y,
88  std::ostream* err = 0) {
89  if (x.rows() != y.rows() && err)
90  *err << "warning: covariance of different length chains";
91  using boost::accumulators::accumulator_set;
92  using boost::accumulators::stats;
93  using boost::accumulators::tag::variance;
94  using boost::accumulators::tag::covariance;
95  using boost::accumulators::tag::covariate1;
96 
97  accumulator_set<double, stats<variance,
98  covariance<double, covariate1> > > acc_xy;
99  accumulator_set<double, stats<variance> > acc_y;
100 
101  int M = std::min(x.size(), y.size());
102  for (int i = 0; i < M; i++) {
103  acc_xy(x(i), boost::accumulators::covariate1 = y(i));
104  acc_y(y(i));
105  }
106 
107  double cov = boost::accumulators::covariance(acc_xy);
108  if (cov > -1e-8 && cov < 1e-8)
109  return cov;
110  return cov / std::sqrt(boost::accumulators::variance(acc_xy)
111  * boost::accumulators::variance(acc_y));
112  }
113 
114  static double quantile(const Eigen::VectorXd& x, const double prob) {
115  using boost::accumulators::accumulator_set;
116  using boost::accumulators::left;
117  using boost::accumulators::quantile;
118  using boost::accumulators::quantile_probability;
119  using boost::accumulators::right;
120  using boost::accumulators::stats;
121  using boost::accumulators::tag::tail;
122  using boost::accumulators::tag::tail_quantile;
123  double M = x.rows();
124  // size_t cache_size = std::min(prob, 1-prob)*M + 2;
125  size_t cache_size = M;
126 
127  if (prob < 0.5) {
128  accumulator_set<double, stats<tail_quantile<left> > >
129  acc(tail<left>::cache_size = cache_size);
130  for (int i = 0; i < M; i++)
131  acc(x(i));
132  return quantile(acc, quantile_probability = prob);
133  }
134  accumulator_set<double, stats<tail_quantile<right> > >
135  acc(tail<right>::cache_size = cache_size);
136  for (int i = 0; i < M; i++)
137  acc(x(i));
138  return quantile(acc, quantile_probability = prob);
139  }
140 
141  static Eigen::VectorXd
142  quantiles(const Eigen::VectorXd& x, const Eigen::VectorXd& probs) {
143  using boost::accumulators::accumulator_set;
144  using boost::accumulators::left;
145  using boost::accumulators::quantile_probability;
146  using boost::accumulators::quantile;
147  using boost::accumulators::right;
148  using boost::accumulators::stats;
149  using boost::accumulators::tag::tail;
150  using boost::accumulators::tag::tail_quantile;
151  double M = x.rows();
152 
153  // size_t cache_size = M/2 + 2;
154  size_t cache_size = M; // 2 + 2;
155 
156  accumulator_set<double, stats<tail_quantile<left> > >
157  acc_left(tail<left>::cache_size = cache_size);
158  accumulator_set<double, stats<tail_quantile<right> > >
159  acc_right(tail<right>::cache_size = cache_size);
160 
161  for (int i = 0; i < M; i++) {
162  acc_left(x(i));
163  acc_right(x(i));
164  }
165 
166  Eigen::VectorXd q(probs.size());
167  for (int i = 0; i < probs.size(); i++) {
168  if (probs(i) < 0.5)
169  q(i) = quantile(acc_left,
170  quantile_probability = probs(i));
171  else
172  q(i) = quantile(acc_right,
173  quantile_probability = probs(i));
174  }
175  return q;
176  }
177 
178  static Eigen::VectorXd autocorrelation(const Eigen::VectorXd& x) {
179  using std::vector;
180  using stan::math::index_type;
181  typedef typename index_type<vector<double> >::type idx_t;
182 
183  std::vector<double> ac;
184  std::vector<double> sample(x.size());
185  for (int i = 0; i < x.size(); i++)
186  sample[i] = x(i);
187  stan::math::autocorrelation(sample, ac);
188 
189  Eigen::VectorXd ac2(ac.size());
190  for (idx_t i = 0; i < ac.size(); i++)
191  ac2(i) = ac[i];
192  return ac2;
193  }
194 
195  static Eigen::VectorXd autocovariance(const Eigen::VectorXd& x) {
196  using std::vector;
197  using stan::math::index_type;
198  typedef typename index_type<vector<double> >::type idx_t;
199 
200  std::vector<double> ac;
201  std::vector<double> sample(x.size());
202  for (int i = 0; i < x.size(); i++)
203  sample[i] = x(i);
204  stan::math::autocovariance(sample, ac);
205 
206  Eigen::VectorXd ac2(ac.size());
207  for (idx_t i = 0; i < ac.size(); i++)
208  ac2(i) = ac[i];
209  return ac2;
210  }
211 
227  double effective_sample_size(const Eigen::Matrix<Eigen::VectorXd,
228  Dynamic, 1> &samples) const {
229  int chains = samples.size();
230 
231  // need to generalize to each jagged samples per chain
232  int n_samples = samples(0).size();
233  for (int chain = 1; chain < chains; chain++) {
234  n_samples = std::min(n_samples,
235  static_cast<int>(samples(chain).size()));
236  }
237 
238  Eigen::Matrix<Eigen::VectorXd, Dynamic, 1> acov(chains);
239  for (int chain = 0; chain < chains; chain++) {
240  acov(chain) = autocovariance(samples(chain));
241  }
242 
243  Eigen::VectorXd chain_mean(chains);
244  Eigen::VectorXd chain_var(chains);
245  for (int chain = 0; chain < chains; chain++) {
246  double n_kept_samples = num_kept_samples(chain);
247  chain_mean(chain) = mean(samples(chain));
248  chain_var(chain) = acov(chain)(0)*n_kept_samples/(n_kept_samples-1);
249  }
250 
251  double mean_var = mean(chain_var);
252  double var_plus = mean_var*(n_samples-1)/n_samples;
253  if (chains > 1)
254  var_plus += variance(chain_mean);
255  Eigen::VectorXd rho_hat_t(n_samples);
256  rho_hat_t.setZero();
257  double rho_hat = 0;
258  int max_t = 0;
259  for (int t = 1; (t < n_samples && rho_hat >= 0); t++) {
260  Eigen::VectorXd acov_t(chains);
261  for (int chain = 0; chain < chains; chain++) {
262  acov_t(chain) = acov(chain)(t);
263  }
264  rho_hat = 1 - (mean_var - mean(acov_t)) / var_plus;
265  if (rho_hat >= 0)
266  rho_hat_t(t) = rho_hat;
267  max_t = t;
268  }
269  double ess = chains * n_samples;
270  if (max_t > 1) {
271  ess /= 1 + 2 * rho_hat_t.sum();
272  }
273  return ess;
274  }
275 
289  double
290  split_potential_scale_reduction(
291  const Eigen::Matrix<Eigen::VectorXd,
292  Dynamic, 1> &samples) const {
293  int chains = samples.size();
294  int n_samples = samples(0).size();
295  for (int chain = 1; chain < chains; chain++) {
296  n_samples = std::min(n_samples,
297  static_cast<int>(samples(chain).size()));
298  }
299  if (n_samples % 2 == 1)
300  n_samples--;
301  int n = n_samples / 2;
302 
303  Eigen::VectorXd split_chain_mean(2*chains);
304  Eigen::VectorXd split_chain_var(2*chains);
305 
306  for (int chain = 0; chain < chains; chain++) {
307  split_chain_mean(2*chain) = mean(samples(chain).topRows(n));
308  split_chain_mean(2*chain+1) = mean(samples(chain).bottomRows(n));
309 
310  split_chain_var(2*chain) = variance(samples(chain).topRows(n));
311  split_chain_var(2*chain+1) = variance(samples(chain).bottomRows(n));
312  }
313 
314  double var_between = n * variance(split_chain_mean);
315  double var_within = mean(split_chain_var);
316 
317  // rewrote [(n-1)*W/n + B/n]/W as (n-1+ B/W)/n
318  return sqrt((var_between/var_within + n-1)/n);
319  }
320 
321  public:
322  explicit chains(const Eigen::Matrix<std::string, Dynamic, 1>& param_names)
323  : param_names_(param_names) { }
324 
325  explicit chains(const std::vector<std::string>& param_names)
326  : param_names_(param_names.size()) {
327  for (size_t i = 0; i < param_names.size(); i++)
328  param_names_(i) = param_names[i];
329  }
330 
331  explicit chains(const stan::io::stan_csv& stan_csv)
332  : param_names_(stan_csv.header) {
333  if (stan_csv.samples.rows() > 0)
334  add(stan_csv);
335  }
336 
337  inline int num_chains() const {
338  return samples_.size();
339  }
340 
341  inline int num_params() const {
342  return param_names_.size();
343  }
344 
345  const Eigen::Matrix<std::string, Dynamic, 1>& param_names() const {
346  return param_names_;
347  }
348 
349  const std::string& param_name(int j) const {
350  return param_names_(j);
351  }
352 
353  int index(const std::string& name) const {
354  int index = -1;
355  for (int i = 0; i < param_names_.size(); i++)
356  if (param_names_(i) == name)
357  return i;
358  return index;
359  }
360 
361  void set_warmup(const int chain, const int warmup) {
362  warmup_(chain) = warmup;
363  }
364 
365  void set_warmup(const int warmup) {
366  warmup_.setConstant(warmup);
367  }
368 
369  const Eigen::VectorXi& warmup() const {
370  return warmup_;
371  }
372 
373  int warmup(const int chain) const {
374  return warmup_(chain);
375  }
376 
377  int num_samples(const int chain) const {
378  return samples_(chain).rows();
379  }
380 
381  int num_samples() const {
382  int n = 0;
383  for (int chain = 0; chain < num_chains(); chain++)
384  n += num_samples(chain);
385  return n;
386  }
387 
388  int num_kept_samples(const int chain) const {
389  return num_samples(chain) - warmup(chain);
390  }
391 
392  int num_kept_samples() const {
393  int n = 0;
394  for (int chain = 0; chain < num_chains(); chain++)
395  n += num_kept_samples(chain);
396  return n;
397  }
398 
399  void add(const int chain,
400  const Eigen::MatrixXd& sample) {
401  if (sample.cols() != num_params())
402  throw std::invalid_argument("add(chain, sample): number of columns"
403  " in sample does not match chains");
404  if (num_chains() == 0 || chain >= num_chains()) {
405  int n = num_chains();
406 
407  // Need this block for Windows. conservativeResize
408  // does not keep the references.
409  Eigen::Matrix<Eigen::MatrixXd, Dynamic, 1>
410  samples_copy(num_chains());
411  Eigen::VectorXi warmup_copy(num_chains());
412  for (int i = 0; i < n; i++) {
413  samples_copy(i) = samples_(i);
414  warmup_copy(i) = warmup_(i);
415  }
416 
417  samples_.resize(chain+1);
418  warmup_.resize(chain+1);
419  for (int i = 0; i < n; i++) {
420  samples_(i) = samples_copy(i);
421  warmup_(i) = warmup_copy(i);
422  }
423  for (int i = n; i < chain+1; i++) {
424  samples_(i) = Eigen::MatrixXd(0, num_params());
425  warmup_(i) = 0;
426  }
427  }
428  int row = samples_(chain).rows();
429  Eigen::MatrixXd new_samples(row+sample.rows(), num_params());
430  new_samples << samples_(chain), sample;
431  samples_(chain) = new_samples;
432  }
433 
434  void add(const Eigen::MatrixXd& sample) {
435  if (sample.rows() == 0)
436  return;
437  if (sample.cols() != num_params())
438  throw std::invalid_argument("add(sample): number of columns in"
439  " sample does not match chains");
440  add(num_chains(), sample);
441  }
442 
450  void add(const std::vector<std::vector<double> >& sample) {
451  int n_row = sample.size();
452  if (n_row == 0)
453  return;
454  int n_col = sample[0].size();
455  Eigen::MatrixXd sample_copy(n_row, n_col);
456  for (int i = 0; i < n_row; i++) {
457  sample_copy.row(i)
458  = Eigen::VectorXd::Map(&sample[i][0], sample[0].size());
459  }
460  add(sample_copy);
461  }
462 
463  void add(const stan::io::stan_csv& stan_csv) {
464  if (stan_csv.header.size() != num_params())
465  throw std::invalid_argument("add(stan_csv): number of columns in"
466  " sample does not match chains");
467  if (!param_names_.cwiseEqual(stan_csv.header).all()) {
468  throw std::invalid_argument("add(stan_csv): header does not match"
469  " chain's header");
470  }
471  add(stan_csv.samples);
472  if (stan_csv.metadata.save_warmup)
473  set_warmup(num_chains()-1, stan_csv.metadata.num_warmup);
474  }
475 
476  Eigen::VectorXd samples(const int chain, const int index) const {
477  return samples_(chain).col(index).bottomRows(num_kept_samples(chain));
478  }
479 
480  Eigen::VectorXd samples(const int index) const {
481  Eigen::VectorXd s(num_kept_samples());
482  int start = 0;
483  for (int chain = 0; chain < num_chains(); chain++) {
484  int n = num_kept_samples(chain);
485  s.middleRows(start, n) = samples_(chain).col(index).bottomRows(n);
486  start += n;
487  }
488  return s;
489  }
490 
491  Eigen::VectorXd samples(const int chain, const std::string& name) const {
492  return samples(chain, index(name));
493  }
494 
495  Eigen::VectorXd samples(const std::string& name) const {
496  return samples(index(name));
497  }
498 
499  double mean(const int chain, const int index) const {
500  return mean(samples(chain, index));
501  }
502 
503  double mean(const int index) const {
504  return mean(samples(index));
505  }
506 
507  double mean(const int chain, const std::string& name) const {
508  return mean(chain, index(name));
509  }
510 
511  double mean(const std::string& name) const {
512  return mean(index(name));
513  }
514 
515  double sd(const int chain, const int index) const {
516  return sd(samples(chain, index));
517  }
518 
519  double sd(const int index) const {
520  return sd(samples(index));
521  }
522 
523  double sd(const int chain, const std::string& name) const {
524  return sd(chain, index(name));
525  }
526 
527  double sd(const std::string& name) const {
528  return sd(index(name));
529  }
530 
531  double variance(const int chain, const int index) const {
532  return variance(samples(chain, index));
533  }
534 
535  double variance(const int index) const {
536  return variance(samples(index));
537  }
538 
539  double variance(const int chain, const std::string& name) const {
540  return variance(chain, index(name));
541  }
542 
543  double variance(const std::string& name) const {
544  return variance(index(name));
545  }
546 
547  double
548  covariance(const int chain, const int index1, const int index2) const {
549  return covariance(samples(chain, index1), samples(chain, index2));
550  }
551 
552  double covariance(const int index1, const int index2) const {
553  return covariance(samples(index1), samples(index2));
554  }
555 
556  double covariance(const int chain, const std::string& name1,
557  const std::string& name2) const {
558  return covariance(chain, index(name1), index(name2));
559  }
560 
561  double
562  covariance(const std::string& name1, const std::string& name2) const {
563  return covariance(index(name1), index(name2));
564  }
565 
566  double
567  correlation(const int chain, const int index1, const int index2) const {
568  return correlation(samples(chain, index1), samples(chain, index2));
569  }
570 
571  double correlation(const int index1, const int index2) const {
572  return correlation(samples(index1), samples(index2));
573  }
574 
575  double correlation(const int chain, const std::string& name1,
576  const std::string& name2) const {
577  return correlation(chain, index(name1), index(name2));
578  }
579 
580  double
581  correlation(const std::string& name1, const std::string& name2) const {
582  return correlation(index(name1), index(name2));
583  }
584 
585  double
586  quantile(const int chain, const int index, const double prob) const {
587  return quantile(samples(chain, index), prob);
588  }
589 
590  double quantile(const int index, const double prob) const {
591  return quantile(samples(index), prob);
592  }
593 
594  double quantile(int chain, const std::string& name, double prob) const {
595  return quantile(chain, index(name), prob);
596  }
597 
598  double quantile(const std::string& name, const double prob) const {
599  return quantile(index(name), prob);
600  }
601 
602  Eigen::VectorXd
603  quantiles(int chain, int index, const Eigen::VectorXd& probs) const {
604  return quantiles(samples(chain, index), probs);
605  }
606 
607  Eigen::VectorXd quantiles(int index, const Eigen::VectorXd& probs) const {
608  return quantiles(samples(index), probs);
609  }
610 
611  Eigen::VectorXd
612  quantiles(int chain, const std::string& name,
613  const Eigen::VectorXd& probs) const {
614  return quantiles(chain, index(name), probs);
615  }
616 
617  Eigen::VectorXd
618  quantiles(const std::string& name, const Eigen::VectorXd& probs) const {
619  return quantiles(index(name), probs);
620  }
621 
622  Eigen::Vector2d central_interval(int chain, int index,
623  double prob) const {
624  double low_prob = (1-prob)/2;
625  double high_prob = 1-low_prob;
626 
627  Eigen::Vector2d interval;
628  interval
629  << quantile(chain, index, low_prob),
630  quantile(chain, index, high_prob);
631  return interval;
632  }
633 
634  Eigen::Vector2d central_interval(int index, double prob) const {
635  double low_prob = (1-prob)/2;
636  double high_prob = 1-low_prob;
637 
638  Eigen::Vector2d interval;
639  interval << quantile(index, low_prob), quantile(index, high_prob);
640  return interval;
641  }
642 
643  Eigen::Vector2d
644  central_interval(int chain, const std::string& name,
645  double prob) const {
646  return central_interval(chain, index(name), prob);
647  }
648 
649  Eigen::Vector2d central_interval(const std::string& name,
650  double prob) const {
651  return central_interval(index(name), prob);
652  }
653 
654  Eigen::VectorXd autocorrelation(const int chain, const int index) const {
655  return autocorrelation(samples(chain, index));
656  }
657 
658  Eigen::VectorXd autocorrelation(int chain,
659  const std::string& name) const {
660  return autocorrelation(chain, index(name));
661  }
662 
663  Eigen::VectorXd autocovariance(const int chain, const int index) const {
664  return autocovariance(samples(chain, index));
665  }
666 
667  Eigen::VectorXd autocovariance(int chain, const std::string& name) const {
668  return autocovariance(chain, index(name));
669  }
670 
671  // FIXME: reimplement using autocorrelation.
672  double effective_sample_size(const int index) const {
673  Eigen::Matrix<Eigen::VectorXd, Dynamic, 1>
674  samples(num_chains());
675  for (int chain = 0; chain < num_chains(); chain++) {
676  samples(chain) = this->samples(chain, index);
677  }
678  return effective_sample_size(samples);
679  }
680 
681  double effective_sample_size(const std::string& name) const {
682  return effective_sample_size(index(name));
683  }
684 
685  double split_potential_scale_reduction(const int index) const {
686  Eigen::Matrix<Eigen::VectorXd, Dynamic, 1>
687  samples(num_chains());
688  for (int chain = 0; chain < num_chains(); chain++) {
689  samples(chain) = this->samples(chain, index);
690  }
691  return split_potential_scale_reduction(samples);
692  }
693 
694  double split_potential_scale_reduction(const std::string& name) const {
695  return split_potential_scale_reduction(index(name));
696  }
697  };
698 
699  }
700 }
701 
702 #endif
double covariance(const int chain, const int index1, const int index2) const
Definition: chains.hpp:548
double variance(const std::string &name) const
Definition: chains.hpp:543
int num_kept_samples() const
Definition: chains.hpp:392
double effective_sample_size(const std::string &name) const
Definition: chains.hpp:681
int index(const std::string &name) const
Definition: chains.hpp:353
double correlation(const int chain, const int index1, const int index2) const
Definition: chains.hpp:567
double correlation(const int chain, const std::string &name1, const std::string &name2) const
Definition: chains.hpp:575
int num_kept_samples(const int chain) const
Definition: chains.hpp:388
void set_warmup(const int warmup)
Definition: chains.hpp:365
const Eigen::VectorXi & warmup() const
Definition: chains.hpp:369
Eigen::MatrixXd samples
int num_chains() const
Definition: chains.hpp:337
double covariance(const int chain, const std::string &name1, const std::string &name2) const
Definition: chains.hpp:556
Eigen::VectorXd samples(const std::string &name) const
Definition: chains.hpp:495
Eigen::VectorXd autocorrelation(const int chain, const int index) const
Definition: chains.hpp:654
double quantile(const int chain, const int index, const double prob) const
Definition: chains.hpp:586
void sample(stan::mcmc::base_mcmc *sampler, int num_warmup, int num_samples, int num_thin, int refresh, bool save, stan::services::sample::mcmc_writer< Model, SampleRecorder, DiagnosticRecorder, MessageRecorder > &mcmc_writer, stan::mcmc::sample &init_s, Model &model, RNG &base_rng, const std::string &prefix, const std::string &suffix, std::ostream &o, StartTransitionCallback &callback, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
Definition: sample.hpp:17
Probability, optimization and sampling library.
Eigen::VectorXd samples(const int chain, const std::string &name) const
Definition: chains.hpp:491
int num_samples(const int chain) const
Definition: chains.hpp:377
Eigen::VectorXd quantiles(int chain, int index, const Eigen::VectorXd &probs) const
Definition: chains.hpp:603
chains(const stan::io::stan_csv &stan_csv)
Definition: chains.hpp:331
Eigen::Vector2d central_interval(const std::string &name, double prob) const
Definition: chains.hpp:649
void add(const stan::io::stan_csv &stan_csv)
Definition: chains.hpp:463
Eigen::VectorXd quantiles(int chain, const std::string &name, const Eigen::VectorXd &probs) const
Definition: chains.hpp:612
double mean(const int index) const
Definition: chains.hpp:503
double sd(const int index) const
Definition: chains.hpp:519
double correlation(const std::string &name1, const std::string &name2) const
Definition: chains.hpp:581
Eigen::Vector2d central_interval(int chain, const std::string &name, double prob) const
Definition: chains.hpp:644
Eigen::VectorXd autocorrelation(int chain, const std::string &name) const
Definition: chains.hpp:658
double covariance(const std::string &name1, const std::string &name2) const
Definition: chains.hpp:562
double split_potential_scale_reduction(const int index) const
Definition: chains.hpp:685
Eigen::VectorXd samples(const int chain, const int index) const
Definition: chains.hpp:476
Eigen::VectorXd samples(const int index) const
Definition: chains.hpp:480
chains(const Eigen::Matrix< std::string, Dynamic, 1 > &param_names)
Definition: chains.hpp:322
double sd(const std::string &name) const
Definition: chains.hpp:527
const Eigen::Matrix< std::string, Dynamic, 1 > & param_names() const
Definition: chains.hpp:345
double correlation(const int index1, const int index2) const
Definition: chains.hpp:571
double variance(const int index) const
Definition: chains.hpp:535
void add(const std::vector< std::vector< double > > &sample)
Convert a vector of vector<double> to Eigen::MatrixXd.
Definition: chains.hpp:450
const std::string & param_name(int j) const
Definition: chains.hpp:349
Eigen::VectorXd quantiles(const std::string &name, const Eigen::VectorXd &probs) const
Definition: chains.hpp:618
double sd(const int chain, const int index) const
Definition: chains.hpp:515
double quantile(int chain, const std::string &name, double prob) const
Definition: chains.hpp:594
chains(const std::vector< std::string > &param_names)
Definition: chains.hpp:325
int warmup(const int chain) const
Definition: chains.hpp:373
int num_samples() const
Definition: chains.hpp:381
int num_params() const
Definition: chains.hpp:341
Eigen::VectorXd autocovariance(int chain, const std::string &name) const
Definition: chains.hpp:667
double variance(const int chain, const int index) const
Definition: chains.hpp:531
double effective_sample_size(const int index) const
Definition: chains.hpp:672
void add(const Eigen::MatrixXd &sample)
Definition: chains.hpp:434
An mcmc::chains object stores parameter names and dimensionalities along with samples from multiple c...
Definition: chains.hpp:46
double covariance(const int index1, const int index2) const
Definition: chains.hpp:552
double quantile(const int index, const double prob) const
Definition: chains.hpp:590
double split_potential_scale_reduction(const std::string &name) const
Definition: chains.hpp:694
void add(const int chain, const Eigen::MatrixXd &sample)
Definition: chains.hpp:399
double mean(const std::string &name) const
Definition: chains.hpp:511
Eigen::VectorXd autocovariance(const int chain, const int index) const
Definition: chains.hpp:663
Eigen::Vector2d central_interval(int chain, int index, double prob) const
Definition: chains.hpp:622
double mean(const int chain, const std::string &name) const
Definition: chains.hpp:507
Eigen::Vector2d central_interval(int index, double prob) const
Definition: chains.hpp:634
double variance(const int chain, const std::string &name) const
Definition: chains.hpp:539
Eigen::VectorXd quantiles(int index, const Eigen::VectorXd &probs) const
Definition: chains.hpp:607
double sd(const int chain, const std::string &name) const
Definition: chains.hpp:523
stan_csv_metadata metadata
double mean(const int chain, const int index) const
Definition: chains.hpp:499
Eigen::Matrix< std::string, Eigen::Dynamic, 1 > header
double quantile(const std::string &name, const double prob) const
Definition: chains.hpp:598
void set_warmup(const int chain, const int warmup)
Definition: chains.hpp:361

     [ Stan Home Page ] © 2011–2016, Stan Development Team.