boost/math/distributions/normal.hpp

   1 //  Copyright John Maddock 2006, 2007.
   2 //  Copyright Paul A. Bristow 2006, 2007.
   3
   4 //  Use, modification and distribution are subject to the
   5 //  Boost Software License, Version 1.0. (See accompanying file
   6 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   7
   8 #ifndef BOOST_STATS_NORMAL_HPP
   9 #define BOOST_STATS_NORMAL_HPP
  10
  11 // http://en.wikipedia.org/wiki/Normal_distribution
  12 // http://www.itl.nist.gov/div898/handbook/eda/section3/eda3661.htm
  13 // Also:
  14 // Weisstein, Eric W. "Normal Distribution."
  15 // From MathWorld--A Wolfram Web Resource.
  16 // http://mathworld.wolfram.com/NormalDistribution.html
  17
  18 #include <boost/math/distributions/fwd.hpp>
  19 #include <boost/math/special_functions/erf.hpp> // for erf/erfc.
  20 #include <boost/math/distributions/complement.hpp>
  21 #include <boost/math/distributions/detail/common_error_handling.hpp>
  22
  23 #include <utility>
  24
  25 namespace boost{ namespace math{
  26
  27 template <class RealType = double, class Policy = policies::policy<> >
  28 class normal_distribution
  29 {
  30 public:
  31    typedef RealType value_type;
  32    typedef Policy policy_type;
  33
  34    normal_distribution(RealType mean = 0, RealType sd = 1)
  35       : m_mean(mean), m_sd(sd)
  36    { // Default is a 'standard' normal distribution N01.
  37      static const char* function = "boost::math::normal_distribution<%1%>::normal_distribution";
  38
  39      RealType result;
  40      detail::check_scale(function, sd, &result, Policy());
  41      detail::check_location(function, mean, &result, Policy());
  42    }
  43
  44    RealType mean()const
  45    { // alias for location.
  46       return m_mean;
  47    }
  48
  49    RealType standard_deviation()const
  50    { // alias for scale.
  51       return m_sd;
  52    }
  53
  54    // Synonyms, provided to allow generic use of find_location and find_scale.
  55    RealType location()const
  56    { // location.
  57       return m_mean;
  58    }
  59    RealType scale()const
  60    { // scale.
  61       return m_sd;
  62    }
  63
  64 private:
  65    //
  66    // Data members:
  67    //
  68    RealType m_mean;  // distribution mean or location.
  69    RealType m_sd;    // distribution standard deviation or scale.
  70 }; // class normal_distribution
  71
  72 typedef normal_distribution<double> normal;
  73
  74 template <class RealType, class Policy>
  75 inline const std::pair<RealType, RealType> range(const normal_distribution<RealType, Policy>& /*dist*/)
  76 { // Range of permissible values for random variable x.
  77    using boost::math::tools::max_value;
  78    return std::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.
  79 }
  80
  81 template <class RealType, class Policy>
  82 inline const std::pair<RealType, RealType> support(const normal_distribution<RealType, Policy>& /*dist*/)
  83 { // Range of supported values for random variable x.
  84    // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
  85
  86    using boost::math::tools::max_value;
  87    return std::pair<RealType, RealType>(-max_value<RealType>(),  max_value<RealType>()); // - to + max value.
  88 }
  89
  90 template <class RealType, class Policy>
  91 inline RealType pdf(const normal_distribution<RealType, Policy>& dist, const RealType& x)
  92 {
  93    BOOST_MATH_STD_USING  // for ADL of std functions
  94
  95    RealType sd = dist.standard_deviation();
  96    RealType mean = dist.mean();
  97
  98    static const char* function = "boost::math::pdf(const normal_distribution<%1%>&, %1%)";
  99    if((boost::math::isinf)(x))
 100    {
 101      return 0; // pdf + and - infinity is zero.
 102    }
 103    // Below produces MSVC 4127 warnings, so the above used instead.
 104    //if(std::numeric_limits<RealType>::has_infinity && abs(x) == std::numeric_limits<RealType>::infinity())
 105    //{ // pdf + and - infinity is zero.
 106    //  return 0;
 107    //}
 108
 109    RealType result;
 110    if(false == detail::check_scale(function, sd, &result, Policy()))
 111    {
 112       return result;
 113    }
 114    if(false == detail::check_location(function, mean, &result, Policy()))
 115    {
 116       return result;
 117    }
 118    if(false == detail::check_x(function, x, &result, Policy()))
 119    {
 120       return result;
 121    }
 122
 123    RealType exponent = x - mean;
 124    exponent *= -exponent;
 125    exponent /= 2 * sd * sd;
 126
 127    result = exp(exponent);
 128    result /= sd * sqrt(2 * constants::pi<RealType>());
 129
 130    return result;
 131 } // pdf
 132
 133 template <class RealType, class Policy>
 134 inline RealType cdf(const normal_distribution<RealType, Policy>& dist, const RealType& x)
 135 {
 136    BOOST_MATH_STD_USING  // for ADL of std functions
 137
 138    RealType sd = dist.standard_deviation();
 139    RealType mean = dist.mean();
 140    static const char* function = "boost::math::cdf(const normal_distribution<%1%>&, %1%)";
 141    RealType result;
 142    if(false == detail::check_scale(function, sd, &result, Policy()))
 143    {
 144       return result;
 145    }
 146    if(false == detail::check_location(function, mean, &result, Policy()))
 147    {
 148       return result;
 149    }
 150    if((boost::math::isinf)(x))
 151    {
 152      if(x < 0) return 0; // -infinity
 153      return 1; // + infinity
 154    }
 155    // These produce MSVC 4127 warnings, so the above used instead.
 156    //if(std::numeric_limits<RealType>::has_infinity && x == std::numeric_limits<RealType>::infinity())
 157    //{ // cdf +infinity is unity.
 158    //  return 1;
 159    //}
 160    //if(std::numeric_limits<RealType>::has_infinity && x == -std::numeric_limits<RealType>::infinity())
 161    //{ // cdf -infinity is zero.
 162    //  return 0;
 163    //}
 164    if(false == detail::check_x(function, x, &result, Policy()))
 165    {
 166      return result;
 167    }
 168    RealType diff = (x - mean) / (sd * constants::root_two<RealType>());
 169    result = boost::math::erfc(-diff, Policy()) / 2;
 170    return result;
 171 } // cdf
 172
 173 template <class RealType, class Policy>
 174 inline RealType quantile(const normal_distribution<RealType, Policy>& dist, const RealType& p)
 175 {
 176    BOOST_MATH_STD_USING  // for ADL of std functions
 177
 178    RealType sd = dist.standard_deviation();
 179    RealType mean = dist.mean();
 180    static const char* function = "boost::math::quantile(const normal_distribution<%1%>&, %1%)";
 181
 182    RealType result;
 183    if(false == detail::check_scale(function, sd, &result, Policy()))
 184       return result;
 185    if(false == detail::check_location(function, mean, &result, Policy()))
 186       return result;
 187    if(false == detail::check_probability(function, p, &result, Policy()))
 188       return result;
 189
 190    result= boost::math::erfc_inv(2 * p, Policy());
 191    result = -result;
 192    result *= sd * constants::root_two<RealType>();
 193    result += mean;
 194    return result;
 195 } // quantile
 196
 197 template <class RealType, class Policy>
 198 inline RealType cdf(const complemented2_type<normal_distribution<RealType, Policy>, RealType>& c)
 199 {
 200    BOOST_MATH_STD_USING  // for ADL of std functions
 201
 202    RealType sd = c.dist.standard_deviation();
 203    RealType mean = c.dist.mean();
 204    RealType x = c.param;
 205    static const char* function = "boost::math::cdf(const complement(normal_distribution<%1%>&), %1%)";
 206
 207    if((boost::math::isinf)(x))
 208    {
 209      if(x < 0) return 1; // cdf complement -infinity is unity.
 210      return 0; // cdf complement +infinity is zero
 211    }
 212    // These produce MSVC 4127 warnings, so the above used instead.
 213    //if(std::numeric_limits<RealType>::has_infinity && x == std::numeric_limits<RealType>::infinity())
 214    //{ // cdf complement +infinity is zero.
 215    //  return 0;
 216    //}
 217    //if(std::numeric_limits<RealType>::has_infinity && x == -std::numeric_limits<RealType>::infinity())
 218    //{ // cdf complement -infinity is unity.
 219    //  return 1;
 220    //}
 221    RealType result;
 222    if(false == detail::check_scale(function, sd, &result, Policy()))
 223       return result;
 224    if(false == detail::check_location(function, mean, &result, Policy()))
 225       return result;
 226    if(false == detail::check_x(function, x, &result, Policy()))
 227       return result;
 228
 229    RealType diff = (x - mean) / (sd * constants::root_two<RealType>());
 230    result = boost::math::erfc(diff, Policy()) / 2;
 231    return result;
 232 } // cdf complement
 233
 234 template <class RealType, class Policy>
 235 inline RealType quantile(const complemented2_type<normal_distribution<RealType, Policy>, RealType>& c)
 236 {
 237    BOOST_MATH_STD_USING  // for ADL of std functions
 238
 239    RealType sd = c.dist.standard_deviation();
 240    RealType mean = c.dist.mean();
 241    static const char* function = "boost::math::quantile(const complement(normal_distribution<%1%>&), %1%)";
 242    RealType result;
 243    if(false == detail::check_scale(function, sd, &result, Policy()))
 244       return result;
 245    if(false == detail::check_location(function, mean, &result, Policy()))
 246       return result;
 247    RealType q = c.param;
 248    if(false == detail::check_probability(function, q, &result, Policy()))
 249       return result;
 250    result = boost::math::erfc_inv(2 * q, Policy());
 251    result *= sd * constants::root_two<RealType>();
 252    result += mean;
 253    return result;
 254 } // quantile
 255
 256 template <class RealType, class Policy>
 257 inline RealType mean(const normal_distribution<RealType, Policy>& dist)
 258 {
 259    return dist.mean();
 260 }
 261
 262 template <class RealType, class Policy>
 263 inline RealType standard_deviation(const normal_distribution<RealType, Policy>& dist)
 264 {
 265    return dist.standard_deviation();
 266 }
 267
 268 template <class RealType, class Policy>
 269 inline RealType mode(const normal_distribution<RealType, Policy>& dist)
 270 {
 271    return dist.mean();
 272 }
 273
 274 template <class RealType, class Policy>
 275 inline RealType median(const normal_distribution<RealType, Policy>& dist)
 276 {
 277    return dist.mean();
 278 }
 279
 280 template <class RealType, class Policy>
 281 inline RealType skewness(const normal_distribution<RealType, Policy>& /*dist*/)
 282 {
 283    return 0;
 284 }
 285
 286 template <class RealType, class Policy>
 287 inline RealType kurtosis(const normal_distribution<RealType, Policy>& /*dist*/)
 288 {
 289    return 3;
 290 }
 291
 292 template <class RealType, class Policy>
 293 inline RealType kurtosis_excess(const normal_distribution<RealType, Policy>& /*dist*/)
 294 {
 295    return 0;
 296 }
 297
 298 } // namespace math
 299 } // namespace boost
 300
 301 // This include must be at the end, *after* the accessors
 302 // for this distribution have been defined, in order to
 303 // keep compilers that support two-phase lookup happy.
 304 #include <boost/math/distributions/detail/derived_accessors.hpp>
 305
 306 #endif // BOOST_STATS_NORMAL_HPP
 307
 308