ceph/src/boost/libs/math/include/boost/math/distributions/hypergeometric.hpp

   1 // Copyright 2008 Gautam Sewani
   2 // Copyright 2008 John Maddock
   3 //
   4 // Use, modification and distribution are subject to the
   5 // Boost Software License, Version 1.0.
   6 // (See accompanying file LICENSE_1_0.txt
   7 // or copy at http://www.boost.org/LICENSE_1_0.txt)
   8
   9 #ifndef BOOST_MATH_DISTRIBUTIONS_HYPERGEOMETRIC_HPP
  10 #define BOOST_MATH_DISTRIBUTIONS_HYPERGEOMETRIC_HPP
  11
  12 #include <boost/math/distributions/detail/common_error_handling.hpp>
  13 #include <boost/math/distributions/complement.hpp>
  14 #include <boost/math/distributions/detail/hypergeometric_pdf.hpp>
  15 #include <boost/math/distributions/detail/hypergeometric_cdf.hpp>
  16 #include <boost/math/distributions/detail/hypergeometric_quantile.hpp>
  17 #include <boost/math/special_functions/fpclassify.hpp>
  18
  19
  20 namespace boost { namespace math {
  21
  22    template <class RealType = double, class Policy = policies::policy<> >
  23    class hypergeometric_distribution
  24    {
  25    public:
  26       typedef RealType value_type;
  27       typedef Policy policy_type;
  28
  29       hypergeometric_distribution(unsigned r, unsigned n, unsigned N) // Constructor.
  30          : m_n(n), m_N(N), m_r(r)
  31       {
  32          static const char* function = "boost::math::hypergeometric_distribution<%1%>::hypergeometric_distribution";
  33          RealType ret;
  34          check_params(function, &ret);
  35       }
  36       // Accessor functions.
  37       unsigned total()const
  38       {
  39          return m_N;
  40       }
  41
  42       unsigned defective()const
  43       {
  44          return m_r;
  45       }
  46
  47       unsigned sample_count()const
  48       {
  49          return m_n;
  50       }
  51
  52       bool check_params(const char* function, RealType* result)const
  53       {
  54          if(m_r > m_N)
  55          {
  56             *result = boost::math::policies::raise_domain_error<RealType>(
  57                function, "Parameter r out of range: must be <= N but got %1%", static_cast<RealType>(m_r), Policy());
  58             return false;
  59          }
  60          if(m_n > m_N)
  61          {
  62             *result = boost::math::policies::raise_domain_error<RealType>(
  63                function, "Parameter n out of range: must be <= N but got %1%", static_cast<RealType>(m_n), Policy());
  64             return false;
  65          }
  66          return true;
  67       }
  68       bool check_x(unsigned x, const char* function, RealType* result)const
  69       {
  70          if(x < static_cast<unsigned>((std::max)(0, (int)(m_n + m_r) - (int)(m_N))))
  71          {
  72             *result = boost::math::policies::raise_domain_error<RealType>(
  73                function, "Random variable out of range: must be > 0 and > m + r - N but got %1%", static_cast<RealType>(x), Policy());
  74             return false;
  75          }
  76          if(x > (std::min)(m_r, m_n))
  77          {
  78             *result = boost::math::policies::raise_domain_error<RealType>(
  79                function, "Random variable out of range: must be less than both n and r but got %1%", static_cast<RealType>(x), Policy());
  80             return false;
  81          }
  82          return true;
  83       }
  84
  85    private:
  86       // Data members:
  87       unsigned m_n;  // number of items picked
  88       unsigned m_N; // number of "total" items
  89       unsigned m_r; // number of "defective" items
  90
  91    }; // class hypergeometric_distribution
  92
  93    typedef hypergeometric_distribution<double> hypergeometric;
  94
  95    template <class RealType, class Policy>
  96    inline const std::pair<unsigned, unsigned> range(const hypergeometric_distribution<RealType, Policy>& dist)
  97    { // Range of permissible values for random variable x.
  98 #ifdef BOOST_MSVC
  99 #  pragma warning(push)
 100 #  pragma warning(disable:4267)
 101 #endif
 102       unsigned r = dist.defective();
 103       unsigned n = dist.sample_count();
 104       unsigned N = dist.total();
 105       unsigned l = static_cast<unsigned>((std::max)(0, (int)(n + r) - (int)(N)));
 106       unsigned u = (std::min)(r, n);
 107       return std::pair<unsigned, unsigned>(l, u);
 108 #ifdef BOOST_MSVC
 109 #  pragma warning(pop)
 110 #endif
 111    }
 112
 113    template <class RealType, class Policy>
 114    inline const std::pair<unsigned, unsigned> support(const hypergeometric_distribution<RealType, Policy>& d)
 115    {
 116       return range(d);
 117    }
 118
 119    template <class RealType, class Policy>
 120    inline RealType pdf(const hypergeometric_distribution<RealType, Policy>& dist, const unsigned& x)
 121    {
 122       static const char* function = "boost::math::pdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 123       RealType result = 0;
 124       if(!dist.check_params(function, &result))
 125          return result;
 126       if(!dist.check_x(x, function, &result))
 127          return result;
 128
 129       return boost::math::detail::hypergeometric_pdf<RealType>(
 130          x, dist.defective(), dist.sample_count(), dist.total(), Policy());
 131    }
 132
 133    template <class RealType, class Policy, class U>
 134    inline RealType pdf(const hypergeometric_distribution<RealType, Policy>& dist, const U& x)
 135    {
 136       BOOST_MATH_STD_USING
 137       static const char* function = "boost::math::pdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 138       RealType r = static_cast<RealType>(x);
 139       unsigned u = itrunc(r, typename policies::normalise<Policy, policies::rounding_error<policies::ignore_error> >::type());
 140       if(u != r)
 141       {
 142          return boost::math::policies::raise_domain_error<RealType>(
 143             function, "Random variable out of range: must be an integer but got %1%", r, Policy());
 144       }
 145       return pdf(dist, u);
 146    }
 147
 148    template <class RealType, class Policy>
 149    inline RealType cdf(const hypergeometric_distribution<RealType, Policy>& dist, const unsigned& x)
 150    {
 151       static const char* function = "boost::math::cdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 152       RealType result = 0;
 153       if(!dist.check_params(function, &result))
 154          return result;
 155       if(!dist.check_x(x, function, &result))
 156          return result;
 157
 158       return boost::math::detail::hypergeometric_cdf<RealType>(
 159          x, dist.defective(), dist.sample_count(), dist.total(), false, Policy());
 160    }
 161
 162    template <class RealType, class Policy, class U>
 163    inline RealType cdf(const hypergeometric_distribution<RealType, Policy>& dist, const U& x)
 164    {
 165       BOOST_MATH_STD_USING
 166       static const char* function = "boost::math::cdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 167       RealType r = static_cast<RealType>(x);
 168       unsigned u = itrunc(r, typename policies::normalise<Policy, policies::rounding_error<policies::ignore_error> >::type());
 169       if(u != r)
 170       {
 171          return boost::math::policies::raise_domain_error<RealType>(
 172             function, "Random variable out of range: must be an integer but got %1%", r, Policy());
 173       }
 174       return cdf(dist, u);
 175    }
 176
 177    template <class RealType, class Policy>
 178    inline RealType cdf(const complemented2_type<hypergeometric_distribution<RealType, Policy>, unsigned>& c)
 179    {
 180       static const char* function = "boost::math::cdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 181       RealType result = 0;
 182       if(!c.dist.check_params(function, &result))
 183          return result;
 184       if(!c.dist.check_x(c.param, function, &result))
 185          return result;
 186
 187       return boost::math::detail::hypergeometric_cdf<RealType>(
 188          c.param, c.dist.defective(), c.dist.sample_count(), c.dist.total(), true, Policy());
 189    }
 190
 191    template <class RealType, class Policy, class U>
 192    inline RealType cdf(const complemented2_type<hypergeometric_distribution<RealType, Policy>, U>& c)
 193    {
 194       BOOST_MATH_STD_USING
 195       static const char* function = "boost::math::cdf(const hypergeometric_distribution<%1%>&, const %1%&)";
 196       RealType r = static_cast<RealType>(c.param);
 197       unsigned u = itrunc(r, typename policies::normalise<Policy, policies::rounding_error<policies::ignore_error> >::type());
 198       if(u != r)
 199       {
 200          return boost::math::policies::raise_domain_error<RealType>(
 201             function, "Random variable out of range: must be an integer but got %1%", r, Policy());
 202       }
 203       return cdf(complement(c.dist, u));
 204    }
 205
 206    template <class RealType, class Policy>
 207    inline RealType quantile(const hypergeometric_distribution<RealType, Policy>& dist, const RealType& p)
 208    {
 209       BOOST_MATH_STD_USING // for ADL of std functions
 210
 211          // Checking function argument
 212          RealType result = 0;
 213       const char* function = "boost::math::quantile(const hypergeometric_distribution<%1%>&, %1%)";
 214       if (false == dist.check_params(function, &result)) return result;
 215       if(false == detail::check_probability(function, p, &result, Policy())) return result;
 216
 217       return static_cast<RealType>(detail::hypergeometric_quantile(p, RealType(1 - p), dist.defective(), dist.sample_count(), dist.total(), Policy()));
 218    } // quantile
 219
 220    template <class RealType, class Policy>
 221    inline RealType quantile(const complemented2_type<hypergeometric_distribution<RealType, Policy>, RealType>& c)
 222    {
 223       BOOST_MATH_STD_USING // for ADL of std functions
 224
 225       // Checking function argument
 226       RealType result = 0;
 227       const char* function = "quantile(const complemented2_type<hypergeometric_distribution<%1%>, %1%>&)";
 228       if (false == c.dist.check_params(function, &result)) return result;
 229       if(false == detail::check_probability(function, c.param, &result, Policy())) return result;
 230
 231       return static_cast<RealType>(detail::hypergeometric_quantile(RealType(1 - c.param), c.param, c.dist.defective(), c.dist.sample_count(), c.dist.total(), Policy()));
 232    } // quantile
 233
 234    template <class RealType, class Policy>
 235    inline RealType mean(const hypergeometric_distribution<RealType, Policy>& dist)
 236    {
 237       return static_cast<RealType>(dist.defective() * dist.sample_count()) / dist.total();
 238    } // RealType mean(const hypergeometric_distribution<RealType, Policy>& dist)
 239
 240    template <class RealType, class Policy>
 241    inline RealType variance(const hypergeometric_distribution<RealType, Policy>& dist)
 242    {
 243       RealType r = static_cast<RealType>(dist.defective());
 244       RealType n = static_cast<RealType>(dist.sample_count());
 245       RealType N = static_cast<RealType>(dist.total());
 246       return n * r  * (N - r) * (N - n) / (N * N * (N - 1));
 247    } // RealType variance(const hypergeometric_distribution<RealType, Policy>& dist)
 248
 249    template <class RealType, class Policy>
 250    inline RealType mode(const hypergeometric_distribution<RealType, Policy>& dist)
 251    {
 252       BOOST_MATH_STD_USING
 253       RealType r = static_cast<RealType>(dist.defective());
 254       RealType n = static_cast<RealType>(dist.sample_count());
 255       RealType N = static_cast<RealType>(dist.total());
 256       return floor((r + 1) * (n + 1) / (N + 2));
 257    }
 258
 259    template <class RealType, class Policy>
 260    inline RealType skewness(const hypergeometric_distribution<RealType, Policy>& dist)
 261    {
 262       BOOST_MATH_STD_USING
 263       RealType r = static_cast<RealType>(dist.defective());
 264       RealType n = static_cast<RealType>(dist.sample_count());
 265       RealType N = static_cast<RealType>(dist.total());
 266       return (N - 2 * r) * sqrt(N - 1) * (N - 2 * n) / (sqrt(n * r * (N - r) * (N - n)) * (N - 2));
 267    } // RealType skewness(const hypergeometric_distribution<RealType, Policy>& dist)
 268
 269    template <class RealType, class Policy>
 270    inline RealType kurtosis_excess(const hypergeometric_distribution<RealType, Policy>& dist)
 271    {
 272       RealType r = static_cast<RealType>(dist.defective());
 273       RealType n = static_cast<RealType>(dist.sample_count());
 274       RealType N = static_cast<RealType>(dist.total());
 275       RealType t1 = N * N * (N - 1) / (r * (N - 2) * (N - 3) * (N - r));
 276       RealType t2 = (N * (N + 1) - 6 * N * (N - r)) / (n * (N - n))
 277          + 3 * r * (N - r) * (N + 6) / (N * N) - 6;
 278       return t1 * t2;
 279    } // RealType kurtosis_excess(const hypergeometric_distribution<RealType, Policy>& dist)
 280
 281    template <class RealType, class Policy>
 282    inline RealType kurtosis(const hypergeometric_distribution<RealType, Policy>& dist)
 283    {
 284       return kurtosis_excess(dist) + 3;
 285    } // RealType kurtosis_excess(const hypergeometric_distribution<RealType, Policy>& dist)
 286 }} // namespaces
 287
 288 // This include must be at the end, *after* the accessors
 289 // for this distribution have been defined, in order to
 290 // keep compilers that support two-phase lookup happy.
 291 #include <boost/math/distributions/detail/derived_accessors.hpp>
 292
 293 #endif // include guard