MUQ  0.4.3
muq::SamplingAlgorithms::MarkovChain Class Reference

A class for storing and working with the results of Markov chain Monte Carlo algorithms. More...

#include <MarkovChain.h>

Inheritance diagram for muq::SamplingAlgorithms::MarkovChain:

Detailed Description

A class for storing and working with the results of Markov chain Monte Carlo algorithms.

The MarkovChain class is a child of SampleCollection where the sample weights correspond to the number of consecutive steps taking the same value, and the weights are unnormalized (i.e., do not sum to one). This is a useful class for storing the chain produced by an MCMC algorithm without storing the duplicate points that result from rejected proposals.

Definition at line 19 of file MarkovChain.h.

Public Member Functions

 MarkovChain ()=default
 
virtual ~MarkovChain ()=default
 
virtual Eigen::VectorXd ESS (std::string const &method="Batch") const override
 
virtual Eigen::VectorXd ESS (int blockDim) const override
 
virtual Eigen::VectorXd ESS (int blockDim, std::string const &method) const override
 
virtual Eigen::VectorXd StandardError (std::string const &method="Batch") const override
 
virtual Eigen::VectorXd StandardError (int blockDim) const override
 
virtual Eigen::VectorXd StandardError (int blockDim, std::string const &method) const override
 
Eigen::VectorXd WolffESS (int blockDim) const
 
Eigen::VectorXd WolffError (int blockDim) const
 
virtual std::shared_ptr< SampleCollectionsegment (unsigned int startInd, unsigned int length, unsigned int skipBy=1) const override
 
- Public Member Functions inherited from muq::SamplingAlgorithms::SampleCollection
 SampleCollection ()=default
 
virtual ~SampleCollection ()=default
 
virtual void Add (std::shared_ptr< SamplingState > newSamp)
 
virtual std::shared_ptr< SamplingStateat (unsigned i)
 
virtual const std::shared_ptr< SamplingStateat (unsigned i) const
 
virtual unsigned int size () const
 
virtual std::shared_ptr< SampleCollectionhead (unsigned int N) const
 
virtual std::shared_ptr< SampleCollectiontail (unsigned int N) const
 
virtual const std::shared_ptr< SamplingStateback () const
 
virtual Eigen::VectorXd CentralMoment (unsigned order, Eigen::VectorXd const &mean, int blockNum=-1) const override
 Computes the componentwise central moments (e.g., variance, skewness, kurtosis, etc..) of a specific order given that we already know the mean. More...
 
virtual Eigen::VectorXd CentralMoment (unsigned int order, int blockNum=-1) const override
 
virtual Eigen::VectorXd Mean (int blockInd=-1) const override
 
virtual Eigen::MatrixXd Covariance (int blockInd=-1) const override
 
virtual Eigen::MatrixXd Covariance (Eigen::VectorXd const &mean, int blockInd=-1) const override
 
virtual std::vector< Eigen::MatrixXd > RunningCovariance (int blockInd=-1) const
 
virtual std::vector< Eigen::MatrixXd > RunningCovariance (Eigen::VectorXd const &mean, int blockInd=-1) const
 
Eigen::VectorXd BatchESS (int blockInd=-1, int batchSize=-1, int overlap=-1) const
 
double MultiBatchESS (int blockInd=-1, int batchSize=-1, int overlap=-1) const
 
virtual Eigen::VectorXd BatchError (int blockInd=-1, int batchSize=-1, int overlap=-1) const
 
virtual Eigen::VectorXd MultiBatchError (int blockInd=-1, int batchSize=-1, int overlap=-1) const
 
virtual Eigen::MatrixXd AsMatrix (int blockDim=-1) const
 
virtual Eigen::VectorXd Weights () const
 
virtual Eigen::VectorXd Rhat (int blockDim, unsigned int numSegments=4, boost::property_tree::ptree options=boost::property_tree::ptree()) const
 
virtual Eigen::VectorXd Rhat (unsigned int numSegments=4, boost::property_tree::ptree options=boost::property_tree::ptree()) const
 
- Public Member Functions inherited from muq::SamplingAlgorithms::SampleEstimator
virtual ~SampleEstimator ()=default
 
virtual unsigned int BlockSize (int blockInd) const =0
 
virtual unsigned int NumBlocks () const =0
 
virtual Eigen::VectorXd CentralMoment (unsigned int order, Eigen::VectorXd const &mean, int blockNum=-1) const
 
virtual Eigen::VectorXd StandardizedMoment (unsigned int order, int blockInd=-1) const
 
virtual Eigen::VectorXd StandardizedMoment (unsigned int order, Eigen::VectorXd const &mean, int blockInd=-1) const
 
virtual Eigen::VectorXd StandardizedMoment (unsigned int order, Eigen::VectorXd const &mean, Eigen::VectorXd const &stdDev, int blockInd=-1) const
 
virtual Eigen::VectorXd Variance (int blockInd=-1) const
 
virtual Eigen::VectorXd Variance (Eigen::VectorXd const &mean, int blockInd=-1) const
 
virtual Eigen::VectorXd Skewness (int blockInd=-1) const
 
virtual Eigen::VectorXd Skewness (Eigen::VectorXd const &mean, int blockInd=-1) const
 
virtual Eigen::VectorXd Skewness (Eigen::VectorXd const &mean, Eigen::VectorXd const &stdDev, int blockInd=-1) const
 
virtual Eigen::VectorXd Kurtosis (int blockInd=-1) const
 
virtual Eigen::VectorXd Kurtosis (Eigen::VectorXd const &mean, int blockInd=-1) const
 
virtual Eigen::VectorXd Kurtosis (Eigen::VectorXd const &mean, Eigen::VectorXd const &stdDev, int blockInd=-1) const
 
virtual Eigen::VectorXd ExpectedValue (std::shared_ptr< muq::Modeling::ModPiece > const &f, std::vector< std::string > const &metains=std::vector< std::string >()) const =0
 

Static Public Member Functions

static double SingleComponentWolffESS (Eigen::Ref< const Eigen::VectorXd > const &trace)
 
- Static Public Member Functions inherited from muq::SamplingAlgorithms::SampleCollection
static std::shared_ptr< SampleCollectionFromMatrix (Eigen::MatrixXd const &samps)
 

Constructor & Destructor Documentation

◆ MarkovChain()

muq::SamplingAlgorithms::MarkovChain::MarkovChain ( )
default

◆ ~MarkovChain()

virtual muq::SamplingAlgorithms::MarkovChain::~MarkovChain ( )
virtualdefault

Member Function Documentation

◆ ESS() [1/3]

virtual Eigen::VectorXd muq::SamplingAlgorithms::MarkovChain::ESS ( int  blockDim) const
inlineoverridevirtual
Parameters
[in]blockDimThe block of the SampleState that we want to use in the ESS computation.
Returns
Either a vector of length \(D\) containing an ESS estimate for each component or a vector of length 1 containing a single multivariate ESS estimate.

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 42 of file MarkovChain.h.

References ESS().

Referenced by ESS().

◆ ESS() [2/3]

Eigen::VectorXd MarkovChain::ESS ( int  blockDim,
std::string const &  method 
) const
overridevirtual
Parameters
[in]blockDimThe block of the SampleState that we want to use in the ESS computation.
[in]methodA string specifying the method used to estimate the ESS. Options in SampleCollection are "Batch" or "MultiBatch". The MarkovChain class also supports a "Wolff" method.
Returns
Either a vector of length \(D\) containing an ESS estimate for each component or a vector of length 1 containing a single multivariate ESS estimate.

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 7 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleCollection::BatchESS(), muq::SamplingAlgorithms::SampleCollection::MultiBatchESS(), and WolffESS().

◆ ESS() [3/3]

virtual Eigen::VectorXd muq::SamplingAlgorithms::MarkovChain::ESS ( std::string const &  method = "Batch") const
inlineoverridevirtual

Computes the effective sample size of the Markov chain.

If method=="Wolff", the spectral method described in "Monte Carlo errors with less error" by Ulli Wolff is employed. This returns an ESS for each component of the chain.

If method=="Batch" (default) The overlapping batch method (OBM) described in [Flegal2010] is used. This method is also applied to each component independently, resulting in an ESS estimate for each component.

If method=="MultiBatch", The multivariate method of [Vats2019] is employed. This method takes into account the joint correlation of all components of the chain and returns a single ESS. This approach is preferred in high dimensional settings.

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 41 of file MarkovChain.h.

References ESS().

Referenced by ESS().

◆ segment()

std::shared_ptr< SampleCollection > MarkovChain::segment ( unsigned int  startInd,
unsigned int  length,
unsigned int  skipBy = 1 
) const
overridevirtual

Returns a new sample collection containing a segment of this sample collection. The Equivalent to a python command like list[startInd:startInd+length:skipBy]

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 158 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleCollection::at(), and muq::SamplingAlgorithms::SampleCollection::size().

◆ SingleComponentWolffESS()

double MarkovChain::SingleComponentWolffESS ( Eigen::Ref< const Eigen::VectorXd > const &  trace)
static

Computes the effective sample size given a vector containing a single component of the Markov chain.

Definition at line 58 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleCollection::size().

Referenced by WolffESS().

◆ StandardError() [1/3]

virtual Eigen::VectorXd muq::SamplingAlgorithms::MarkovChain::StandardError ( int  blockDim) const
inlineoverridevirtual

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 60 of file MarkovChain.h.

References StandardError().

Referenced by StandardError().

◆ StandardError() [2/3]

Eigen::VectorXd MarkovChain::StandardError ( int  blockInd,
std::string const &  method 
) const
overridevirtual

Returns an estimate of the the Monte Carlo standard error (MCSE) \(\hat{\sigma}\) for a Monte Carlo estimate of the mean derived using this SampleEstimator. Recall at the MCSE is the standard deviation of the estimator variance employed in the Central Limit Theorem.

Parameters
[in]blockIndSpecifies the block of the sampling state we're interested in. Defaults to -1, which will result in all blocks of the sampling state being concatenated in the MCSE estimate.
[in]methodA string describing what method should be used to estimate the MCSE. Defaults to "Batch". Options are the same as the ESS function: "Batch", "MultiBatch", and "Wolff". The "Wolff" option is only valid for the "MarkovChain" and "MultiIndexEstimator" classes.
Returns
A vector containing either the MCSE for each component.

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 24 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleCollection::BatchError(), muq::SamplingAlgorithms::SampleCollection::MultiBatchError(), and WolffError().

◆ StandardError() [3/3]

virtual Eigen::VectorXd muq::SamplingAlgorithms::MarkovChain::StandardError ( std::string const &  method = "Batch") const
inlineoverridevirtual

Computes the Monte Carlo standard error (MCSE) of the Markov chain.

If method=="Wolff", the spectral method described in "Monte Carlo errors with less error" by Ulli Wolff [Wolff2004] is employed. This returns an MCSE for each component of the chain.

If method=="Batch" (default) The overlapping batch method (OBM) described in [Flegal2010] is used. This method is also applied to each component independently, resulting in an MCSE estimate for each component.

If method=="MultiBatch", The multivariate method of [Vats2019] is employed. This method takes into account the joint correlation of all components of the chain and returns a single generalized MCSE. This approach is preferred in high dimensional settings.

Reimplemented from muq::SamplingAlgorithms::SampleCollection.

Definition at line 59 of file MarkovChain.h.

References StandardError().

Referenced by StandardError().

◆ WolffError()

Eigen::VectorXd MarkovChain::WolffError ( int  blockDim) const

Computes the MCSE based on the effective sample size returned by WolffESS

Parameters
[in]blockIndSpecifies the block of the sampling state we're interested in. Defaults to -1, which will result in all blocks of the sampling state being concatenated in the MCSE estimate.
Returns
A vector containing the MCSE for each component of the chain.

Definition at line 40 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleEstimator::Variance(), and WolffESS().

Referenced by StandardError().

◆ WolffESS()

Eigen::VectorXd MarkovChain::WolffESS ( int  blockDim) const

Computes the effective sample size using the spectral method of [Wolff2004]

Parameters
[in]blockIndSpecifies the block of the sampling state we're interested in. Defaults to -1, which will result in all blocks of the sampling state being concatenated in the ESS estimate.
Returns
A vector containing the ESS for each component of the chain.

Definition at line 47 of file MarkovChain.cpp.

References muq::SamplingAlgorithms::SampleCollection::AsMatrix(), and SingleComponentWolffESS().

Referenced by ESS(), and WolffError().

Member Data Documentation

◆ meta

std::vector<std::unordered_map<std::string, boost::any> > muq::SamplingAlgorithms::MarkovChain::meta
private

Definition at line 85 of file MarkovChain.h.


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