boost/random/negative_binomial_distribution.hpp

0001 /* boost random/negative_binomial_distribution.hpp header file
0002  *
0003  * Copyright Steven Watanabe 2010
0004  * Distributed under the Boost Software License, Version 1.0. (See
0005  * accompanying file LICENSE_1_0.txt or copy at
0006  * http://www.boost.org/LICENSE_1_0.txt)
0007  *
0008  * See http://www.boost.org for most recent version including documentation.
0009  *
0010  * $Id$
0011  */
0012
0013 #ifndef BOOST_RANDOM_NEGATIVE_BINOMIAL_DISTRIBUTION_HPP_INCLUDED
0014 #define BOOST_RANDOM_NEGATIVE_BINOMIAL_DISTRIBUTION_HPP_INCLUDED
0015
0016 #include <iosfwd>
0017
0018 #include <boost/limits.hpp>
0019 #include <boost/random/detail/config.hpp>
0020 #include <boost/random/gamma_distribution.hpp>
0021 #include <boost/random/poisson_distribution.hpp>
0022
0023 namespace boost {
0024 namespace random {
0025
0026 /**
0027  * The negative binomial distribution is an integer valued
0028  * distribution with two parameters, @c k and @c p.  The
0029  * distribution produces non-negative values.
0030  *
0031  * The distribution function is
0032  * \f$\displaystyle P(i) = {k+i-1\choose i}p^k(1-p)^i\f$.
0033  *
0034  * This implementation uses a gamma-poisson mixture.
0035  */
0036 template<class IntType = int, class RealType = double>
0037 class negative_binomial_distribution {
0038 public:
0039     typedef IntType result_type;
0040     typedef RealType input_type;
0041
0042     class param_type {
0043     public:
0044         typedef negative_binomial_distribution distribution_type;
0045         /**
0046          * Construct a param_type object.  @c k and @c p
0047          * are the parameters of the distribution.
0048          *
0049          * Requires: k >=0 && 0 <= p <= 1
0050          */
0051         explicit param_type(IntType k_arg = 1, RealType p_arg = RealType (0.5))
0052           : _k(k_arg), _p(p_arg)
0053         {}
0054         /** Returns the @c k parameter of the distribution. */
0055         IntType k() const { return _k; }
0056         /** Returns the @c p parameter of the distribution. */
0057         RealType p() const { return _p; }
0058 #ifndef BOOST_RANDOM_NO_STREAM_OPERATORS
0059         /** Writes the parameters of the distribution to a @c std::ostream. */
0060         template<class CharT, class Traits>
0061         friend std::basic_ostream<CharT,Traits>&
0062         operator<<(std::basic_ostream<CharT,Traits>& os,
0063                    const param_type& parm)
0064         {
0065             os << parm._p << " " << parm._k;
0066             return os;
0067         }
0068
0069         /** Reads the parameters of the distribution from a @c std::istream. */
0070         template<class CharT, class Traits>
0071         friend std::basic_istream<CharT,Traits>&
0072         operator>>(std::basic_istream<CharT,Traits>& is, param_type& parm)
0073         {
0074             is >> parm._p >> std::ws >> parm._k;
0075             return is;
0076         }
0077 #endif
0078         /** Returns true if the parameters have the same values. */
0079         friend bool operator==(const param_type& lhs, const param_type& rhs)
0080         {
0081             return lhs._k == rhs._k && lhs._p == rhs._p;
0082         }
0083         /** Returns true if the parameters have different values. */
0084         friend bool operator!=(const param_type& lhs, const param_type& rhs)
0085         {
0086             return !(lhs == rhs);
0087         }
0088     private:
0089         IntType _k;
0090         RealType _p;
0091     };
0092
0093     /**
0094      * Construct a @c negative_binomial_distribution object. @c k and @c p
0095      * are the parameters of the distribution.
0096      *
0097      * Requires: k >=0 && 0 <= p <= 1
0098      */
0099     explicit negative_binomial_distribution(IntType k_arg = 1,
0100                                             RealType p_arg = RealType(0.5))
0101       : _k(k_arg), _p(p_arg)
0102     {}
0103
0104     /**
0105      * Construct an @c negative_binomial_distribution object from the
0106      * parameters.
0107      */
0108     explicit negative_binomial_distribution(const param_type& parm)
0109       : _k(parm.k()), _p(parm.p())
0110     {}
0111
0112     /**
0113      * Returns a random variate distributed according to the
0114      * negative binomial distribution.
0115      */
0116     template<class URNG>
0117     IntType operator()(URNG& urng) const
0118     {
0119         gamma_distribution<RealType> gamma(_k, (1-_p)/_p);
0120         poisson_distribution<IntType, RealType> poisson(gamma(urng));
0121         return poisson(urng);
0122     }
0123
0124     /**
0125      * Returns a random variate distributed according to the negative
0126      * binomial distribution with parameters specified by @c param.
0127      */
0128     template<class URNG>
0129     IntType operator()(URNG& urng, const param_type& parm) const
0130     {
0131         return negative_binomial_distribution(parm)(urng);
0132     }
0133
0134     /** Returns the @c k parameter of the distribution. */
0135     IntType k() const { return _k; }
0136     /** Returns the @c p parameter of the distribution. */
0137     RealType p() const { return _p; }
0138
0139     /** Returns the smallest value that the distribution can produce. */
0140     IntType min BOOST_PREVENT_MACRO_SUBSTITUTION() const { return 0; }
0141     /** Returns the largest value that the distribution can produce. */
0142     IntType max BOOST_PREVENT_MACRO_SUBSTITUTION() const
0143     { return (std::numeric_limits<IntType>::max)(); }
0144
0145     /** Returns the parameters of the distribution. */
0146     param_type param() const { return param_type(_k, _p); }
0147     /** Sets parameters of the distribution. */
0148     void param(const param_type& parm)
0149     {
0150         _k = parm.k();
0151         _p = parm.p();
0152     }
0153
0154     /**
0155      * Effects: Subsequent uses of the distribution do not depend
0156      * on values produced by any engine prior to invoking reset.
0157      */
0158     void reset() { }
0159
0160 #ifndef BOOST_RANDOM_NO_STREAM_OPERATORS
0161     /** Writes the parameters of the distribution to a @c std::ostream. */
0162     template<class CharT, class Traits>
0163     friend std::basic_ostream<CharT,Traits>&
0164     operator<<(std::basic_ostream<CharT,Traits>& os,
0165                const negative_binomial_distribution& bd)
0166     {
0167         os << bd.param();
0168         return os;
0169     }
0170
0171     /** Reads the parameters of the distribution from a @c std::istream. */
0172     template<class CharT, class Traits>
0173     friend std::basic_istream<CharT,Traits>&
0174     operator>>(std::basic_istream<CharT,Traits>& is,
0175                negative_binomial_distribution& bd)
0176     {
0177         bd.read(is);
0178         return is;
0179     }
0180 #endif
0181
0182     /** Returns true if the two distributions will produce the same
0183         sequence of values, given equal generators. */
0184     friend bool operator==(const negative_binomial_distribution& lhs,
0185                            const negative_binomial_distribution& rhs)
0186     {
0187         return lhs._k == rhs._k && lhs._p == rhs._p;
0188     }
0189     /** Returns true if the two distributions could produce different
0190         sequences of values, given equal generators. */
0191     friend bool operator!=(const negative_binomial_distribution& lhs,
0192                            const negative_binomial_distribution& rhs)
0193     {
0194         return !(lhs == rhs);
0195     }
0196
0197 private:
0198
0199     /// @cond \show_private
0200
0201     template<class CharT, class Traits>
0202     void read(std::basic_istream<CharT, Traits>& is) {
0203         param_type parm;
0204         if(is >> parm) {
0205             param(parm);
0206         }
0207     }
0208
0209     // parameters
0210     IntType _k;
0211     RealType _p;
0212
0213     /// @endcond
0214 };
0215
0216 }
0217
0218 }
0219
0220 #endif