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

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