Negative Binomial Distribution: degenerate cases

[HS <tan@…>]

In the following calculations (please see code included below), the 4 results cannot be all correct at the same time.

Also, if one will to uncomment the last block of code, one would find that its calculation would fail.

Thoughts:

• Now when the probability of "success" p0 (you used p in your documentation) is 1, shouldn't the value of the corresponding random variable be simply 0? If so, shouldn't the quantiles at p=0 and p=1 be both at 0?

• Similarly, when p0 is 0, I would think it is quite "natural" to set/define the corresponding random variable to be one with mass 1 at infinity ... no need to work so hard to search for the zero as the error to the final calculation shows.

Thank you.
HS

P.S. Please copy the included code to a test.cpp and then build with g++ -o test test.cpp.

P.ps. If the evaluation policy will to change to ignore the error, it appears that the calculation never complete. Can something be done to prevent this? Suggestions?

#define BOOST_MATH_DOMAIN_ERROR_POLICY ignore_error
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
// #define BOOST_MATH_EVALUATION_ERROR_POLICY ignore_error
#define BOOST_MATH_DISCRETE_QUANTILE_POLICY integer_round_up

#include <cstdlib>
#include <cstdio>

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

int main()
{
   double _r0 = 3.0;
   double _p0 = 1.0;

   boost::math::negative_binomial_distribution<> _dist( _r0, _p0 );

   std::printf(
         "\n"
         "   r = %g\n"
         "   p = %g\n\n",
         _r0,
         _p0 );

   std::printf(
         "   Quantile(%g;upperTail=false) = %g\n",
         1.0,
         boost::math::quantile( _dist, 1.0 ) );

   std::printf(
         "   Quantile(%g;upperTail=true) = %g\n\n",
         1.0,
         boost::math::quantile(
            boost::math::complement(
               _dist,
               1.0 ) ) );

   std::printf(
         "   Quantile(%g;upperTail=false) = %g\n",
         0.0,
         boost::math::quantile( _dist, 0.0 ) );

   std::printf(
         "   Quantile(%g;upperTail=true) = %g\n\n",
         0.0,
         boost::math::quantile(
            boost::math::complement(
               _dist,
               0.0 ) ) );

   /*
   {
      _p0 = 0.0;

      boost::math::negative_binomial_distribution<> _dist( _r0, _p0 );

      std::printf(
            "\n"
            "   r = %g\n"
            "   p = %g\n\n",
            _r0,
            _p0 );

      std::printf(
            "   Quantile(%g;upperTail=false) = %g\n",
            0.32,
            boost::math::quantile( _dist, 0.32 ) );
   }
   */

   return EXIT_SUCCESS;
}

[John Maddock]

See ​https://github.com/boostorg/math/commit/0c01f682eb9c562d3cca6372c0aa8e68951decd1

I haven't made the quantile for probability 1 zero in the corner case you gave: rather it's always at infinity. I think this is more consistent, though there are an infinite number of "correct" answers for this case ;-)

[HS <tan@…>]

Replying to johnmaddock:

See ​https://github.com/boostorg/math/commit/0c01f682eb9c562d3cca6372c0aa8e68951decd1

I haven't made the quantile for probability 1 zero in the corner case you gave: rather it's always at infinity. I think this is more consistent, though there are an infinite number of "correct" answers for this case ;-)

Dear John,

Thank you (and your team) for the great work!

I will not quibble with you about this---it is an edge case which is generally not very interesting in the applied world. However, as a teacher (was) of probability classes, the inverse CDF (or quantile) at p of a distribution F, to us, is the infimum of all x such that F(x) >= p; in this degenerate case of a discrete distribution, this is 0. And one will also find this is so on R too:

> qnbinom(1,3,1,lower.tail=FALSE)
[1] 0
> qnbinom(1,3,1,lower.tail=TRUE)
[1] 0
> qnbinom(0,3,1,lower.tail=FALSE)
[1] 0
> qnbinom(0,3,1,lower.tail=TRUE)
[1] 0

Mathmatica 9, on the other hand, says the following:

In[1]:= Quantile[NegativeBinomialDistribution[3,1],0]
Out[1]= 0

In[2]:= Quantile[NegativeBinomialDistribution[3,1],1]
Out[2]= Infinity

But it is very easy for anyone to work around the decision here if he or she disagrees with it. :)

With best regards,
HS