autoppl  v0.8
A C++ template library for probabilistic programming
ppl::math Namespace Reference

Namespaces

 details
 

Classes

struct  WelfordVar
 

Typedefs

using dist_value_t = util::dist_value_t
 

Functions

size_t padded_length (size_t N)
 
template<class T >
auto autocorrelation (const Eigen::MatrixBase< T > &x)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_pdf (const XType &x, const MeanType &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_pdf (const Eigen::MatrixBase< XType > &x, const MeanType &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MeanType > &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> >>
dist_value_t normal_pdf (const Eigen::MatrixBase< XType > &x, const MeanType &mean, const Eigen::MatrixBase< SigmaType > &sigma)
 
template<class XType , class MeanType , class SigmaType >
dist_value_t normal_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MeanType > &mean, const Eigen::MatrixBase< SigmaType > &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_log_pdf (const XType &x, const MeanType &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_log_pdf (const Eigen::MatrixBase< XType > &x, const MeanType &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<SigmaType> >>
dist_value_t normal_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MeanType > &mean, const SigmaType &sigma)
 
template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> >>
dist_value_t normal_log_pdf (const Eigen::MatrixBase< XType > &x, const MeanType &mean, const Eigen::MatrixBase< SigmaType > &sigma)
 
template<class XType , class MeanType , class SigmaType >
dist_value_t normal_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MeanType > &mean, const Eigen::MatrixBase< SigmaType > &sigma)
 
template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<LocType> && std::is_arithmetic_v<ScaleType> >>
dist_value_t cauchy_log_pdf (const XType &x, const LocType &loc, const ScaleType &scale)
 
template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<LocType> && std::is_arithmetic_v<ScaleType> >>
dist_value_t cauchy_log_pdf (const Eigen::MatrixBase< XType > &x, const LocType &loc, const ScaleType &scale)
 
template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<ScaleType> >>
dist_value_t cauchy_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< LocType > &loc, const ScaleType &scale)
 
template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<LocType> >>
dist_value_t cauchy_log_pdf (const Eigen::MatrixBase< XType > &x, const LocType &loc, const Eigen::MatrixBase< ScaleType > &scale)
 
template<class XType , class LocType , class ScaleType >
dist_value_t cauchy_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< LocType > &loc, const Eigen::MatrixBase< ScaleType > &scale)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_pdf (const XType &x, const MinType &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_pdf (const Eigen::MatrixBase< XType > &x, const MinType &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MinType > &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> >>
dist_value_t uniform_pdf (const Eigen::MatrixBase< XType > &x, const MinType &min, const Eigen::MatrixBase< MaxType > &max)
 
template<class XType , class MinType , class MaxType >
dist_value_t uniform_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MinType > &min, const Eigen::MatrixBase< MaxType > &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_log_pdf (const XType &x, const MinType &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_log_pdf (const Eigen::MatrixBase< XType > &x, const MinType &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MaxType> >>
dist_value_t uniform_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MinType > &min, const MaxType &max)
 
template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> >>
dist_value_t uniform_log_pdf (const Eigen::MatrixBase< XType > &x, const MinType &min, const Eigen::MatrixBase< MaxType > &max)
 
template<class XType , class MinType , class MaxType >
dist_value_t uniform_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< MinType > &min, const Eigen::MatrixBase< MaxType > &max)
 
template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<PType> >>
dist_value_t bernoulli_pdf (const XType &x, const PType &p)
 
template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<PType> >>
dist_value_t bernoulli_pdf (const Eigen::MatrixBase< XType > &x, const PType &p)
 
template<class XType , class PType >
dist_value_t bernoulli_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< PType > &p)
 
template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<PType> >>
dist_value_t bernoulli_log_pdf (const XType &x, const PType &p)
 
template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<PType> >>
dist_value_t bernoulli_log_pdf (const Eigen::MatrixBase< XType > &x, const PType &p)
 
template<class XType , class PType >
dist_value_t bernoulli_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< PType > &p)
 
template<class XType , class VType , class NType >
dist_value_t wishart_log_pdf (const Eigen::MatrixBase< XType > &x, const Eigen::MatrixBase< VType > &v, const NType &n)
 
template<class T >
auto ess (const details::vec_cref_t< T > &samples)
 
template<class T >
auto ess (const Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &samples)
 
template<class T >
lse (T x, T y)
 

Variables

constexpr double SQRT_TWO_PI
 
constexpr double LOG_SQRT_TWO_PI
 
template<class T >
constexpr T inf
 
template<class T >
constexpr T neg_inf
 

Typedef Documentation

◆ dist_value_t

Function Documentation

◆ autocorrelation()

template<class T >
auto ppl::math::autocorrelation ( const Eigen::MatrixBase< T > &  x)
inline

Computes autocorrelation of x where each column of x is a component of a process and hence each row is a time point. More mathematically, x(i,...) is the process value at time i.

For more detail, see: https://lingpipe-blog.com/2012/06/08/autocorrelation-fft-kiss-eigen/ https://github.com/stan-dev/math/blob/41e548e19da5675121c245b535d4019c8bbd754b/stan/math/prim/mat/fun/autocorrelation.hpp

Template Parameters
Tunderlying value type (usually double)
Parameters
xprocess matrix

◆ bernoulli_log_pdf() [1/3]

template<class XType , class PType >
dist_value_t ppl::math::bernoulli_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< PType > &  p 
)
inline

◆ bernoulli_log_pdf() [2/3]

template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<PType> >>
dist_value_t ppl::math::bernoulli_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const PType &  p 
)
inline

◆ bernoulli_log_pdf() [3/3]

template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<PType> >>
dist_value_t ppl::math::bernoulli_log_pdf ( const XType &  x,
const PType &  p 
)
inline

◆ bernoulli_pdf() [1/3]

template<class XType , class PType >
dist_value_t ppl::math::bernoulli_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< PType > &  p 
)
inline

◆ bernoulli_pdf() [2/3]

template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<PType> >>
dist_value_t ppl::math::bernoulli_pdf ( const Eigen::MatrixBase< XType > &  x,
const PType &  p 
)
inline

◆ bernoulli_pdf() [3/3]

template<class XType , class PType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<PType> >>
dist_value_t ppl::math::bernoulli_pdf ( const XType &  x,
const PType &  p 
)
inline

Bernoulli pdf and log pdf (pmf actually). It is defined to clip when p is out of the range [0,1], i.e. if p < 0, then we take p = 0 and if p > 1, then we take p = 1.

◆ cauchy_log_pdf() [1/5]

template<class XType , class LocType , class ScaleType >
dist_value_t ppl::math::cauchy_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< LocType > &  loc,
const Eigen::MatrixBase< ScaleType > &  scale 
)
inline

◆ cauchy_log_pdf() [2/5]

template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<ScaleType> >>
dist_value_t ppl::math::cauchy_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< LocType > &  loc,
const ScaleType &  scale 
)
inline

◆ cauchy_log_pdf() [3/5]

template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<LocType> >>
dist_value_t ppl::math::cauchy_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const LocType &  loc,
const Eigen::MatrixBase< ScaleType > &  scale 
)
inline

◆ cauchy_log_pdf() [4/5]

template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<LocType> && std::is_arithmetic_v<ScaleType> >>
dist_value_t ppl::math::cauchy_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const LocType &  loc,
const ScaleType &  scale 
)
inline

◆ cauchy_log_pdf() [5/5]

template<class XType , class LocType , class ScaleType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<LocType> && std::is_arithmetic_v<ScaleType> >>
dist_value_t ppl::math::cauchy_log_pdf ( const XType &  x,
const LocType &  loc,
const ScaleType &  scale 
)
inline

◆ ess() [1/2]

template<class T >
auto ppl::math::ess ( const details::vec_cref_t< T > &  samples)
inline

Computes the effective sample size (ESS) for a given vector of samples (matrices). Every element of the vector is a matrix of samples for each chain. Every matrix contains the samples as rows, i.e. every row is a sample of an p-dimensional vector, where p is the number of columns of the matrix (number of parameters).

The algorithm assumes that every sample matrix has the same dimensions.

Template Parameters
Tunderlying Eigen expression type
Parameters
samplesvector of samples
Returns
a vector of ESS for each component If number of samples is 1 or less, or there are 0 components, or number of chains is 0, return an empty vector. In either case, the dimension of the return vector is same as the number of components.

◆ ess() [2/2]

template<class T >
auto ppl::math::ess ( const Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &  samples)
inline

Computes the effective sample size (ESS) for a given sample matrix. This is an overload for when there is only chain and can supply a single matrix instead. See above overload for more details.

Template Parameters
Tunderlying Eigen expression type
Parameters
samplessample matrix (1 chain)
Returns
a vector of ESS for each component

◆ lse()

template<class T >
T ppl::math::lse ( x,
y 
)
inline

LogSumExp taken from wikipedia: log(e^x + e^y)

◆ normal_log_pdf() [1/5]

template<class XType , class MeanType , class SigmaType >
dist_value_t ppl::math::normal_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MeanType > &  mean,
const Eigen::MatrixBase< SigmaType > &  sigma 
)
inline

◆ normal_log_pdf() [2/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MeanType > &  mean,
const SigmaType &  sigma 
)
inline

◆ normal_log_pdf() [3/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> >>
dist_value_t ppl::math::normal_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const MeanType &  mean,
const Eigen::MatrixBase< SigmaType > &  sigma 
)
inline

◆ normal_log_pdf() [4/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const MeanType &  mean,
const SigmaType &  sigma 
)
inline

◆ normal_log_pdf() [5/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_log_pdf ( const XType &  x,
const MeanType &  mean,
const SigmaType &  sigma 
)
inline

◆ normal_pdf() [1/5]

template<class XType , class MeanType , class SigmaType >
dist_value_t ppl::math::normal_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MeanType > &  mean,
const Eigen::MatrixBase< SigmaType > &  sigma 
)
inline

◆ normal_pdf() [2/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MeanType > &  mean,
const SigmaType &  sigma 
)
inline

◆ normal_pdf() [3/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> >>
dist_value_t ppl::math::normal_pdf ( const Eigen::MatrixBase< XType > &  x,
const MeanType &  mean,
const Eigen::MatrixBase< SigmaType > &  sigma 
)
inline

◆ normal_pdf() [4/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_pdf ( const Eigen::MatrixBase< XType > &  x,
const MeanType &  mean,
const SigmaType &  sigma 
)
inline

◆ normal_pdf() [5/5]

template<class XType , class MeanType , class SigmaType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MeanType> && std::is_arithmetic_v<SigmaType> >>
dist_value_t ppl::math::normal_pdf ( const XType &  x,
const MeanType &  mean,
const SigmaType &  sigma 
)
inline

◆ padded_length()

size_t ppl::math::padded_length ( size_t  N)
inline

◆ uniform_log_pdf() [1/5]

template<class XType , class MinType , class MaxType >
dist_value_t ppl::math::uniform_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MinType > &  min,
const Eigen::MatrixBase< MaxType > &  max 
)
inline

◆ uniform_log_pdf() [2/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MinType > &  min,
const MaxType &  max 
)
inline

◆ uniform_log_pdf() [3/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> >>
dist_value_t ppl::math::uniform_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const MinType &  min,
const Eigen::MatrixBase< MaxType > &  max 
)
inline

◆ uniform_log_pdf() [4/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const MinType &  min,
const MaxType &  max 
)
inline

◆ uniform_log_pdf() [5/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_log_pdf ( const XType &  x,
const MinType &  min,
const MaxType &  max 
)
inline

◆ uniform_pdf() [1/5]

template<class XType , class MinType , class MaxType >
dist_value_t ppl::math::uniform_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MinType > &  min,
const Eigen::MatrixBase< MaxType > &  max 
)
inline

◆ uniform_pdf() [2/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< MinType > &  min,
const MaxType &  max 
)
inline

◆ uniform_pdf() [3/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> >>
dist_value_t ppl::math::uniform_pdf ( const Eigen::MatrixBase< XType > &  x,
const MinType &  min,
const Eigen::MatrixBase< MaxType > &  max 
)
inline

◆ uniform_pdf() [4/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_pdf ( const Eigen::MatrixBase< XType > &  x,
const MinType &  min,
const MaxType &  max 
)
inline

◆ uniform_pdf() [5/5]

template<class XType , class MinType , class MaxType , class = std::enable_if_t< std::is_arithmetic_v<XType> && std::is_arithmetic_v<MinType> && std::is_arithmetic_v<MaxType> >>
dist_value_t ppl::math::uniform_pdf ( const XType &  x,
const MinType &  min,
const MaxType &  max 
)
inline

◆ wishart_log_pdf()

template<class XType , class VType , class NType >
dist_value_t ppl::math::wishart_log_pdf ( const Eigen::MatrixBase< XType > &  x,
const Eigen::MatrixBase< VType > &  v,
const NType &  n 
)
inline

Variable Documentation

◆ inf

template<class T >
constexpr T ppl::math::inf
inlineconstexpr
Initial value:
=
std::numeric_limits<T>::is_iec559 ?
std::numeric_limits<T>::infinity() :
std::numeric_limits<T>::max()

◆ LOG_SQRT_TWO_PI

constexpr double ppl::math::LOG_SQRT_TWO_PI
inlineconstexpr
Initial value:
=
0.918938533204672741780329736405617

◆ neg_inf

template<class T >
constexpr T ppl::math::neg_inf
inlineconstexpr
Initial value:
=
std::numeric_limits<T>::is_iec559 ?
-std::numeric_limits<T>::infinity() :
std::numeric_limits<T>::lowest()

◆ SQRT_TWO_PI

constexpr double ppl::math::SQRT_TWO_PI
inlineconstexpr
Initial value:
=
2.506628274631000502415765284811045