26 #ifndef O2SCL_PROB_DENS_FUNC_H 27 #define O2SCL_PROB_DENS_FUNC_H 29 #include <gsl/gsl_rng.h> 30 #include <gsl/gsl_randist.h> 31 #include <gsl/gsl_cdf.h> 34 #include <gsl/gsl_poly.h> 38 #include <boost/numeric/ublas/vector.hpp> 39 #include <boost/numeric/ublas/operation.hpp> 41 #include <o2scl/hist.h> 42 #include <o2scl/rng_gsl.h> 43 #include <o2scl/search_vec.h> 44 #include <o2scl/cholesky.h> 46 #include <o2scl/vec_stats.h> 48 #ifndef DOXYGEN_NO_O2NS 76 virtual double pdf(
double x)
const {
88 virtual double cdf(
double x)
const {
149 O2SCL_ERR2(
"Tried to create a Gaussian dist. with sigma",
150 "<0 in prob_dens_gaussian::prob_dens_gaussian().",
194 "in prob_dens_gaussian::prob_dens_gaussian().",
209 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
218 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
221 return cent_+gsl_ran_gaussian(&r,sigma_);
225 virtual double pdf(
double x)
const {
227 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
230 return gsl_ran_gaussian_pdf(x-cent_,sigma_);
236 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
239 return log(gsl_ran_gaussian_pdf(x-cent_,sigma_));
243 virtual double cdf(
double x)
const {
245 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
248 return gsl_cdf_gaussian_P(x-cent_,sigma_);
254 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
257 if (in_cdf<0.0 || in_cdf>1.0) {
258 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
259 "prob_dens_gaussian::invert_cdf().",
exc_einval);
261 return gsl_cdf_gaussian_Pinv(in_cdf,sigma_)+cent_;
267 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
285 virtual double lower_limit()
const=0;
288 virtual double upper_limit()
const=0;
348 if (
this==&pdg)
return *
this;
377 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
386 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
395 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
398 return gsl_ran_flat(&r,ll,ul);
402 virtual double pdf(
double x)
const {
404 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
407 if (x<ll || x>ul)
return 0.0;
408 return gsl_ran_flat_pdf(x,ll,ul);
414 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
417 if (x<ll || x>ul)
return 0.0;
418 return log(gsl_ran_flat_pdf(x,ll,ul));
422 virtual double cdf(
double x)
const {
424 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
427 if (x<ll)
return 0.0;
428 if (x>ul)
return 1.0;
429 return gsl_cdf_flat_P(x,ll,ul);
435 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
438 if (in_cdf<0.0 || in_cdf>1.0) {
439 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
440 "prob_dens_uniform::invert_cdf().",
exc_einval);
442 return gsl_cdf_flat_Pinv(in_cdf,ll,ul);
448 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
510 O2SCL_ERR2(
"Tried to create log normal dist. with mu or sigma",
511 "<0 in prob_dens_lognormal::prob_dens_lognormal().",
530 if (
this==&pdg)
return *
this;
540 O2SCL_ERR2(
"Tried to set mu or sigma negative",
541 "in prob_dens_lognormal::prob_dens_lognormal().",
555 return gsl_ran_lognormal(&r,mu_,sigma_);
559 virtual double pdf(
double x)
const {
563 return gsl_ran_lognormal_pdf(x,mu_,sigma_);
571 return log(gsl_ran_lognormal_pdf(x,mu_,sigma_));
575 virtual double cdf(
double x)
const {
579 return gsl_cdf_lognormal_P(x,mu_,sigma_);
584 if (in_cdf<0.0 || in_cdf>1.0) {
585 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
586 "prob_dens_lognormal::invert_cdf().",
exc_einval);
588 return gsl_cdf_lognormal_Pinv(in_cdf,mu_,sigma_);
594 O2SCL_ERR2(
"Parameters not set in prob_dens_lognormal::",
654 virtual double lower_limit()
const;
657 virtual double upper_limit()
const;
660 virtual double pdf(
double x)
const;
663 virtual double log_pdf(
double x)
const;
666 virtual double cdf(
double x)
const;
683 template<
class vec_t=boost::numeric::ublas::vector<
double> >
689 virtual size_t dim()
const {
695 virtual double pdf(
const vec_t &x)
const {
701 virtual double log_pdf(
const vec_t &x)
const {
703 if (!std::isfinite(val) || val<0.0) {
707 double val2=log(
pdf(x));
708 if (!std::isfinite(val2)) {
709 std::cout << val <<
" " << val2 << std::endl;
727 template<
class vec_t=boost::numeric::ublas::vector<
double> >
733 std::vector<prob_dens_func>
list;
757 virtual size_t dim()
const {
762 virtual double pdf(
const vec_t &x)
const {
764 for(
size_t i=0;i<list.size();i++) ret*=list[i].
pdf(x[i]);
769 virtual double log_pdf(
const vec_t &x)
const {
771 for(
size_t i=0;i<list.size();i++) ret*=list[i].
pdf(x[i]);
777 for(
size_t i=0;i<list.size();i++) x[i]=list[i]();
835 template<
class vec_t=boost::numeric::ublas::vector<
double> >
888 void set(
double x_cent,
double y_cent,
double x_std,
double y_std,
890 if (fabs(covar)>=1.0) {
891 O2SCL_ERR2(
"Covariance cannot have magnitude equal or larger than ",
892 "1 in prob_dens_mdim_biv_gaussian::set().",
905 virtual double pdf(
const vec_t &v)
const {
906 double x=v[0], y=v[1];
907 double arg=-((x-x0)*(x-x0)/sig_x/sig_x+
908 (y-y0)*(y-y0)/sig_y/sig_y-
909 2.0*rho*(x-x0)*(y-y0)/sig_x/sig_y)/2.0/(1.0-rho*rho);
922 double arg=-2.0*log(1.0-integral);
926 sig_y/sqrt(1.0-rho*rho);
932 virtual void contour(
double level,
double theta, vec_t &x) {
934 O2SCL_ERR2(
"Cannot produce contours for negative values in ",
935 "prob_dens_mdim_biv_gaussian::contour().",
940 O2SCL_ERR2(
"Cannot produce contours larger than maximum in ",
941 "prob_dens_mdim_biv_gaussian::contour().",
945 sqrt(1.0-rho*rho))*2.0*(1.0-rho*rho);
946 double r2=arg/(cos(theta)*cos(theta)/sig_x/sig_x+
947 sin(theta)*sin(theta)/sig_y/sig_y-
948 2.0*rho/sig_x/sig_y*cos(theta)*sin(theta));
949 x[0]=sqrt(r2)*cos(theta)+x0;
950 x[1]=sqrt(r2)*sin(theta)+y0;
985 template<
class vec_t=boost::numeric::ublas::vector<
double>,
986 class mat_t=boost::numeric::ublas::matrix<
double> >
1083 template<
class mat2_t,
class vec2_t,
1086 const vec2_t &vals) {
1088 vec_t peak2(p_mdim);
1089 mat_t covar(p_mdim,p_mdim);
1093 for(
size_t i=0;i<p_mdim;i++) {
1094 const mat2_col_t col(pts,i);
1095 peak2[i]=o2scl::wvector_mean<mat2_col_t>(n_pts,col,vals);
1097 covar(i,i)=o2scl::wvector_stddev<mat2_col_t>(n_pts,col,vals);
1098 covar(i,i)*=covar(i,i);
1101 for(
size_t i=0;i<p_mdim;i++) {
1102 mat2_col_t col_i(pts,i);
1103 for(
size_t j=i+1;j<p_mdim;j++) {
1104 const mat2_col_t col_j(pts,j);
1110 set(p_mdim,peak2,covar);
1116 set(p_ndim,p_peak,covar);
1124 void set(
size_t p_ndim, vec_t &p_peak, mat_t &covar) {
1126 O2SCL_ERR(
"Zero dimension in prob_dens_mdim_gaussian::set().",
1132 for(
size_t i=0;i<ndim;i++) peak[i]=p_peak[i];
1142 o2scl_linalg::cholesky_invert<mat_t>(ndim,covar_inv);
1146 double sqrt_det=1.0;
1147 for(
size_t i=0;i<ndim;i++) {
1148 if (!std::isfinite(chol(i,i))) {
1149 O2SCL_ERR2(
"An entry of the Cholesky decomposition was not finite ",
1152 sqrt_det*=chol(i,i);
1153 for(
size_t j=0;j<ndim;j++) {
1154 if (i<j) chol(i,j)=0.0;
1160 if (!std::isfinite(norm)) {
1169 void set_alt(
size_t p_ndim, vec_t &p_peak, mat_t &p_chol,
1170 mat_t &p_covar_inv,
double p_norm) {
1174 covar_inv=p_covar_inv;
1196 template<
class vec_vec_t,
class mat_col_t,
class func_t>
1198 vec_vec_t &x, vec_t &y, func_t &fcovar) {
1202 mat_t KXsX(n_dim,n_init);
1203 for(
size_t irow=n_init;irow<n_dim+n_init;irow++) {
1204 for(
size_t icol=0;icol<n_init;icol++) {
1205 KXsX(irow-n_init,icol)=fcovar(x[irow],x[icol]);
1209 mat_t KXXs=boost::numeric::ublas::trans(KXsX);
1211 mat_t KXX(n_init,n_init);
1212 for(
size_t irow=0;irow<n_init;irow++) {
1213 for(
size_t icol=0;icol<n_init;icol++) {
1215 KXX(irow,icol)=KXX(icol,irow);
1217 KXX(irow,icol)=fcovar(x[irow],x[icol]);
1222 mat_t KXsXs(n_dim,n_dim);
1223 for(
size_t irow=n_init;irow<n_dim+n_init;irow++) {
1224 for(
size_t icol=n_init;icol<n_dim+n_init;icol++) {
1226 KXsXs(irow-n_init,icol-n_init)=KXsXs(icol-n_init,irow-n_init);
1228 KXsXs(irow-n_init,icol-n_init)=fcovar(x[irow],x[icol]);
1234 mat_t inv_KXX(n_init,n_init);
1240 "prob_dens_mdim_gaussian::set_gproc().",
1243 o2scl_linalg::LU_invert<mat_t,mat_t,mat_col_t>(n_init,KXX,p,inv_KXX);
1246 vec_t prod(n_init), mean(n_dim);
1247 boost::numeric::ublas::axpy_prod(inv_KXX,y,prod,
true);
1248 boost::numeric::ublas::axpy_prod(KXsX,prod,mean,
true);
1251 mat_t covar(n_dim,n_dim), prod2(n_init,n_dim), prod3(n_dim,n_dim);
1252 boost::numeric::ublas::axpy_prod(inv_KXX,KXXs,prod2,
true);
1253 boost::numeric::ublas::axpy_prod(KXsX,prod2,prod3,
true);
1257 this->
set(n_dim,mean,covar);
1262 virtual double pdf(
const vec_t &x)
const {
1264 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1268 for(
size_t i=0;i<ndim;i++) {
1271 vtmp=prod(covar_inv,q);
1272 ret*=exp(-0.5*inner_prod(q,vtmp));
1279 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1282 double ret=log(norm);
1283 for(
size_t i=0;i<ndim;i++) q[i]=x[i]-peak[i];
1284 vtmp=prod(covar_inv,q);
1285 ret+=-0.5*inner_prod(q,vtmp);
1292 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1295 for(
size_t i=0;i<ndim;i++) q[i]=pdg();
1297 for(
size_t i=0;i<ndim;i++) x[i]=peak[i]+vtmp[i];
1312 template<
class vec_t=boost::numeric::ublas::vector<
double>,
1313 class mat_t=boost::numeric::ublas::matrix<
double> >
1343 vec_t &p_low, vec_t &p_high) {
1344 set(p_ndim,p_peak,covar,p_low,p_high);
1354 void set(
size_t p_ndim, vec_t &p_peak, mat_t &covar,
1355 vec_t &p_low, vec_t &p_high) {
1365 virtual double pdf(
const vec_t &x)
const {
1366 for(
size_t i=0;i<this->ndim;i++) {
1368 O2SCL_ERR(
"Parameter too small in pdf().",
1372 O2SCL_ERR(
"Parameter too large in pdf().",
1383 for(
size_t i=0;i<this->ndim;i++) {
1385 O2SCL_ERR(
"Parameter too small in log_pdf().",
1389 O2SCL_ERR(
"Parameter too large in log_pdf().",
1401 while (done==
false) {
1405 for(
size_t i=0;i<this->ndim;i++) {
1411 }
else if (x[i]>high[i]) {
1420 "prob_dens_mdim_bound_gaussian::operator().",
1444 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1457 virtual double pdf(
const vec_t &x_B,
const vec_t &x_A)
const=0;
1462 virtual double log_pdf(
const vec_t &x_B,
const vec_t &x_A)
const=0;
1465 virtual void operator()(
const vec_t &x_B, vec_t &x_A)
const=0;
1506 double val1=
log_pdf(x_prime,x);
1507 double val2=
log_pdf(x,x_prime);
1508 if (!std::isfinite(val1)) {
1509 std::cout <<
"val1: " << val1 << std::endl;
1510 O2SCL_ERR(
"Log pdf not finite in log_metrop_hast 1.",
1513 if (!std::isfinite(val2)) {
1514 std::cout <<
"val2: " << val2 << std::endl;
1515 O2SCL_ERR(
"Log pdf not finite in log_metrop_hast 2.",
1565 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1597 for(
size_t i=0;i<sz;i++) {
1600 if (!std::isfinite(low[i]) || !std::isfinite(high[i])) {
1601 O2SCL_ERR2(
"Limit not finite in prob_cond_mdim_fixed_step::",
1606 if (low[i]>high[i]) {
1652 (vec_t &step, vec_t &low, vec_t &high) {
1653 if (step.size()!=low.size()) {
1654 O2SCL_ERR2(
"Vectors 'step' and 'low' mismatched in ",
1655 "prob_cond_mdim_fixed_step constructor.",
1658 if (step.size()!=high.size()) {
1659 O2SCL_ERR2(
"Vectors 'step' and 'high' mismatched in ",
1660 "prob_cond_mdim_fixed_step constructor.",
1663 set_internal(step.size(),step,low,high);
1668 virtual int set(vec_t &step, vec_t &low, vec_t &high) {
1669 if (step.size()!=low.size()) {
1670 O2SCL_ERR2(
"Vectors 'step' and 'low' mismatched in ",
1671 "prob_cond_mdim_fixed_step::set().",
1674 if (step.size()!=high.size()) {
1675 O2SCL_ERR2(
"Vectors 'step' and 'high' mismatched in ",
1676 "prob_cond_mdim_fixed_step::set().",
1679 set_internal(step.size(),step,low,high);
1685 return u_step.size();
1691 virtual double pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1693 for(
size_t i=0;i<u_step.size();i++) {
1695 if (fabs(x_A[i]-x_B[i])>u_step[i])
return 0.0;
1697 if (x_B[i]-u_step[i]<u_low[i]) {
1698 if (x_B[i]+u_step[i]>u_high[i]) {
1700 vol1*=u_high[i]-u_low[i];
1703 vol1*=x_B[i]+u_step[i]-u_low[i];
1706 if (x_B[i]+u_step[i]>u_high[i]) {
1708 vol1*=u_high[i]-(x_B[i]-u_step[i]);
1712 vol1*=2.0*u_step[i];
1722 virtual double log_pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1723 return log(
pdf(x_B,x_A));
1728 size_t nv=u_step.size();
1729 for(
size_t i=0;i<nv;i++) {
1730 if (x_B[i]<u_low[i] || x_B[i]>u_high[i]) {
1731 std::cout <<
"Input out of bounds in fixed_step::operator(): " 1732 << i <<
" " << x_B[i] <<
" " 1733 << u_low[i] <<
" " << u_high[i] << std::endl;
1734 O2SCL_ERR(
"Input out of bounds in fixed_step::operator().",
1738 x_A[i]=x_B[i]+u_step[i]*(rg.
random()*2.0-1.0);
1739 }
while (x_A[i]<u_low[i] || x_A[i]>u_high[i]);
1757 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1794 virtual double pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1795 return base.
pdf(x_A);
1801 virtual double log_pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1823 template<
class vec_t=boost::numeric::ublas::vector<
double>,
1824 class mat_t=boost::numeric::ublas::matrix<
double> >
1901 void set(
size_t p_ndim, mat_t &covar) {
1903 O2SCL_ERR(
"Zero dimension in prob_cond_mdim_gaussian::set().",
1917 o2scl_linalg::cholesky_invert<mat_t>(ndim,covar_inv);
1920 double sqrt_det=1.0;
1921 for(
size_t i=0;i<ndim;i++) {
1922 if (!std::isfinite(chol(i,i))) {
1923 O2SCL_ERR2(
"An entry of the Cholesky decomposition was not finite ",
1926 sqrt_det*=chol(i,i);
1927 for(
size_t j=0;j<ndim;j++) {
1928 if (i<j) chol(i,j)=0.0;
1934 if (!std::isfinite(norm)) {
1943 virtual double pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1945 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1949 for(
size_t i=0;i<ndim;i++) q[i]=x_A[i]-x_B[i];
1950 vtmp=prod(covar_inv,q);
1951 ret*=exp(-0.5*inner_prod(q,vtmp));
1958 virtual double log_pdf(
const vec_t &x_B,
const vec_t &x_A)
const {
1960 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1963 double ret=log(norm);
1964 for(
size_t i=0;i<ndim;i++) q[i]=x_A[i]-x_B[i];
1965 vtmp=prod(covar_inv,q);
1966 ret+=-0.5*inner_prod(q,vtmp);
1980 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1983 for (
size_t i=0;i<ndim;i++) q[i]=pdg();
1985 for (
size_t i=0;i<ndim;i++) x_A[i]=x_B[i]+vtmp[i];
1996 #ifdef O2SCL_NEVER_DEFINED 2003 template<
class vec_t=boost::numeric::ublas::vector<
double>,
2004 class mat_t=boost::numeric::ublas::matrix<
double>,
2005 class mat_col_t=boost::numeric::ublas::matrix_column<mat_t> >
2006 class prob_dens_mdim_gproc :
2014 #ifndef DOXYGEN_NO_O2NS vec_t peak
Location of the peak.
double sigma_
Width parameter.
prob_dens_mdim_bound_gaussian(size_t p_ndim, vec_t &p_peak, mat_t &covar, vec_t &p_low, vec_t &p_high)
Create a distribution with the specified peak, covariance matrix, lower limits, and upper limits...
virtual double entropy() const
Entropy of the distribution ( )
void set_gproc(size_t n_dim, size_t n_init, vec_vec_t &x, vec_t &y, func_t &fcovar)
Given a data set and a covariance function, construct probability distribution based on a Gaussian pr...
prob_cond_mdim_gaussian & operator=(const prob_cond_mdim_gaussian &pcmg)
Copy constructor with operator=.
virtual double pdf(const vec_t &x_B, const vec_t &x_A) const
The conditional probability of x_A given x_B, i.e. .
rng_gsl rg
Internal random number generator.
Gaussian distribution bounded by a hypercube.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
prob_cond_mdim_indep & operator=(const prob_cond_mdim_indep &pcmi)
Copy constructor with operator=.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
prob_cond_mdim_gaussian(size_t p_ndim, mat_t &covar)
Create a distribution from the covariance matrix.
double mean()
Get the center.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
virtual size_t dim() const
The dimensionality.
int set_internal(size_t sz, vec_t &step, vec_t &low, vec_t &high)
Internal set function.
virtual double pdf(const vec_t &v) const
Compute the normalized probability density.
virtual double pdf(const vec_t &x) const
The normalized density.
virtual size_t dim() const
The dimensionality.
Conditional probability for a random walk inside a hypercube.
void set_seed(unsigned long int s)
Set the seed.
prob_dens_gaussian & operator=(const prob_dens_gaussian &pdg)
Copy constructor with operator=.
prob_cond_mdim_gaussian(const prob_cond_mdim_gaussian &pcmg)
Copy constructor.
virtual double log_pdf(const vec_t &x_B, const vec_t &x_A) const
The log of the conditional probability of x_A given x_B i.e. .
o2scl::prob_dens_gaussian pdg
Standard normal.
prob_dens_mdim_gaussian(size_t p_mdim, size_t n_pts, const mat2_t &pts, const vec2_t &vals)
Create a distribution from a set of samples from a multidimensional Gaussian.
Generic object which represents a column of a const matrix.
void set_sigma(double sigma)
Set the Gaussian width (must be positive)
invalid argument supplied by user
prob_dens_mdim_gaussian(size_t p_ndim, vec_t &p_peak, mat_t &covar)
Create a distribution from the covariance matrix.
ubvector sum
Normalized partial sum of histogram bins.
A multi-dimensional conditional probability density function independent of the input.
size_t n
Number of original histogram bins.
const vec_t & get_peak()
Get the peak location.
vec_t vtmp
Temporary storage 2.
virtual double entropy() const
Entropy of the distribution ( )
const ubvector & bin_ranges()
Get reference to bin ranges.
virtual double level_fixed_integral(double integral)
Return the contour level corresponding to a fixed integral.
prob_dens_lognormal & operator=(const prob_dens_lognormal &pdg)
Copy constructor with operator=.
virtual void operator()(vec_t &x) const
Sample the distribution.
A class for representing permutations.
virtual double log_pdf(const vec_t &x) const
Compute the natural log of the probability density function (arbitrary normalization) ...
double sig_x
The x standard deviation.
prob_dens_lognormal()
Create a blank lognormal distribution.
int LU_decomp(const size_t N, mat_t &A, o2scl::permutation &p, int &signum)
Compute the LU decomposition of the matrix A.
o2scl::rng_gsl r
Base GSL random number generator.
o2scl::prob_dens_gaussian pdg
Standard normal.
virtual double pdf(const vec_t &x_B, const vec_t &x_A) const
The conditional probability of x_A given x_B, i.e. .
double random()
Return a random number in .
double x0
The x coordinate of the centroid.
prob_cond_mdim_fixed_step & operator=(const prob_cond_mdim_fixed_step &pcmfs)
Copy constructor with operator=.
A one-dimensional Gaussian probability density.
A multi-dimensional conditional probability density function.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
rng_gsl r
The GSL random number generator.
prob_dens_gaussian()
Create a standard normal distribution.
mat_t chol
Cholesky decomposition.
virtual size_t dim() const
The dimensionality.
virtual void operator()(const vec_t &x_B, vec_t &x_A) const
Sample the distribution.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
mat_t covar_inv
Inverse of the covariance matrix.
virtual double invert_cdf(double cdf) const
The inverse cumulative distribution function.
void set_mu_sigma(double mu, double sigma)
Set the maximum and width of the lognormal distribution.
A one-dimensional histogram class.
prob_dens_gaussian(const prob_dens_gaussian &pdg)
Copy constructor.
virtual double pdf(const vec_t &x) const
The normalized density.
const double & get_norm()
Get the normalization.
prob_dens_lognormal(const prob_dens_lognormal &pdg)
Copy constructor.
requested feature not (yet) implemented
std::vector< double > u_step
Step sizes.
vec_t vtmp
Temporary storage 2.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
virtual size_t dim() const
The dimensionality.
prob_cond_mdim_gaussian()
Create an empty distribution.
prob_dens_gaussian(double cent, double sigma)
Create a Gaussian distribution with width sigma.
A multi-dimensional Gaussian conditional probability density function.
virtual void operator()(const vec_t &x_B, vec_t &x_A) const
Sample the distribution.
double wvector_covariance(size_t n, const vec_t &data1, const vec2_t &data2, const vec3_t &weights)
The weighted covariance of two vectors.
virtual double log_pdf(const vec_t &x_B, const vec_t &x_A) const
The log of the conditional probability of x_A given x_B i.e. .
A multidimensional distribution formed by the product of several one-dimensional distributions.
double norm
Normalization factor.
prob_dens_mdim_factor(const prob_dens_mdim_factor &pdmf)
Copy constructor.
prob_dens_mdim_gaussian(const prob_dens_mdim_gaussian &pdmg)
Copy constructor.
A one-dimensional probability density function.
virtual double pdf(double x) const
The normalized density.
virtual double entropy() const
Entropy of the distribution ( )
A one-dimensional probability density over a finite range.
virtual double entropy() const
Entropy of the distribution ( )
virtual void operator()(vec_t &x) const
Sample the distribution.
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
const mat_t & get_chol()
Get the Cholesky decomposition.
A multi-dimensional probability density function.
const mat_t & get_covar_inv()
Get the inverse of the covariance matrix.
void set_seed(unsigned long int s)
Set the random number generator seed.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
virtual void operator()(const vec_t &x_B, vec_t &x_A) const
Sample the distribution.
int cholesky_decomp(const size_t M, mat_t &A, bool err_on_fail=true)
Compute the in-place Cholesky decomposition of a symmetric positive-definite square matrix...
ubvector range
Vector specifying original histogram bins.
size_t samp_max
Maximum number of samples.
A bivariate gaussian probability distribution.
double cent_
Central value.
virtual double pdf(const vec_t &x) const
The normalized density.
virtual size_t dim() const
Return the dimensionality.
search_vec< ubvector > sv
Search through the partial sums.
virtual double log_pdf(double x) const
The log of the normalized density.
virtual double operator()() const
Sample from the specified density.
std::vector< double > u_high
Upper limits.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
virtual double operator()() const
Sample from the specified density.
double sigma_
Width parameter.
vec_t q
Temporary storage 1.
virtual double log_pdf(double x) const
The log of the normalized density.
virtual double pdf(const vec_t &x_B, const vec_t &x_A) const
The conditional probability of x_A given x_B, i.e. .
double rho
The covariance.
A one-dimensional probability density over the positive real numbers.
prob_dens_mdim_biv_gaussian(const prob_dens_mdim_biv_gaussian &pdmbg)
Copy constructor.
vec_t q
Temporary storage 1.
virtual double log_metrop_hast(const vec_t &x, vec_t &x_prime) const
Sample the distribution and return the log of the Metropolis-Hastings ratio.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
A multi-dimensional Gaussian probability density function using a Cholesky decomposition.
virtual size_t dim() const
The dimensionality.
prob_cond_mdim_indep(prob_dens_mdim< vec_t > &out)
Create a conditional probability distribution based on the specified probability distribution.
Random number generator (GSL)
virtual void operator()(vec_t &x) const
Sample the distribution.
Searching class for monotonic data with caching.
prob_dens_lognormal(double mu, double sigma)
Create lognormal distribution with mean parameter mu and width parameter sigma.
double norm
Normalization factor, .
virtual void contour(double level, double theta, vec_t &x)
Return a point on the contour for a specified level given an angle.
rng_gsl rng
Random number generator.
Probability density function based on a histogram.
void set_alt(size_t p_ndim, vec_t &p_peak, mat_t &p_chol, mat_t &p_covar_inv, double p_norm)
Alternate set function for use when covariance matrix has already been decomposed and inverted...
size_t ndim
Number of dimensions.
std::vector< prob_dens_func > list
Vector of one-dimensional distributions.
virtual double pdf(double x) const
The normalized density.
virtual double log_pdf(const vec_t &x_B, const vec_t &x_A) const
The log of the conditional probability of x_A given x_B i.e. .
int diagonal_has_zero(const size_t N, mat_t &A)
Return 1 if at least one of the elements in the diagonal is zero.
virtual size_t dim() const
Return the dimensionality.
virtual double operator()() const
Sample from the specified density.
mat_t covar_inv
Inverse of the covariance matrix.
virtual double pdf(const vec_t &x) const
Compute the probability density function (arbitrary normalization)
void set_seed(unsigned long int s)
Set the seed.
prob_dens_mdim_gaussian & operator=(const prob_dens_mdim_gaussian &pdmg)
Copy constructor with operator=.
prob_cond_mdim_fixed_step(const prob_cond_mdim_fixed_step &pcmfs)
Copy constructor.
prob_cond_mdim_indep(const prob_cond_mdim_indep &pcmi)
Copy constructor.
void set_seed(unsigned long int s)
Set the seed.
std::vector< double > u_low
Lower limits.
virtual void operator()(vec_t &x) const
Sample the distribution.
size_t ndim
Number of dimensions.
virtual double pdf(double x) const
The normalized density.
virtual double log_pdf(double x) const
The log of the normalized density.
prob_dens_mdim_gaussian()
Create an empty distribution.
void set_seed(unsigned long int s)
Set the seed.
const ubvector & partial_sums()
Get reference to partial sums.
double sig_y
The y standard deviation.
void set(size_t p_ndim, vec_t &p_peak, mat_t &covar)
Set the peak and covariance matrix for the distribution.
Lognormal density function.
prob_dens_mdim_bound_gaussian()
Create an empty distribution.
double stddev()
Get the Gaussian width.
mat_t chol
Cholesky decomposition.
void set_center(double cent)
Set the center.
double y0
The y coordinate of the centroid.