2
0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00

Add logpdf to chi squared distribution

[ci skip]
This commit is contained in:
Matt Borland
2022-02-24 12:19:15 +01:00
parent 6f1cda5bb3
commit 93593fb8f2
2 changed files with 45 additions and 0 deletions

View File

@@ -131,6 +131,47 @@ RealType pdf(const chi_squared_distribution<RealType, Policy>& dist, const RealT
return gamma_p_derivative(degrees_of_freedom / 2, chi_square / 2, Policy()) / 2;
} // pdf
template <class RealType, class Policy>
RealType logpdf(const chi_squared_distribution<RealType, Policy>& dist, const RealType& chi_square)
{
BOOST_MATH_STD_USING // for ADL of std functions
RealType k = dist.degrees_of_freedom();
// Error check:
RealType error_result;
static const char* function = "boost::math::logpdf(const chi_squared_distribution<%1%>&, %1%)";
if(false == detail::check_df(function, k, &error_result, Policy()))
{
return error_result;
}
if((chi_square < 0) || !(boost::math::isfinite)(chi_square))
{
return policies::raise_domain_error<RealType>(
function, "Chi Square parameter was %1%, but must be > 0 !", chi_square, Policy());
}
if(chi_square == 0)
{
// Handle special cases:
if(k < 2)
{
return policies::raise_overflow_error<RealType>(function, 0, Policy());
}
else if(k == 2)
{
return -boost::math::constants::ln_two<RealType>();
}
else
{
return -std::numeric_limits<RealType>::infinity();
}
}
return log(pdf(dist, chi_square));
} // logpdf
template <class RealType, class Policy>
inline RealType cdf(const chi_squared_distribution<RealType, Policy>& dist, const RealType& chi_square)
{

View File

@@ -36,6 +36,8 @@ using std::cout;
using std::endl;
#include <limits>
using std::numeric_limits;
#include <cmath>
using std::log;
template <class RealType>
RealType naive_pdf(RealType df, RealType x)
@@ -59,6 +61,8 @@ void test_spot(
cdf(dist, cs), P, tol);
BOOST_CHECK_CLOSE(
pdf(dist, cs), naive_pdf(dist.degrees_of_freedom(), cs), tol);
BOOST_CHECK_CLOSE(
logpdf(dist, cs), log(pdf(dist, cs)), tol);
if((P < 0.99) && (Q < 0.99))
{
//