MUQ  0.4.3
DistributionWrapper.cpp
Go to the documentation of this file.
1 #include "AllClassWrappers.h"
2 
4 
15 
16 #include <pybind11/pybind11.h>
17 #include <pybind11/stl.h>
18 #include <pybind11/eigen.h>
19 
20 #include <string>
21 
22 #include <functional>
23 #include <vector>
24 
25 using namespace muq::Modeling::PythonBindings;
26 using namespace muq::Modeling;
27 namespace py = pybind11;
28 
29 class PyDistributionTramp : public PyDistribution {
30 public:
31  // Inherit the constructors
33 
34  virtual inline Eigen::VectorXd SampleImpl(std::vector<Eigen::VectorXd> const& inputs) override {
35  PYBIND11_OVERLOAD_PURE(Eigen::VectorXd, PyDistribution, SampleImpl, inputs);
36  }
37 
38  virtual inline double LogDensityImpl(std::vector<Eigen::VectorXd> const& inputs) override {
39  PYBIND11_OVERLOAD_PURE(double, PyDistribution, LogDensityImpl, inputs);
40  }
41 private:
42 };
43 
44 class PyGaussianTramp : public PyGaussianBase {
45 public:
46  /* Inherit the constructors */
47  using PyGaussianBase::PyGaussianBase;
48 
49  /* Trampoline (need one for each virtual function) */
50  unsigned int Dimension() const override {
51  PYBIND11_OVERLOAD(unsigned int,GaussianBase,Dimension);
52  }
53 
54  Eigen::MatrixXd ApplyCovariance(Eigen::Ref<const Eigen::MatrixXd> const& x) const override{
55  PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd, PyGaussianBase, ApplyCovariance, x);
56  }
57 
58  Eigen::MatrixXd ApplyPrecision(Eigen::Ref<const Eigen::MatrixXd> const& x) const override{
59  PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd, PyGaussianBase, ApplyPrecision, x);
60  }
61 
62  virtual Eigen::MatrixXd ApplyCovSqrt(Eigen::Ref<const Eigen::MatrixXd> const& x) const override{
63  PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd, PyGaussianBase, ApplyCovSqrt, x);
64  }
65 
66  Eigen::MatrixXd ApplyPrecSqrt(Eigen::Ref<const Eigen::MatrixXd> const& x) const override{
67  PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd, PyGaussianBase, ApplyPrecSqrt, x);
68  }
69 
70  Eigen::VectorXd SampleImpl(ref_vector<Eigen::VectorXd> const& params) override{
71  PYBIND11_OVERLOAD(Eigen::VectorXd, PyGaussianBase, SampleImpl, params);
72  }
73 
74  void ResetHyperparameters(ref_vector<Eigen::VectorXd> const& params) override{
75  PYBIND11_OVERLOAD(void, PyGaussianBase, ResetHyperparameters, params);
76  }
77 
78  double LogDeterminant() const override{
79  PYBIND11_OVERLOAD(double, PyGaussianBase, LogDeterminant);
80  }
81 
82 };
83 
84 class Publicist : public PyDistribution {
85 public:
86  // Expose protected functions
89 };
90 
92 {
93  py::class_<Distribution, std::shared_ptr<Distribution>> dist(m, "Distribution");
94  dist
95  .def("LogDensity", (double (Distribution::*)(std::vector<Eigen::VectorXd> const&)) &Distribution::LogDensity)
96  .def("LogDensity", (double (Distribution::*)(Eigen::VectorXd const&)) &Distribution::LogDensity)
97  .def("GradLogDensity", (Eigen::VectorXd (Distribution::*)(unsigned int, std::vector<Eigen::VectorXd> const&)) &Distribution::GradLogDensity)
98  .def("ApplyLogDensityHessian", (Eigen::VectorXd (Distribution::*)(unsigned int, unsigned int, std::vector<Eigen::VectorXd> const&, Eigen::VectorXd const&)) &Distribution::ApplyLogDensityHessian)
99  .def("Sample", (Eigen::VectorXd (Distribution::*)(std::vector<Eigen::VectorXd> const&)) &Distribution::Sample)
100  .def("Sample", (Eigen::VectorXd (Distribution::*)()) &Distribution::Sample)
101  .def("AsDensity", &Distribution::AsDensity)
102  .def("AsVariable", &Distribution::AsVariable)
103  .def_readonly("varSize", &Distribution::varSize)
104  .def_readonly("hyperSizes", &Distribution::hyperSizes);
105 
106 
107  py::class_<PyDistribution, PyDistributionTramp, Distribution, std::shared_ptr<PyDistribution>> pydist(m, "PyDistribution");
108  pydist.def(py::init<unsigned int>());
109  pydist.def(py::init<unsigned int, Eigen::VectorXi const&>());
110  pydist.def("SampleImpl", (Eigen::VectorXd (PyDistribution::*)(std::vector<Eigen::VectorXd> const&)) &Publicist::SampleImpl);
111  pydist.def("LogDensityImpl", (double (PyDistribution::*)(std::vector<Eigen::VectorXd> const&)) &Publicist::LogDensityImpl);
112 
113  py::class_<DensityBase, Distribution, ModPiece, std::shared_ptr<DensityBase>> densBase(m, "DensityBase");
114  densBase
115  .def(py::init<Eigen::VectorXi>());
116 
117  py::class_<Density, DensityBase, std::shared_ptr<Density>> dens(m, "Density");
118  dens
119  .def("GetDistribution", &Density::GetDistribution);
120 
121  py::class_<DensityProduct, DensityBase, std::shared_ptr<DensityProduct>> densProd(m,"DensityProduct");
122  densProd
123  .def(py::init<int>());
124 
125  py::class_<RandomVariable, Distribution, ModPiece, std::shared_ptr<RandomVariable>> rv(m, "RandomVariable");
126  rv
127  .def(py::init<std::shared_ptr<Distribution>>());
128 
129  py::class_<UniformBox, Distribution, std::shared_ptr<UniformBox>> uBox(m, "UniformBox");
130  uBox
131  .def(py::init<Eigen::MatrixXd const&>());
132 
133  py::class_<GaussianBase, Distribution, std::shared_ptr<GaussianBase>>(m,"GaussianBase")
134  .def("Dimension", &GaussianBase::Dimension)
135  .def("GetMean", &GaussianBase::GetMean)
136  .def("SetMean", &GaussianBase::SetMean);
137 
138  py::class_<Gaussian, GaussianBase, std::shared_ptr<Gaussian>> gauss(m,"Gaussian");
139  gauss
140  .def(py::init<unsigned int>())
141  .def(py::init<unsigned int, Gaussian::InputMask>())
142  .def(py::init<Eigen::VectorXd const&>())
143  .def(py::init<Eigen::VectorXd const&, Gaussian::InputMask>())
144  .def(py::init<Eigen::VectorXd const&, Eigen::MatrixXd const&>())
145  .def(py::init<Eigen::VectorXd const&, Eigen::MatrixXd const&, Gaussian::Mode>())
146  .def(py::init<Eigen::VectorXd const&, Eigen::MatrixXd const&, Gaussian::Mode, Gaussian::InputMask>())
147  .def("GetMode", &Gaussian::GetMode)
148  .def("GetCovariance", &Gaussian::GetCovariance)
149  .def("GetPrecision", &Gaussian::GetPrecision)
150  .def("ApplyPrecision", &Gaussian::ApplyPrecision)
151  .def("ApplyCovariance", &Gaussian::ApplyCovariance)
152  .def("ApplyCovSqrt", &Gaussian::ApplyCovSqrt)
153  .def("ApplyPrecSqrt", &Gaussian::ApplyPrecSqrt)
154  .def("SetCovariance", &Gaussian::SetCovariance)
155  .def("SetPrecision", &Gaussian::SetPrecision)
156  .def("Condition", &Gaussian::Condition);
157 
158  py::class_<PyGaussianBase, PyGaussianTramp, GaussianBase, Distribution, std::shared_ptr<PyGaussianBase>>(m, "PyGaussianBase")
159  .def(py::init<unsigned int>())
160  .def(py::init<unsigned int, Eigen::VectorXi>())
161  .def(py::init<Eigen::VectorXd>())
162  .def(py::init<Eigen::VectorXd, Eigen::VectorXi>())
163  .def("Dimension", &PyGaussianBase::Dimension)
164  .def("ApplyCovariance", (Eigen::MatrixXd (PyGaussianBase::*)(Eigen::MatrixXd const&) const) &PyGaussianBase::ApplyCovariance)
165  .def("ApplyPrecision", (Eigen::MatrixXd (PyGaussianBase::*)(Eigen::MatrixXd const&) const) &PyGaussianBase::ApplyPrecision)
166  .def("ApplyCovSqrt",(Eigen::MatrixXd (PyGaussianBase::*)(Eigen::MatrixXd const&) const) &PyGaussianBase::ApplyCovSqrt)
167  .def("ApplyPrecSqrt", (Eigen::MatrixXd (PyGaussianBase::*)(Eigen::MatrixXd const&) const) &PyGaussianBase::ApplyPrecSqrt)
168  .def("LogDeterminant", &PyGaussianBase::LogDeterminant)
169  .def("SetMean", &PyGaussianBase::SetMean)
170  .def("GetMean", &PyGaussianBase::GetMean)
171  .def("SampleImpl", (Eigen::VectorXd (PyGaussianBase::*)(ref_vector<Eigen::VectorXd> const&)) &PyGaussianBase::SampleImpl);
172 
173 
174  py::enum_<Gaussian::Mode>(gauss, "Mode")
175  .value("Covariance", Gaussian::Mode::Covariance)
176  .value("Precision", Gaussian::Mode::Precision)
177  .export_values();
178 
179  py::enum_<Gaussian::ExtraInputs>(gauss, "ExtraInputs", py::arithmetic())
180  .value("None", Gaussian::ExtraInputs::None)
181  .value("Mean", Gaussian::ExtraInputs::Mean)
182  .value("DiagCovariance", Gaussian::ExtraInputs::DiagCovariance)
183  .value("DiagPrecision", Gaussian::ExtraInputs::DiagPrecision)
184  .value("FullCovariance", Gaussian::ExtraInputs::FullCovariance)
185  .value("FullPrecision", Gaussian::ExtraInputs::FullPrecision)
186  .export_values();
187 
188 
189  py::class_<InverseGamma, Distribution, std::shared_ptr<InverseGamma>> ig(m, "InverseGamma");
190  ig
191  .def(py::init<double,double>())
192  .def(py::init<Eigen::VectorXd const&,Eigen::VectorXd const&>())
193  .def_readonly("alpha", &InverseGamma::alpha)
194  .def_readonly("beta", &InverseGamma::beta);
195 
196  py::class_<Gamma, Distribution, std::shared_ptr<Gamma>>(m, "Gamma")
197  .def(py::init<double,double>())
198  .def(py::init<Eigen::VectorXd const&,Eigen::VectorXd const&>())
199  .def_readonly("alpha", &Gamma::alpha)
200  .def_readonly("beta", &Gamma::beta)
201  .def_static("FromMoments", &Gamma::FromMoments);
202 
203  py::class_<MixtureDistribution, Distribution, std::shared_ptr<MixtureDistribution>>(m, "MixtureDistribution")
204  .def(py::init<std::vector<std::shared_ptr<Distribution>> const&, Eigen::VectorXd const&>())
205  .def("Components", &MixtureDistribution::Components)
206  .def("Probabilities", &MixtureDistribution::Probabilities);
207 }
virtual std::shared_ptr< Distribution > GetDistribution()
Definition: Density.h:52
std::shared_ptr< Density > AsDensity()
Returns a density built from this distribution.
Definition: Distribution.cpp:9
virtual Eigen::VectorXd ApplyLogDensityHessian(unsigned int const inWrt1, unsigned int const inWrt2, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec)
virtual double LogDensity(ref_vector< Eigen::VectorXd > const &inputs)
Evaluate the log-density.
Eigen::VectorXd Sample()
Sample the distribution with no inputs.
virtual Eigen::VectorXd GradLogDensity(unsigned int wrt, std::vector< Eigen::VectorXd > const &inputs)
Definition: Distribution.h:40
std::shared_ptr< RandomVariable > AsVariable()
Returns a random variable built from this distribution.
const unsigned int varSize
Definition: Distribution.h:145
const Eigen::VectorXi hyperSizes
Definition: Distribution.h:146
const Eigen::VectorXd beta
Definition: Gamma.h:36
const Eigen::VectorXd alpha
Definition: Gamma.h:35
static std::shared_ptr< Gamma > FromMoments(Eigen::VectorXd const &mean, Eigen::VectorXd const &var)
Definition: Gamma.cpp:87
Defines an abstract Gaussian class.@seealso Gaussian.
Definition: GaussianBase.h:15
virtual Eigen::VectorXd const & GetMean() const
Definition: GaussianBase.h:91
virtual unsigned int Dimension() const
virtual double LogDeterminant() const
Definition: GaussianBase.h:105
virtual Eigen::MatrixXd ApplyPrecSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const =0
virtual Eigen::MatrixXd ApplyCovariance(Eigen::Ref< const Eigen::MatrixXd > const &x) const =0
virtual Eigen::MatrixXd ApplyPrecision(Eigen::Ref< const Eigen::MatrixXd > const &x) const =0
virtual void SetMean(Eigen::VectorXd const &newMu)
virtual Eigen::MatrixXd ApplyCovSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const =0
Eigen::MatrixXd GetPrecision() const
Definition: Gaussian.cpp:108
void SetCovariance(Eigen::MatrixXd const &newCov)
Set the covariance matrix.
Definition: Gaussian.cpp:182
void SetPrecision(Eigen::MatrixXd const &newPrec)
Set the precision matrix.
Definition: Gaussian.cpp:199
Mode GetMode() const
Definition: Gaussian.h:61
virtual Eigen::MatrixXd ApplyCovSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
Definition: Gaussian.cpp:215
std::shared_ptr< Gaussian > Condition(Eigen::MatrixXd const &obsMat, Eigen::VectorXd const &data, Eigen::MatrixXd const &obsCov) const
Returns a new Gaussian distribution conditioned on a linear observation.
Definition: Gaussian.cpp:282
virtual Eigen::MatrixXd ApplyPrecSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
Definition: Gaussian.cpp:231
virtual Eigen::MatrixXd ApplyCovariance(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
Definition: Gaussian.cpp:248
virtual Eigen::MatrixXd ApplyPrecision(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
Definition: Gaussian.cpp:264
Eigen::MatrixXd GetCovariance() const
Get the covariance.
Definition: Gaussian.cpp:91
const Eigen::VectorXd beta
Definition: InverseGamma.h:31
const Eigen::VectorXd alpha
Definition: InverseGamma.h:30
std::vector< std::shared_ptr< Distribution > > Components()
virtual Eigen::VectorXd SampleImpl(std::vector< Eigen::VectorXd > const &inputs)=0
virtual double LogDensityImpl(ref_vector< Eigen::VectorXd > const &inputs) override
Implement the log-density.
PyDistribution(unsigned int varSizeIn, Eigen::VectorXi const &hyperSizesIn=Eigen::VectorXi())
virtual Eigen::VectorXd SampleImpl(ref_vector< Eigen::VectorXd > const &inputs) override
Sample the distribution.
void DistributionWrapper(pybind11::module &m)
std::vector< std::reference_wrapper< const T > > ref_vector
A vector of references to something ...
Definition: WorkPiece.h:37