MUQ  0.4.3
SampleWrapper.cpp
Go to the documentation of this file.
1 #include "AllClassWrappers.h"
2 
9 
11 
13 
14 #include <pybind11/pybind11.h>
15 #include <pybind11/stl.h>
16 #include <pybind11/eigen.h>
17 
18 #include <string>
19 
20 #include <functional>
21 #include <vector>
22 
23 using namespace muq::SamplingAlgorithms;
24 using namespace muq::Utilities;
25 
26 namespace py = pybind11;
27 
28 void PythonBindings::SampleWrapper(py::module &m)
29 {
30  py::class_<SamplingStateIdentity, std::shared_ptr<SamplingStateIdentity>> ssID(m, "SamplingStateIdentity");
31  ssID
32  .def(py::init<int>())
33  .def_readonly("blockInd", &SamplingStateIdentity::blockInd);
34 
35  // py::class_<SamplingStatePartialMoment, std::shared_ptr<SamplingStatePartialMoment>> ssParMom(m, "SamplingStatePartialMoment");
36  // ssParMom
37  // .def(py::init<int, int, Eigen::VectorXd const&>())
38  // .def_readonly("blockInd", &SamplingStatePartialMoment::blockInd)
39  // .def_readonly("momentPower", &SamplingStatePartialMoment::momentPower);
40  // //.def_readonly("mu", &SamplingStatePartialMoment::mu);
41 
42  py::class_<SampleEstimator, std::shared_ptr<SampleEstimator>>(m,"SampleEstimator")
43  .def("CentralMoment", (Eigen::VectorXd (SampleEstimator::*)(unsigned int, int) const) &SampleEstimator::CentralMoment, py::arg("order"), py::arg("blockDim") = -1)
44  .def("CentralMoment", (Eigen::VectorXd (SampleEstimator::*)(unsigned int, Eigen::VectorXd const&, int) const) &SampleEstimator::CentralMoment, py::arg("order"), py::arg("mean"), py::arg("blockDim") = -1)
45  .def("Mean", &SampleEstimator::Mean, py::arg("blockDim") = -1)
46  .def("Variance", (Eigen::VectorXd (SampleEstimator::*)(int) const) &SampleEstimator::Variance, py::arg("blockDim") = -1)
47  .def("Variance", (Eigen::VectorXd (SampleEstimator::*)(Eigen::VectorXd const&, int) const) &SampleEstimator::Variance, py::arg("mean"), py::arg("blockDim") = -1)
48  .def("Covariance", (Eigen::MatrixXd (SampleEstimator::*)(int) const) &SampleEstimator::Covariance, py::arg("blockDim") = -1)
49  .def("Covariance", (Eigen::MatrixXd (SampleEstimator::*)(Eigen::VectorXd const&, int) const) &SampleEstimator::Covariance, py::arg("mean"), py::arg("blockDim") = -1)
50  .def("StandardizedMoment", (Eigen::VectorXd (SampleEstimator::*)(unsigned int, int) const) &SampleEstimator::CentralMoment, py::arg("order"), py::arg("blockDim") = -1)
51  .def("StandardizedMoment", (Eigen::VectorXd (SampleEstimator::*)(unsigned int, Eigen::VectorXd const&, int) const) &SampleEstimator::StandardizedMoment, py::arg("order"), py::arg("mean"), py::arg("blockDim") = -1)
52  .def("StandardizedMoment", (Eigen::VectorXd (SampleEstimator::*)(unsigned int, Eigen::VectorXd const&, Eigen::VectorXd const&, int) const) &SampleEstimator::StandardizedMoment, py::arg("order"), py::arg("mean"),py::arg("stdDev"), py::arg("blockDim") = -1)
53  .def("Skewness", (Eigen::VectorXd (SampleEstimator::*)(int) const) &SampleEstimator::Skewness, py::arg("blockDim") = -1)
54  .def("Skewness", (Eigen::VectorXd (SampleEstimator::*)(Eigen::VectorXd const&, int) const) &SampleEstimator::Skewness, py::arg("mean"), py::arg("blockDim") = -1)
55  .def("Skewness", (Eigen::VectorXd (SampleEstimator::*)(Eigen::VectorXd const&, Eigen::VectorXd const&, int) const) &SampleEstimator::Skewness, py::arg("mean"),py::arg("stdDev"), py::arg("blockDim") = -1)
56  .def("Kurtosis", (Eigen::VectorXd (SampleEstimator::*)(int) const) &SampleEstimator::Kurtosis, py::arg("blockDim") = -1)
57  .def("Kurtosis", (Eigen::VectorXd (SampleEstimator::*)(Eigen::VectorXd const&, int) const) &SampleEstimator::Kurtosis, py::arg("mean"), py::arg("blockDim") = -1)
58  .def("Kurtosis", (Eigen::VectorXd (SampleEstimator::*)(Eigen::VectorXd const&, Eigen::VectorXd const&, int) const) &SampleEstimator::Kurtosis, py::arg("mean"),py::arg("stdDev"), py::arg("blockDim") = -1)
59  .def("ExpectedValue", &SampleEstimator::ExpectedValue, py::arg("f"), py::arg("metasIn")=std::vector<std::string>())
60  .def("ESS", (Eigen::VectorXd (SampleEstimator::*)(int, std::string const&) const) &SampleEstimator::ESS, py::arg("blockDim")=-1, py::arg("method")="Batch")
61  .def("StandardError", (Eigen::VectorXd (SampleEstimator::*)(int, std::string const&) const) &SampleEstimator::StandardError, py::arg("blockDim")=-1, py::arg("method")="Batch");
62 
63 
64  py::class_<MultiIndexEstimator, SampleEstimator, std::shared_ptr<MultiIndexEstimator>>(m, "MultiIndexEstimator")
65  .def(py::init<std::vector<std::shared_ptr<MIMCMCBox>>>());
66 
67  py::class_<SampleCollection, SampleEstimator, std::shared_ptr<SampleCollection>> sampColl(m, "SampleCollection");
68  sampColl
69  .def(py::init<>())
70  .def("__getitem__", (const std::shared_ptr<SamplingState> (SampleCollection::*)(unsigned) const) &SampleCollection::at)
71 // .def("at", &SampleCollection::at)
72  .def("size", &SampleCollection::size)
73  .def("RunningCovariance", (std::vector<Eigen::MatrixXd> (SampleCollection::*)(Eigen::VectorXd const&, int) const) &SampleCollection::RunningCovariance, py::arg("mean"), py::arg("blockDim") = -1)
74  .def("RunningCovariance", (std::vector<Eigen::MatrixXd> (SampleCollection::*)(int) const) &SampleCollection::RunningCovariance, py::arg("blockDim") = -1)
75  .def("Add", &SampleCollection::Add)
76  .def("Weights", &SampleCollection::Weights)
77  .def("AsMatrix", &SampleCollection::AsMatrix, py::arg("blockDim")=-1)
78  .def_static("FromMatrix", &SampleCollection::FromMatrix)
79  .def("GetMeta", (Eigen::MatrixXd (SampleCollection::*)(std::string const&) const) &SampleCollection::GetMeta)
80  .def("ListMeta", &SampleCollection::ListMeta, py::arg("requireAll")=true)
81  .def("WriteToFile", (void (SampleCollection::*)(std::string const&, std::string const&) const) &SampleCollection::WriteToFile, py::arg("filename"), py::arg("dataset") = "/")
82  .def("head", &SampleCollection::head)
83  .def("tail", &SampleCollection::tail)
84  .def("segment", &SampleCollection::segment, py::arg("startInd"),py::arg("length"),py::arg("skipBy")=1)
85  .def("BatchESS", &SampleCollection::BatchESS, py::arg("blockDim")=-1, py::arg("batchSize")=-1, py::arg("overlap")=-1)
86  .def("MultiBatchESS", &SampleCollection::MultiBatchESS, py::arg("blockDim")=-1, py::arg("batchSize")=-1, py::arg("overlap")=-1)
87  .def("BatchError", &SampleCollection::BatchError, py::arg("blockDim")=-1, py::arg("batchSize")=-1, py::arg("overlap")=-1)
88  .def("MultiBatchError", &SampleCollection::MultiBatchError, py::arg("blockDim")=-1, py::arg("batchSize")=-1, py::arg("overlap")=-1);
89 
90  py::class_<MarkovChain, SampleCollection, SampleEstimator, std::shared_ptr<MarkovChain>>(m,"MarkovChain")
91  .def(py::init<>())
92  .def("WolfESS", &MarkovChain::WolffESS, py::arg("blockInd")=-1)
93  .def("WolfError", &MarkovChain::WolffError, py::arg("blockInd")=-1)
94  .def_static("SingleComponentWolffESS", &MarkovChain::SingleComponentWolffESS);
95 
96  m.def_submodule("Diagnostics")
97  .def("Rhat", [](std::vector<std::shared_ptr<SampleCollection>> const& collections){return Diagnostics::Rhat(collections);})
98  .def("Rhat", [](std::vector<std::shared_ptr<SampleCollection>> const& collections, py::dict opts){return Diagnostics::Rhat(collections, ConvertDictToPtree(opts));})
99  .def("Rhat", [](std::vector<std::shared_ptr<MarkovChain>> const& collections){return Diagnostics::Rhat(collections);})
100  .def("Rhat", [](std::vector<std::shared_ptr<MarkovChain>> const& collections, py::dict opts){return Diagnostics::Rhat(collections, ConvertDictToPtree(opts));})
101  .def("Rhat", [](std::vector<std::shared_ptr<MultiIndexEstimator>> const& collections){return Diagnostics::Rhat(collections);})
102  .def("Rhat", [](std::vector<std::shared_ptr<MultiIndexEstimator>> const& collections, py::dict opts){return Diagnostics::Rhat(collections, ConvertDictToPtree(opts));})
103  .def("BasicRhat", &Diagnostics::BasicRhat)
104  .def("BasicMPSRF", &Diagnostics::BasicMPSRF)
105  .def("SplitChains", [](std::vector<std::shared_ptr<SampleCollection>> const& origChains, unsigned int numSegments){return Diagnostics::SplitChains(origChains,numSegments);})
106  .def("TransformChains", &Diagnostics::TransformChains)
107  .def("ComputeRanks", &Diagnostics::ComputeRanks);
108 
109 
110  py::class_<SamplingState, std::shared_ptr<SamplingState>> sampState(m, "SamplingState");
111  sampState
112  .def(py::init<Eigen::VectorXd const&>())
113  .def(py::init<Eigen::VectorXd const&, double>())
114  .def(py::init<std::vector<Eigen::VectorXd> const&>())
115  .def(py::init<std::vector<Eigen::VectorXd> const&, double>())
116  .def_readonly("weight", &SamplingState::weight)
117  .def_readonly("state", &SamplingState::state)
118  .def("HasMeta", &SamplingState::HasMeta)
119  .def("GetMeta", [](std::shared_ptr<SamplingState> self, std::string const& metaKey)
120  -> boost::any& {
121  return self->meta.at(metaKey);
122  })
123  .def("GetMetaSamplingState", [](std::shared_ptr<SamplingState> self, std::string const& metaKey)
124  -> std::shared_ptr<SamplingState> {
125  return muq::Utilities::AnyCast(self->meta.at(metaKey));
126  })
127  .def("TotalDim", &SamplingState::TotalDim)
128  .def("ToVector", &SamplingState::ToVector,py::arg("blockInd")=-1);
129 }
Eigen::VectorXd WolffError(int blockDim) const
Definition: MarkovChain.cpp:40
static double SingleComponentWolffESS(Eigen::Ref< const Eigen::VectorXd > const &trace)
Definition: MarkovChain.cpp:58
Eigen::VectorXd WolffESS(int blockDim) const
Definition: MarkovChain.cpp:47
A class to hold and analyze a collection of SamplingState objects.
virtual std::shared_ptr< SamplingState > at(unsigned i)
std::unordered_map< std::string, Eigen::MatrixXd > GetMeta() const
double MultiBatchESS(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 std::vector< Eigen::MatrixXd > RunningCovariance(int blockInd=-1) const
virtual std::shared_ptr< SampleCollection > tail(unsigned int N) const
virtual std::shared_ptr< SampleCollection > segment(unsigned int startInd, unsigned int length, unsigned int skipBy=1) const
virtual void Add(std::shared_ptr< SamplingState > newSamp)
virtual Eigen::MatrixXd AsMatrix(int blockDim=-1) const
virtual std::shared_ptr< SampleCollection > head(unsigned int N) const
Eigen::VectorXd BatchESS(int blockInd=-1, int batchSize=-1, int overlap=-1) const
virtual Eigen::VectorXd Weights() const
static std::shared_ptr< SampleCollection > FromMatrix(Eigen::MatrixXd const &samps)
virtual Eigen::VectorXd BatchError(int blockInd=-1, int batchSize=-1, int overlap=-1) const
Abstract base class for computing sample-based approximations of expectations.
virtual Eigen::MatrixXd Covariance(int blockInd=-1) const
virtual Eigen::VectorXd Skewness(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
virtual Eigen::VectorXd StandardError(int blockInd, std::string const &method) const =0
virtual Eigen::VectorXd CentralMoment(unsigned int order, int blockNum=-1) const
virtual Eigen::VectorXd StandardizedMoment(unsigned int order, int blockInd=-1) const
virtual Eigen::VectorXd Mean(int blockInd=-1) const
virtual Eigen::VectorXd Kurtosis(int blockInd=-1) const
virtual Eigen::VectorXd Variance(int blockInd=-1) const
virtual Eigen::VectorXd ESS(int blockInd, std::string const &method) const =0
int TotalDim() const
The total number of parameters in the state, i.e., the sum of state[i].size()
double weight
The weight of this state.
Definition: SamplingState.h:51
std::vector< Eigen::VectorXd > state
The state variables.
Definition: SamplingState.h:48
bool HasMeta(std::string const &metaKey)
Checks to see if the meta map contains a particular key.
Eigen::VectorXd ToVector(int blockInd=-1) const
Class for easily casting boost::any's in assignment operations.
Definition: AnyHelpers.h:34
void SampleWrapper(pybind11::module &m)
boost::property_tree::ptree ConvertDictToPtree(pybind11::dict dict)