1 #include "AllClassWrappers.h"
16 #include <pybind11/pybind11.h>
17 #include <pybind11/stl.h>
18 #include <pybind11/eigen.h>
34 virtual inline Eigen::VectorXd SampleImpl(std::vector<Eigen::VectorXd>
const& inputs)
override {
35 PYBIND11_OVERLOAD_PURE(Eigen::VectorXd,
PyDistribution, SampleImpl, inputs);
38 virtual inline double LogDensityImpl(std::vector<Eigen::VectorXd>
const& inputs)
override {
39 PYBIND11_OVERLOAD_PURE(
double,
PyDistribution, LogDensityImpl, inputs);
47 using PyGaussianBase::PyGaussianBase;
50 unsigned int Dimension()
const override {
54 Eigen::MatrixXd ApplyCovariance(Eigen::Ref<const Eigen::MatrixXd>
const& x)
const override{
55 PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd,
PyGaussianBase, ApplyCovariance, x);
58 Eigen::MatrixXd ApplyPrecision(Eigen::Ref<const Eigen::MatrixXd>
const& x)
const override{
59 PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd,
PyGaussianBase, ApplyPrecision, x);
62 virtual Eigen::MatrixXd ApplyCovSqrt(Eigen::Ref<const Eigen::MatrixXd>
const& x)
const override{
63 PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd,
PyGaussianBase, ApplyCovSqrt, x);
66 Eigen::MatrixXd ApplyPrecSqrt(Eigen::Ref<const Eigen::MatrixXd>
const& x)
const override{
67 PYBIND11_OVERLOAD_PURE(Eigen::MatrixXd,
PyGaussianBase, ApplyPrecSqrt, x);
71 PYBIND11_OVERLOAD(Eigen::VectorXd,
PyGaussianBase, SampleImpl, params);
75 PYBIND11_OVERLOAD(
void,
PyGaussianBase, ResetHyperparameters, params);
78 double LogDeterminant()
const override{
93 py::class_<Distribution, std::shared_ptr<Distribution>> dist(m,
"Distribution");
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);
113 py::class_<DensityBase, Distribution, ModPiece, std::shared_ptr<DensityBase>> densBase(m,
"DensityBase");
115 .def(py::init<Eigen::VectorXi>());
117 py::class_<Density, DensityBase, std::shared_ptr<Density>> dens(m,
"Density");
121 py::class_<DensityProduct, DensityBase, std::shared_ptr<DensityProduct>> densProd(m,
"DensityProduct");
123 .def(py::init<int>());
125 py::class_<RandomVariable, Distribution, ModPiece, std::shared_ptr<RandomVariable>> rv(m,
"RandomVariable");
127 .def(py::init<std::shared_ptr<Distribution>>());
129 py::class_<UniformBox, Distribution, std::shared_ptr<UniformBox>> uBox(m,
"UniformBox");
131 .def(py::init<Eigen::MatrixXd const&>());
133 py::class_<GaussianBase, Distribution, std::shared_ptr<GaussianBase>>(m,
"GaussianBase")
138 py::class_<Gaussian, GaussianBase, std::shared_ptr<Gaussian>> gauss(m,
"Gaussian");
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>())
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>())
174 py::enum_<Gaussian::Mode>(gauss,
"Mode")
175 .value(
"Covariance", Gaussian::Mode::Covariance)
176 .value(
"Precision", Gaussian::Mode::Precision)
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)
189 py::class_<InverseGamma, Distribution, std::shared_ptr<InverseGamma>> ig(m,
"InverseGamma");
191 .def(py::init<double,double>())
192 .def(py::init<Eigen::VectorXd const&,Eigen::VectorXd const&>())
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&>())
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&>())
virtual std::shared_ptr< Distribution > GetDistribution()
std::shared_ptr< Density > AsDensity()
Returns a density built from this distribution.
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)
std::shared_ptr< RandomVariable > AsVariable()
Returns a random variable built from this distribution.
const unsigned int varSize
const Eigen::VectorXi hyperSizes
const Eigen::VectorXd beta
const Eigen::VectorXd alpha
static std::shared_ptr< Gamma > FromMoments(Eigen::VectorXd const &mean, Eigen::VectorXd const &var)
Defines an abstract Gaussian class.@seealso Gaussian.
virtual Eigen::VectorXd const & GetMean() const
virtual unsigned int Dimension() const
virtual double LogDeterminant() const
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
void SetCovariance(Eigen::MatrixXd const &newCov)
Set the covariance matrix.
void SetPrecision(Eigen::MatrixXd const &newPrec)
Set the precision matrix.
virtual Eigen::MatrixXd ApplyCovSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
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.
virtual Eigen::MatrixXd ApplyPrecSqrt(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
virtual Eigen::MatrixXd ApplyCovariance(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
virtual Eigen::MatrixXd ApplyPrecision(Eigen::Ref< const Eigen::MatrixXd > const &x) const override
Eigen::MatrixXd GetCovariance() const
Get the covariance.
const Eigen::VectorXd beta
const Eigen::VectorXd alpha
std::vector< std::shared_ptr< Distribution > > Components()
Eigen::VectorXd Probabilities()
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 ...