Skip to content

Commit

Permalink
Merge pull request #946 from mborland/logcdf
Browse files Browse the repository at this point in the history
Add logcdf to distributions
  • Loading branch information
mborland committed Apr 13, 2023
2 parents 9a3b8bc + 3886633 commit 93448ac
Show file tree
Hide file tree
Showing 17 changed files with 842 additions and 36 deletions.
15 changes: 15 additions & 0 deletions include/boost/math/distributions/detail/derived_accessors.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -122,6 +122,13 @@ inline typename Distribution::value_type cdf(const Distribution& dist, const Rea
typedef typename Distribution::value_type value_type;
return cdf(dist, static_cast<value_type>(x));
}
template <class Distribution, class Realtype>
inline typename Distribution::value_type logcdf(const Distribution& dist, const Realtype& x)
{
using std::log;
using value_type = typename Distribution::value_type;
return log(cdf(dist, static_cast<value_type>(x)));
}
template <class Distribution, class RealType>
inline typename Distribution::value_type quantile(const Distribution& dist, const RealType& x)
{
Expand All @@ -143,6 +150,14 @@ inline typename Distribution::value_type cdf(const complemented2_type<Distributi
return cdf(complement(c.dist, static_cast<value_type>(c.param)));
}

template <class Distribution, class RealType>
inline typename Distribution::value_type logcdf(const complemented2_type<Distribution, RealType>& c)
{
using std::log;
typedef typename Distribution::value_type value_type;
return log(cdf(complement(c.dist, static_cast<value_type>(c.param))));
}

template <class Distribution, class RealType>
inline typename Distribution::value_type quantile(const complemented2_type<Distribution, RealType>& c)
{
Expand Down
39 changes: 39 additions & 0 deletions include/boost/math/distributions/exponential.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -163,6 +163,24 @@ inline RealType cdf(const exponential_distribution<RealType, Policy>& dist, cons
return result;
} // cdf

template <class RealType, class Policy>
inline RealType logcdf(const exponential_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // for ADL of std functions

static const char* function = "boost::math::logcdf(const exponential_distribution<%1%>&, %1%)";

RealType result = 0;
RealType lambda = dist.lambda();
if(0 == detail::verify_lambda(function, lambda, &result, Policy()))
return result;
if(0 == detail::verify_exp_x(function, x, &result, Policy()))
return result;
result = boost::math::log1p(-exp(-x * lambda), Policy());

return result;
} // cdf

template <class RealType, class Policy>
inline RealType quantile(const exponential_distribution<RealType, Policy>& dist, const RealType& p)
{
Expand Down Expand Up @@ -207,6 +225,27 @@ inline RealType cdf(const complemented2_type<exponential_distribution<RealType,
return result;
}

template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<exponential_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING // for ADL of std functions

static const char* function = "boost::math::logcdf(const exponential_distribution<%1%>&, %1%)";

RealType result = 0;
RealType lambda = c.dist.lambda();
if(0 == detail::verify_lambda(function, lambda, &result, Policy()))
return result;
if(0 == detail::verify_exp_x(function, c.param, &result, Policy()))
return result;
// Workaround for VC11/12 bug:
if (c.param >= tools::max_value<RealType>())
return 0;
result = -c.param * lambda;

return result;
}

template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<exponential_distribution<RealType, Policy>, RealType>& c)
{
Expand Down
50 changes: 50 additions & 0 deletions include/boost/math/distributions/extreme_value.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -174,6 +174,32 @@ inline RealType cdf(const extreme_value_distribution<RealType, Policy>& dist, co
return result;
} // cdf

template <class RealType, class Policy>
inline RealType logcdf(const extreme_value_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // for ADL of std functions

static const char* function = "boost::math::logcdf(const extreme_value_distribution<%1%>&, %1%)";

if((boost::math::isinf)(x))
return x < 0 ? 0.0f : 1.0f;
RealType a = dist.location();
RealType b = dist.scale();
RealType result = 0;
if(0 == detail::verify_scale_b(function, b, &result, Policy()))
return result;
if(0 == detail::check_finite(function, a, &result, Policy()))
return result;
if(0 == detail::check_finite(function, a, &result, Policy()))
return result;
if(0 == detail::check_x("boost::math::logcdf(const extreme_value_distribution<%1%>&, %1%)", x, &result, Policy()))
return result;

result = -exp((a-x)/b);

return result;
} // logcdf

template <class RealType, class Policy>
RealType quantile(const extreme_value_distribution<RealType, Policy>& dist, const RealType& p)
{
Expand Down Expand Up @@ -225,6 +251,30 @@ inline RealType cdf(const complemented2_type<extreme_value_distribution<RealType
return result;
}

template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<extreme_value_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING // for ADL of std functions

static const char* function = "boost::math::logcdf(const extreme_value_distribution<%1%>&, %1%)";

if((boost::math::isinf)(c.param))
return c.param < 0 ? 1.0f : 0.0f;
RealType a = c.dist.location();
RealType b = c.dist.scale();
RealType result = 0;
if(0 == detail::verify_scale_b(function, b, &result, Policy()))
return result;
if(0 == detail::check_finite(function, a, &result, Policy()))
return result;
if(0 == detail::check_x(function, c.param, &result, Policy()))
return result;

result = log1p(-exp(-exp((a-c.param)/b)), Policy());

return result;
}

template <class RealType, class Policy>
RealType quantile(const complemented2_type<extreme_value_distribution<RealType, Policy>, RealType>& c)
{
Expand Down
62 changes: 59 additions & 3 deletions include/boost/math/distributions/geometric.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -43,9 +43,11 @@
#include <boost/math/special_functions/fpclassify.hpp> // isnan.
#include <boost/math/tools/roots.hpp> // for root finding.
#include <boost/math/distributions/detail/inv_discrete_quantile.hpp>
#include <boost/math/special_functions/log1p.hpp>

#include <limits> // using std::numeric_limits;
#include <utility>
#include <cmath>

#if defined (BOOST_MSVC)
# pragma warning(push)
Expand Down Expand Up @@ -378,9 +380,39 @@ namespace boost
return probability;
} // cdf Cumulative Distribution Function geometric.

template <class RealType, class Policy>
inline RealType cdf(const complemented2_type<geometric_distribution<RealType, Policy>, RealType>& c)
{ // Complemented Cumulative Distribution Function geometric.
template <class RealType, class Policy>
inline RealType logcdf(const geometric_distribution<RealType, Policy>& dist, const RealType& k)
{ // Cumulative Distribution Function of geometric.
using std::pow;
static const char* function = "boost::math::logcdf(const geometric_distribution<%1%>&, %1%)";

// k argument may be integral, signed, or unsigned, or floating point.
// If necessary, it has already been promoted from an integral type.
RealType p = dist.success_fraction();
// Error check:
RealType result = 0;
if(false == geometric_detail::check_dist_and_k(
function,
p,
k,
&result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}
if(k == 0)
{
return log(p); // success_fraction
}
//RealType q = 1 - p; // Bad for small p
//RealType probability = 1 - std::pow(q, k+1);

RealType z = boost::math::log1p(-p, Policy()) * (k + 1);
return log1p(-exp(z), Policy());
} // logcdf Cumulative Distribution Function geometric.

template <class RealType, class Policy>
inline RealType cdf(const complemented2_type<geometric_distribution<RealType, Policy>, RealType>& c)
{ // Complemented Cumulative Distribution Function geometric.
BOOST_MATH_STD_USING
static const char* function = "boost::math::cdf(const geometric_distribution<%1%>&, %1%)";
// k argument may be integral, signed, or unsigned, or floating point.
Expand All @@ -403,6 +435,30 @@ namespace boost
return probability;
} // cdf Complemented Cumulative Distribution Function geometric.

template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<geometric_distribution<RealType, Policy>, RealType>& c)
{ // Complemented Cumulative Distribution Function geometric.
BOOST_MATH_STD_USING
static const char* function = "boost::math::logcdf(const geometric_distribution<%1%>&, %1%)";
// k argument may be integral, signed, or unsigned, or floating point.
// If necessary, it has already been promoted from an integral type.
RealType const& k = c.param;
geometric_distribution<RealType, Policy> const& dist = c.dist;
RealType p = dist.success_fraction();
// Error check:
RealType result = 0;
if(false == geometric_detail::check_dist_and_k(
function,
p,
k,
&result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}

return boost::math::log1p(-p, Policy()) * (k+1);
} // logcdf Complemented Cumulative Distribution Function geometric.

template <class RealType, class Policy>
inline RealType quantile(const geometric_distribution<RealType, Policy>& dist, const RealType& x)
{ // Quantile, percentile/100 or Percent Point geometric function.
Expand Down
85 changes: 85 additions & 0 deletions include/boost/math/distributions/laplace.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
#ifndef BOOST_STATS_LAPLACE_HPP
#define BOOST_STATS_LAPLACE_HPP

#include <boost/math/special_functions/log1p.hpp>
#include <boost/math/distributions/detail/common_error_handling.hpp>
#include <boost/math/distributions/complement.hpp>
#include <boost/math/constants/constants.hpp>
Expand Down Expand Up @@ -226,6 +227,50 @@ inline RealType cdf(const laplace_distribution<RealType, Policy>& dist, const Re
return result;
} // cdf

template <class RealType, class Policy>
inline RealType logcdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // For ADL of std functions.

RealType result = 0;
// Checking function argument.
const char* function = "boost::math::logcdf(const laplace_distribution<%1%>&, %1%)";
// Check scale and location.
if (false == dist.check_parameters(function, &result))
{
return result;
}

// Special cdf values:
if((boost::math::isinf)(x))
{
if(x < 0)
{
return 0; // -infinity.
}
return 1; // + infinity.
}

if (false == detail::check_x(function, x, &result, Policy()))
{
return result;
}

// General cdf values
RealType scale( dist.scale() );
RealType location( dist.location() );

if (x < location)
{
result = ((x - location) / scale) - boost::math::constants::ln_two<RealType>();
}
else
{
result = log1p(-exp((location - x) / scale) / 2);
}

return result;
} // logcdf

template <class RealType, class Policy>
inline RealType quantile(const laplace_distribution<RealType, Policy>& dist, const RealType& p)
Expand Down Expand Up @@ -302,6 +347,46 @@ inline RealType cdf(const complemented2_type<laplace_distribution<RealType, Poli
return result;
} // cdf complement

template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
{
// Calculate complement of logcdf.
BOOST_MATH_STD_USING // for ADL of std functions

RealType scale = c.dist.scale();
RealType location = c.dist.location();
RealType x = c.param;
RealType result = 0;

// Checking function argument.
const char* function = "boost::math::logcdf(const complemented2_type<laplace_distribution<%1%>, %1%>&)";

// Check scale and location.
if (false == c.dist.check_parameters(function, &result)) return result;

// Special cdf values.
if((boost::math::isinf)(x))
{
if(x < 0)
{
return 1; // cdf complement -infinity is unity.
}

return 0; // cdf complement +infinity is zero.
}
if(false == detail::check_x(function, x, &result, Policy()))return result;

// Cdf interval value.
if (-x < -location)
{
result = (-x+location)/scale - boost::math::constants::ln_two<RealType>();
}
else
{
result = log1p(-exp( (-location+x)/scale )/2, Policy());
}
return result;
} // cdf complement

template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
Expand Down
Loading

0 comments on commit 93448ac

Please sign in to comment.