boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00
Code Issues Projects Releases Wiki Activity

Add logcdf specialization to pareto distribution

Browse Source
This commit is contained in:
Matt Borland
2023-02-11 17:49:29 -08:00
parent e387f9c429
commit 614445505a
2 changed files with 70 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files

48
include/boost/math/distributions/pareto.hpp
View File

@@ -1,5 +1,6 @@
// Copyright John Maddock 2007.
// Copyright Paul A. Bristow 2007, 2009
// Copyright Matt Borland 2023.
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0. (See accompanying file
// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
@@ -20,8 +21,11 @@
#include <boost/math/distributions/complement.hpp>
#include <boost/math/distributions/detail/common_error_handling.hpp>
#include <boost/math/special_functions/powm1.hpp>
#include <boost/math/special_functions/log1p.hpp>
#include <utility> // for BOOST_CURRENT_VALUE?
#include <limits>
#include <cmath>
namespace boost
{
@@ -225,6 +229,28 @@ namespace boost
return result;
} // cdf
template <class RealType, class Policy>
inline RealType logcdf(const pareto_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // for ADL of std function pow.
static const char* function = "boost::math::logcdf(const pareto_distribution<%1%>&, %1%)";
RealType scale = dist.scale();
RealType shape = dist.shape();
RealType result = 0;
if(false == (detail::check_pareto_x(function, x, &result, Policy())
&& detail::check_pareto(function, scale, shape, &result, Policy())))
return result;
if (x <= scale)
{ // regardless of shape, cdf is zero.
return -std::numeric_limits<RealType>::infinity();
}
result = log1p(-pow(scale/x, shape), Policy());
return result;
} // logcdf
template <class RealType, class Policy>
inline RealType quantile(const pareto_distribution<RealType, Policy>& dist, const RealType& p)
{
@@ -274,6 +300,28 @@ namespace boost
return result;
} // cdf complement
template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<pareto_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING // for ADL of std function pow.
static const char* function = "boost::math::logcdf(const pareto_distribution<%1%>&, %1%)";
RealType result = 0;
RealType x = c.param;
RealType scale = c.dist.scale();
RealType shape = c.dist.shape();
if(false == (detail::check_pareto_x(function, x, &result, Policy())
&& detail::check_pareto(function, scale, shape, &result, Policy())))
return result;
if (x <= scale)
{ // regardless of shape, cdf is zero, and complement is unity.
return 0;
}
result = log(pow((scale/x), shape));
return result;
} // logcdf complement
template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<pareto_distribution<RealType, Policy>, RealType>& c)
{

22
test/test_pareto.cpp
View File

@@ -1,5 +1,6 @@
// Copyright Paul A. Bristow 2007, 2009.
// Copyright John Maddock 2006.
// Copyright Matt Borland 2023.
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0.
// (See accompanying file LICENSE_1_0.txt
@@ -40,10 +41,18 @@
using std::setprecision;
#include <limits>
using std::numeric_limits;
#include <type_traits>
template <class RealType>
void check_pareto(RealType scale, RealType shape, RealType x, RealType p, RealType q, RealType tol)
{
RealType logtol = tol * 10;
BOOST_IF_CONSTEXPR (std::is_same<RealType, long double>::value ||
std::is_same<RealType, boost::math::concepts::real_concept>::value)
{
logtol *= 100;
}
BOOST_CHECK_CLOSE_FRACTION(
::boost::math::cdf(
pareto_distribution<RealType>(scale, shape), // distribution.
@@ -57,6 +66,19 @@
x)), // random variable.
q, // probability complement.
tol); // tolerance eps.
BOOST_CHECK_CLOSE_FRACTION(
::boost::math::logcdf(
pareto_distribution<RealType>(scale, shape), // distribution.
x), // random variable.
log(p), // probability.
logtol); // tolerance eps.
BOOST_CHECK_CLOSE_FRACTION(
::boost::math::logcdf(
complement(
pareto_distribution<RealType>(scale, shape), // distribution.
x)), // random variable.
log(q), // probability complement.
logtol); // tolerance eps.
BOOST_CHECK_CLOSE_FRACTION(
::boost::math::quantile(
pareto_distribution<RealType>(scale, shape), // distribution.