369 lines · plain
1// Copyright John Maddock 2006, 2007.2// Copyright Paul A. Bristow 2006, 2007.3// Copyright Matt Borland 2024.4// Use, modification and distribution are subject to the5// Boost Software License, Version 1.0. (See accompanying file6// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)7 8#ifndef BOOST_STATS_NORMAL_HPP9#define BOOST_STATS_NORMAL_HPP10 11// http://en.wikipedia.org/wiki/Normal_distribution12// http://www.itl.nist.gov/div898/handbook/eda/section3/eda3661.htm13// Also:14// Weisstein, Eric W. "Normal Distribution."15// From MathWorld--A Wolfram Web Resource.16// http://mathworld.wolfram.com/NormalDistribution.html17 18#include <boost/math/tools/config.hpp>19#include <boost/math/tools/tuple.hpp>20#include <boost/math/tools/numeric_limits.hpp>21#include <boost/math/tools/promotion.hpp>22#include <boost/math/distributions/fwd.hpp>23#include <boost/math/special_functions/erf.hpp> // for erf/erfc.24#include <boost/math/distributions/complement.hpp>25#include <boost/math/distributions/detail/common_error_handling.hpp>26#include <boost/math/constants/constants.hpp>27#include <boost/math/policies/policy.hpp>28 29namespace boost{ namespace math{30 31template <class RealType = double, class Policy = policies::policy<> >32class normal_distribution33{34public:35 using value_type = RealType;36 using policy_type = Policy;37 38 BOOST_MATH_GPU_ENABLED explicit normal_distribution(RealType l_mean = 0, RealType sd = 1)39 : m_mean(l_mean), m_sd(sd)40 { // Default is a 'standard' normal distribution N01.41 constexpr auto function = "boost::math::normal_distribution<%1%>::normal_distribution";42 43 RealType result;44 detail::check_scale(function, sd, &result, Policy());45 detail::check_location(function, l_mean, &result, Policy());46 }47 48 BOOST_MATH_GPU_ENABLED RealType mean()const49 { // alias for location.50 return m_mean;51 }52 53 BOOST_MATH_GPU_ENABLED RealType standard_deviation()const54 { // alias for scale.55 return m_sd;56 }57 58 // Synonyms, provided to allow generic use of find_location and find_scale.59 BOOST_MATH_GPU_ENABLED RealType location()const60 { // location.61 return m_mean;62 }63 BOOST_MATH_GPU_ENABLED RealType scale()const64 { // scale.65 return m_sd;66 }67 68private:69 //70 // Data members:71 //72 RealType m_mean; // distribution mean or location.73 RealType m_sd; // distribution standard deviation or scale.74}; // class normal_distribution75 76using normal = normal_distribution<double>;77 78//79// Deduction guides, note we don't check the 80// value of __cpp_deduction_guides, just assume81// they work as advertised, even if this is pre-final C++17.82//83#ifdef __cpp_deduction_guides84 85template <class RealType>86normal_distribution(RealType, RealType)->normal_distribution<typename boost::math::tools::promote_args<RealType>::type>;87template <class RealType>88normal_distribution(RealType)->normal_distribution<typename boost::math::tools::promote_args<RealType>::type>;89 90#endif91 92#ifdef _MSC_VER93#pragma warning(push)94#pragma warning(disable:4127)95#endif96 97template <class RealType, class Policy>98BOOST_MATH_GPU_ENABLED inline boost::math::pair<RealType, RealType> range(const normal_distribution<RealType, Policy>& /*dist*/)99{ // Range of permissible values for random variable x.100 BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<RealType>::has_infinity)101 { 102 return boost::math::pair<RealType, RealType>(-boost::math::numeric_limits<RealType>::infinity(), boost::math::numeric_limits<RealType>::infinity()); // - to + infinity.103 }104 else105 { // Can only use max_value.106 using boost::math::tools::max_value;107 return boost::math::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.108 }109}110 111template <class RealType, class Policy>112BOOST_MATH_GPU_ENABLED inline boost::math::pair<RealType, RealType> support(const normal_distribution<RealType, Policy>& /*dist*/)113{ // This is range values for random variable x where cdf rises from 0 to 1, and outside it, the pdf is zero.114 BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<RealType>::has_infinity)115 { 116 return boost::math::pair<RealType, RealType>(-boost::math::numeric_limits<RealType>::infinity(), boost::math::numeric_limits<RealType>::infinity()); // - to + infinity.117 }118 else119 { // Can only use max_value.120 using boost::math::tools::max_value;121 return boost::math::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.122 }123}124 125#ifdef _MSC_VER126#pragma warning(pop)127#endif128 129template <class RealType, class Policy>130BOOST_MATH_GPU_ENABLED inline RealType pdf(const normal_distribution<RealType, Policy>& dist, const RealType& x)131{132 BOOST_MATH_STD_USING // for ADL of std functions133 134 RealType sd = dist.standard_deviation();135 RealType mean = dist.mean();136 137 constexpr auto function = "boost::math::pdf(const normal_distribution<%1%>&, %1%)";138 139 RealType result = 0;140 if(false == detail::check_scale(function, sd, &result, Policy()))141 {142 return result;143 }144 if(false == detail::check_location(function, mean, &result, Policy()))145 {146 return result;147 }148 if((boost::math::isinf)(x))149 {150 return 0; // pdf + and - infinity is zero.151 }152 if(false == detail::check_x(function, x, &result, Policy()))153 {154 return result;155 }156 157 RealType exponent = x - mean;158 exponent *= -exponent;159 exponent /= 2 * sd * sd;160 161 result = exp(exponent);162 result /= sd * sqrt(2 * constants::pi<RealType>());163 164 return result;165} // pdf166 167template <class RealType, class Policy>168BOOST_MATH_GPU_ENABLED inline RealType logpdf(const normal_distribution<RealType, Policy>& dist, const RealType& x)169{170 BOOST_MATH_STD_USING // for ADL of std functions171 172 const RealType sd = dist.standard_deviation();173 const RealType mean = dist.mean();174 175 constexpr auto function = "boost::math::logpdf(const normal_distribution<%1%>&, %1%)";176 177 RealType result = -boost::math::numeric_limits<RealType>::infinity();178 if(false == detail::check_scale(function, sd, &result, Policy()))179 {180 return result;181 }182 if(false == detail::check_location(function, mean, &result, Policy()))183 {184 return result;185 }186 if((boost::math::isinf)(x))187 {188 return result; // pdf + and - infinity is zero so logpdf is -inf189 }190 if(false == detail::check_x(function, x, &result, Policy()))191 {192 return result;193 }194 195 const RealType pi = boost::math::constants::pi<RealType>();196 const RealType half = boost::math::constants::half<RealType>();197 198 result = -log(sd) - half*log(2*pi) - (x-mean)*(x-mean)/(2*sd*sd);199 200 return result;201}202 203template <class RealType, class Policy>204BOOST_MATH_GPU_ENABLED inline RealType cdf(const normal_distribution<RealType, Policy>& dist, const RealType& x)205{206 BOOST_MATH_STD_USING // for ADL of std functions207 208 RealType sd = dist.standard_deviation();209 RealType mean = dist.mean();210 constexpr auto function = "boost::math::cdf(const normal_distribution<%1%>&, %1%)";211 RealType result = 0;212 if(false == detail::check_scale(function, sd, &result, Policy()))213 {214 return result;215 }216 if(false == detail::check_location(function, mean, &result, Policy()))217 {218 return result;219 }220 if((boost::math::isinf)(x))221 {222 if(x < 0) return 0; // -infinity223 return 1; // + infinity224 }225 if(false == detail::check_x(function, x, &result, Policy()))226 {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} // cdf233 234template <class RealType, class Policy>235BOOST_MATH_GPU_ENABLED inline RealType quantile(const normal_distribution<RealType, Policy>& dist, const RealType& p)236{237 BOOST_MATH_STD_USING // for ADL of std functions238 239 RealType sd = dist.standard_deviation();240 RealType mean = dist.mean();241 constexpr auto function = "boost::math::quantile(const normal_distribution<%1%>&, %1%)";242 243 RealType result = 0;244 if(false == detail::check_scale(function, sd, &result, Policy()))245 return result;246 if(false == detail::check_location(function, mean, &result, Policy()))247 return result;248 if(false == detail::check_probability(function, p, &result, Policy()))249 return result;250 251 result= boost::math::erfc_inv(2 * p, Policy());252 result = -result;253 result *= sd * constants::root_two<RealType>();254 result += mean;255 return result;256} // quantile257 258template <class RealType, class Policy>259BOOST_MATH_GPU_ENABLED inline RealType cdf(const complemented2_type<normal_distribution<RealType, Policy>, RealType>& c)260{261 BOOST_MATH_STD_USING // for ADL of std functions262 263 RealType sd = c.dist.standard_deviation();264 RealType mean = c.dist.mean();265 RealType x = c.param;266 constexpr auto function = "boost::math::cdf(const complement(normal_distribution<%1%>&), %1%)";267 268 RealType result = 0;269 if(false == detail::check_scale(function, sd, &result, Policy()))270 return result;271 if(false == detail::check_location(function, mean, &result, Policy()))272 return result;273 if((boost::math::isinf)(x))274 {275 if(x < 0) return 1; // cdf complement -infinity is unity.276 return 0; // cdf complement +infinity is zero277 }278 if(false == detail::check_x(function, x, &result, Policy()))279 return result;280 281 RealType diff = (x - mean) / (sd * constants::root_two<RealType>());282 result = boost::math::erfc(diff, Policy()) / 2;283 return result;284} // cdf complement285 286template <class RealType, class Policy>287BOOST_MATH_GPU_ENABLED inline RealType quantile(const complemented2_type<normal_distribution<RealType, Policy>, RealType>& c)288{289 BOOST_MATH_STD_USING // for ADL of std functions290 291 RealType sd = c.dist.standard_deviation();292 RealType mean = c.dist.mean();293 constexpr auto function = "boost::math::quantile(const complement(normal_distribution<%1%>&), %1%)";294 RealType result = 0;295 if(false == detail::check_scale(function, sd, &result, Policy()))296 return result;297 if(false == detail::check_location(function, mean, &result, Policy()))298 return result;299 RealType q = c.param;300 if(false == detail::check_probability(function, q, &result, Policy()))301 return result;302 result = boost::math::erfc_inv(2 * q, Policy());303 result *= sd * constants::root_two<RealType>();304 result += mean;305 return result;306} // quantile307 308template <class RealType, class Policy>309BOOST_MATH_GPU_ENABLED inline RealType mean(const normal_distribution<RealType, Policy>& dist)310{311 return dist.mean();312}313 314template <class RealType, class Policy>315BOOST_MATH_GPU_ENABLED inline RealType standard_deviation(const normal_distribution<RealType, Policy>& dist)316{317 return dist.standard_deviation();318}319 320template <class RealType, class Policy>321BOOST_MATH_GPU_ENABLED inline RealType mode(const normal_distribution<RealType, Policy>& dist)322{323 return dist.mean();324}325 326template <class RealType, class Policy>327BOOST_MATH_GPU_ENABLED inline RealType median(const normal_distribution<RealType, Policy>& dist)328{329 return dist.mean();330}331 332template <class RealType, class Policy>333BOOST_MATH_GPU_ENABLED inline RealType skewness(const normal_distribution<RealType, Policy>& /*dist*/)334{335 return 0;336}337 338template <class RealType, class Policy>339BOOST_MATH_GPU_ENABLED inline RealType kurtosis(const normal_distribution<RealType, Policy>& /*dist*/)340{341 return 3;342}343 344template <class RealType, class Policy>345BOOST_MATH_GPU_ENABLED inline RealType kurtosis_excess(const normal_distribution<RealType, Policy>& /*dist*/)346{347 return 0;348}349 350template <class RealType, class Policy>351BOOST_MATH_GPU_ENABLED inline RealType entropy(const normal_distribution<RealType, Policy> & dist)352{353 BOOST_MATH_STD_USING354 RealType arg = constants::two_pi<RealType>()*constants::e<RealType>()*dist.standard_deviation()*dist.standard_deviation();355 return log(arg)/2;356}357 358} // namespace math359} // namespace boost360 361// This include must be at the end, *after* the accessors362// for this distribution have been defined, in order to363// keep compilers that support two-phase lookup happy.364#include <boost/math/distributions/detail/derived_accessors.hpp>365 366#endif // BOOST_STATS_NORMAL_HPP367 368 369