Opened 9 years ago

Closed 9 years ago

#9126 closed Bugs (fixed)

Logistic distribution pdf() and cdf(complement()) fail to catch invalid scale and location parameters

Reported by: Paul McClellan <paulm@…> Owned by: Paul A. Bristow
Milestone: Boost 1.55.0 Component: math
Version: Boost 1.54.0 Severity: Problem
Keywords: Cc:

Description

This issue is related to Ticket #9042

In boost/math/distributions/logistic.hpp:

    template <class RealType, class Policy>
    inline RealType pdf(const logistic_distribution<RealType, Policy>& dist, const RealType& x)
    {
       RealType scale = dist.scale();
       RealType location = dist.location();

       static const char* function = "boost::math::pdf(const logistic_distribution<%1%>&, %1%)";
       if((boost::math::isinf)(x))
       {
          return 0; // pdf + and - infinity is zero.
       }

       RealType result = 0;
       if(false == detail::check_scale(function, scale , &result, Policy()))
       {
          return result;
       }
       if(false == detail::check_location(function, location, &result, Policy()))
       {
          return result;
       }
       if(false == detail::check_x(function, x, &result, Policy()))
       {
          return result;
       }

       BOOST_MATH_STD_USING
       RealType exp_term = (location - x) / scale;
       if(fabs(exp_term) > tools::log_max_value<RealType>())
          return 0;
       exp_term = exp(exp_term);
       if((exp_term * scale > 1) && (exp_term > tools::max_value<RealType>() / (scale * exp_term)))
          return 1 / (scale * exp_term);
       return (exp_term) / (scale * (1 + exp_term) * (1 + exp_term));
    } 

The test if((boost::math::isinf)(x)) occurs before check_scale() or check_location() is called, returning 0 for infinite x, even if scale or location are invalid.

The calls to check_scale() and check_location() should be moved up just before the test if((boost::math::isinf)(x)).

See, for example, how the tests are ordered here:

    template <class RealType, class Policy>
    inline RealType cdf(const logistic_distribution<RealType, Policy>& dist, const RealType& x)
    {
       RealType scale = dist.scale();
       RealType location = dist.location();
       RealType result = 0; // of checks.
       static const char* function = "boost::math::cdf(const logistic_distribution<%1%>&, %1%)";
       if(false == detail::check_scale(function, scale, &result, Policy()))
       {
          return result;
       }
       if(false == detail::check_location(function, location, &result, Policy()))
       {
          return result;
       }

       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;
       }
       BOOST_MATH_STD_USING
       RealType power = (location - x) / scale;
       if(power > tools::log_max_value<RealType>())
          return 0;
       if(power < -tools::log_max_value<RealType>())
          return 1;
       return 1 / (1 + exp(power)); 
    } 

The same issue exists with the complementary cdf:

    template <class RealType, class Policy>
    inline RealType cdf(const complemented2_type<logistic_distribution<RealType, Policy>, RealType>& c)
    {
       BOOST_MATH_STD_USING
       RealType location = c.dist.location();
       RealType scale = c.dist.scale();
       RealType x = c.param;
       static const char* function = "boost::math::cdf(const complement(logistic_distribution<%1%>&), %1%)";

       if((boost::math::isinf)(x))
       {
          if(x < 0) return 1; // cdf complement -infinity is unity.
          return 0; // cdf complement +infinity is zero
       }
       RealType result = 0;
       if(false == detail::check_scale(function, scale, &result, Policy()))
          return result;
       if(false == detail::check_location(function, location, &result, Policy()))
          return result;
       if(false == detail::check_x(function, x, &result, Policy()))
          return result;
       RealType power = (x - location) / scale;
       if(power > tools::log_max_value<RealType>())
          return 0;
       if(power < -tools::log_max_value<RealType>())
          return 1;
       return 1 / (1 + exp(power)); 
    } 

Change History (9)

comment:1 by Paul A. Bristow, 9 years ago

Milestone: To Be DeterminedBoost 1.55.0
Owner: changed from John Maddock to Paul A. Bristow
Status: newassigned

Thanks for reporting this which should be fixed for 1.55.

comment:2 by Paul A. Bristow, 9 years ago

Have you an actual use example or test case for this condition (or is this comment just from examination of the code)?

comment:3 by Paul McClellan <paulm@…>, 9 years ago

I have wrappers around boost:

#include <boost/math/distributions/logistic.hpp>

double logistic_pdf( double dX, double dL, double dS )
{
	// logistic distribution object:
	boost::math::logistic logis(dL, dS);
	return pdf(logis, dX);
}

double logistic_utp( double dX, double dL, double dS )
{
	// logistic distribution object:
	boost::math::logistic logis(dL, dS);
	return cdf(complement(logis, dX));
}

These test cases should error, but return 0.0:

logistic_pdf(Pc_PosInf, Pc_PosInf, 1.0);
logistic_pdf(Pc_PosInf, 0.0, Pc_PosInf);
logistic_pdf(Pc_PosInf, 0.0, 0.0);

logistic_utp(Pc_PosInf, Pc_PosInf, 1.0);
logistic_utp(Pc_PosInf, 0.0, Pc_PosInf);
logistic_utp(Pc_PosInf, 0.0, 0.0);

comment:4 by Paul McClellan <paulm@…>, 9 years ago

I should have included:

#define Pc_PosInf (std::numeric_limits<double>::infinity())

comment:5 by Paul A. Bristow, 9 years ago

Fixed in boost trunk at https://svn.boost.org/svn/boost/trunk/boost/math/distributions/logistic.hpp

(you could download this version).

SVN commit 87950

but I am still puzzled why the checks in the distribution object constructor didn't catch these tests.

Did you change the default policy (from 'to throw exception') to ignore errors?

Did you use try'n'catch blocks so that the error messages appear? (as always recommended).

Did you get the same in debug and release (optimised) version?

Should make boost release 1.55.

comment:6 by Paul McClellan <paulm@…>, 9 years ago

I apologize, I should have included my user.hpp settings. I've tried to attach my user.hpp file but Trac rejects it as spam. I've included the relevant portions here:

//
// Default policies follow:
//
// Domain errors:
//
// #define BOOST_MATH_DOMAIN_ERROR_POLICY throw_on_error
#define BOOST_MATH_DOMAIN_ERROR_POLICY ignore_error
//
// Pole errors:
//
// #define BOOST_MATH_POLE_ERROR_POLICY throw_on_error
#define BOOST_MATH_POLE_ERROR_POLICY ignore_error
//
// Overflow Errors:
//
// #define BOOST_MATH_OVERFLOW_ERROR_POLICY throw_on_error
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
//
// Internal Evaluation Errors:
//
// #define BOOST_MATH_EVALUATION_ERROR_POLICY throw_on_error
#define BOOST_MATH_EVALUATION_ERROR_POLICY ignore_error
//
// Underfow:
//
// #define BOOST_MATH_UNDERFLOW_ERROR_POLICY ignore_error
//
// Denorms:
//
// #define BOOST_MATH_DENORM_ERROR_POLICY ignore_error
//
// Max digits to use for internal calculations:
//
// #define BOOST_MATH_DIGITS10_POLICY 0
//
// Promote floats to doubles internally?
//
// #define BOOST_MATH_PROMOTE_FLOAT_POLICY true
//
// Promote doubles to long double internally:
//
// #define BOOST_MATH_PROMOTE_DOUBLE_POLICY true
//
// What do discrete quantiles return?
//
// #define BOOST_MATH_DISCRETE_QUANTILE_POLICY integer_round_outwards
//
// If a function is mathematically undefined
// (for example the Cauchy distribution has no mean),
// then do we stop the code from compiling?
//
// #define BOOST_MATH_ASSERT_UNDEFINED_POLICY true
//
// Maximum series iterstions permitted:
//
// #define BOOST_MATH_MAX_SERIES_ITERATION_POLICY 1000000
//
// Maximum root finding steps permitted:
//
// define BOOST_MATH_MAX_ROOT_ITERATION_POLICY 200

// [PJM] Needed for students_t_ltq(), students_t_utq()
// [PJM] changed with fix to Ticket #8837
//#define BOOST_MATH_ROUNDING_ERROR_POLICY ignore_error

comment:7 by Paul A. Bristow, 9 years ago

Ah ha, as I eventually suspected, you are ignoring the carefully constructed error reporting mechanism :-)

I trust that the code changes will meet your wishes, but one could argue that the more correct result from 'bad' distribution constructor arguments is always NaN (not infinity). You could define your custom policy to achieve just that. Or you could use the default error handling policy and throw'n'catch from the distribution construction - the intended usage, as given in many examples.

comment:8 by John Maddock, 9 years ago

(In [85987]) Merge accumulated patches from Trunk. Refs #8384, Refs #8855, refs #9107, refs #9109, refs #8333, refs #8621, refs #8732, refs #8733, refs #8837, refs #8940, refs #9042, refs #9087, refs #9104, refs #9126.

comment:9 by John Maddock, 9 years ago

Resolution: fixed
Status: assignedclosed
Note: See TracTickets for help on using tickets.