laplace.hpp 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489
  1. // Copyright Thijs van den Berg, 2008.
  2. // Copyright John Maddock 2008.
  3. // Copyright Paul A. Bristow 2008, 2014.
  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. // This module implements the Laplace distribution.
  8. // Weisstein, Eric W. "Laplace Distribution." From MathWorld--A Wolfram Web Resource.
  9. // http://mathworld.wolfram.com/LaplaceDistribution.html
  10. // http://en.wikipedia.org/wiki/Laplace_distribution
  11. //
  12. // Abramowitz and Stegun 1972, p 930
  13. // http://www.math.sfu.ca/~cbm/aands/page_930.htm
  14. #ifndef BOOST_STATS_LAPLACE_HPP
  15. #define BOOST_STATS_LAPLACE_HPP
  16. #include <boost/math/special_functions/log1p.hpp>
  17. #include <boost/math/distributions/detail/common_error_handling.hpp>
  18. #include <boost/math/distributions/complement.hpp>
  19. #include <boost/math/constants/constants.hpp>
  20. #include <limits>
  21. namespace boost{ namespace math{
  22. #ifdef _MSC_VER
  23. # pragma warning(push)
  24. # pragma warning(disable:4127) // conditional expression is constant
  25. #endif
  26. template <class RealType = double, class Policy = policies::policy<> >
  27. class laplace_distribution
  28. {
  29. public:
  30. // ----------------------------------
  31. // public Types
  32. // ----------------------------------
  33. using value_type = RealType;
  34. using policy_type = Policy;
  35. // ----------------------------------
  36. // Constructor(s)
  37. // ----------------------------------
  38. explicit laplace_distribution(RealType l_location = 0, RealType l_scale = 1)
  39. : m_location(l_location), m_scale(l_scale)
  40. {
  41. RealType result;
  42. check_parameters("boost::math::laplace_distribution<%1%>::laplace_distribution()", &result);
  43. }
  44. // ----------------------------------
  45. // Public functions
  46. // ----------------------------------
  47. RealType location() const
  48. {
  49. return m_location;
  50. }
  51. RealType scale() const
  52. {
  53. return m_scale;
  54. }
  55. bool check_parameters(const char* function, RealType* result) const
  56. {
  57. if(false == detail::check_scale(function, m_scale, result, Policy())) return false;
  58. if(false == detail::check_location(function, m_location, result, Policy())) return false;
  59. return true;
  60. }
  61. private:
  62. RealType m_location;
  63. RealType m_scale;
  64. }; // class laplace_distribution
  65. //
  66. // Convenient type synonym for double.
  67. using laplace = laplace_distribution<double>;
  68. #ifdef __cpp_deduction_guides
  69. template <class RealType>
  70. laplace_distribution(RealType)->laplace_distribution<typename boost::math::tools::promote_args<RealType>::type>;
  71. template <class RealType>
  72. laplace_distribution(RealType,RealType)->laplace_distribution<typename boost::math::tools::promote_args<RealType>::type>;
  73. #endif
  74. //
  75. // Non-member functions.
  76. template <class RealType, class Policy>
  77. inline std::pair<RealType, RealType> range(const laplace_distribution<RealType, Policy>&)
  78. {
  79. if (std::numeric_limits<RealType>::has_infinity)
  80. { // Can use infinity.
  81. return std::pair<RealType, RealType>(-std::numeric_limits<RealType>::infinity(), std::numeric_limits<RealType>::infinity()); // - to + infinity.
  82. }
  83. else
  84. { // Can only use max_value.
  85. using boost::math::tools::max_value;
  86. return std::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.
  87. }
  88. }
  89. template <class RealType, class Policy>
  90. inline std::pair<RealType, RealType> support(const laplace_distribution<RealType, Policy>&)
  91. {
  92. if (std::numeric_limits<RealType>::has_infinity)
  93. { // Can Use infinity.
  94. return std::pair<RealType, RealType>(-std::numeric_limits<RealType>::infinity(), std::numeric_limits<RealType>::infinity()); // - to + infinity.
  95. }
  96. else
  97. { // Can only use max_value.
  98. using boost::math::tools::max_value;
  99. return std::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.
  100. }
  101. }
  102. template <class RealType, class Policy>
  103. inline RealType pdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
  104. {
  105. BOOST_MATH_STD_USING // for ADL of std functions
  106. // Checking function argument
  107. RealType result = 0;
  108. const char* function = "boost::math::pdf(const laplace_distribution<%1%>&, %1%))";
  109. // Check scale and location.
  110. if (false == dist.check_parameters(function, &result)) return result;
  111. // Special pdf values.
  112. if((boost::math::isinf)(x))
  113. {
  114. return 0; // pdf + and - infinity is zero.
  115. }
  116. if (false == detail::check_x(function, x, &result, Policy())) return result;
  117. // General case
  118. RealType scale( dist.scale() );
  119. RealType location( dist.location() );
  120. RealType exponent = x - location;
  121. if (exponent>0) exponent = -exponent;
  122. exponent /= scale;
  123. result = exp(exponent);
  124. result /= 2 * scale;
  125. return result;
  126. } // pdf
  127. template <class RealType, class Policy>
  128. inline RealType logpdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
  129. {
  130. BOOST_MATH_STD_USING // for ADL of std functions
  131. // Checking function argument
  132. RealType result = -std::numeric_limits<RealType>::infinity();
  133. const char* function = "boost::math::logpdf(const laplace_distribution<%1%>&, %1%))";
  134. // Check scale and location.
  135. if (false == dist.check_parameters(function, &result))
  136. {
  137. return result;
  138. }
  139. // Special pdf values.
  140. if((boost::math::isinf)(x))
  141. {
  142. return result; // pdf + and - infinity is zero so logpdf is -INF
  143. }
  144. if (false == detail::check_x(function, x, &result, Policy()))
  145. {
  146. return result;
  147. }
  148. const RealType mu = dist.scale();
  149. const RealType b = dist.location();
  150. // if b is 0 avoid divde by 0 error
  151. if(abs(b) < std::numeric_limits<RealType>::epsilon())
  152. {
  153. result = log(pdf(dist, x));
  154. }
  155. else
  156. {
  157. // General case
  158. const RealType log2 = boost::math::constants::ln_two<RealType>();
  159. result = -abs(x-mu)/b - log(b) - log2;
  160. }
  161. return result;
  162. } // logpdf
  163. template <class RealType, class Policy>
  164. inline RealType cdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
  165. {
  166. BOOST_MATH_STD_USING // For ADL of std functions.
  167. RealType result = 0;
  168. // Checking function argument.
  169. const char* function = "boost::math::cdf(const laplace_distribution<%1%>&, %1%)";
  170. // Check scale and location.
  171. if (false == dist.check_parameters(function, &result)) return result;
  172. // Special cdf values:
  173. if((boost::math::isinf)(x))
  174. {
  175. if(x < 0) return 0; // -infinity.
  176. return 1; // + infinity.
  177. }
  178. if (false == detail::check_x(function, x, &result, Policy())) return result;
  179. // General cdf values
  180. RealType scale( dist.scale() );
  181. RealType location( dist.location() );
  182. if (x < location)
  183. {
  184. result = exp( (x-location)/scale )/2;
  185. }
  186. else
  187. {
  188. result = 1 - exp( (location-x)/scale )/2;
  189. }
  190. return result;
  191. } // cdf
  192. template <class RealType, class Policy>
  193. inline RealType logcdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
  194. {
  195. BOOST_MATH_STD_USING // For ADL of std functions.
  196. RealType result = 0;
  197. // Checking function argument.
  198. const char* function = "boost::math::logcdf(const laplace_distribution<%1%>&, %1%)";
  199. // Check scale and location.
  200. if (false == dist.check_parameters(function, &result))
  201. {
  202. return result;
  203. }
  204. // Special cdf values:
  205. if((boost::math::isinf)(x))
  206. {
  207. if(x < 0)
  208. {
  209. return 0; // -infinity.
  210. }
  211. return 1; // + infinity.
  212. }
  213. if (false == detail::check_x(function, x, &result, Policy()))
  214. {
  215. return result;
  216. }
  217. // General cdf values
  218. RealType scale( dist.scale() );
  219. RealType location( dist.location() );
  220. if (x < location)
  221. {
  222. result = ((x - location) / scale) - boost::math::constants::ln_two<RealType>();
  223. }
  224. else
  225. {
  226. result = log1p(-exp((location - x) / scale) / 2);
  227. }
  228. return result;
  229. } // logcdf
  230. template <class RealType, class Policy>
  231. inline RealType quantile(const laplace_distribution<RealType, Policy>& dist, const RealType& p)
  232. {
  233. BOOST_MATH_STD_USING // for ADL of std functions.
  234. // Checking function argument
  235. RealType result = 0;
  236. const char* function = "boost::math::quantile(const laplace_distribution<%1%>&, %1%)";
  237. if (false == dist.check_parameters(function, &result)) return result;
  238. if(false == detail::check_probability(function, p, &result, Policy())) return result;
  239. // Extreme values of p:
  240. if(p == 0)
  241. {
  242. result = policies::raise_overflow_error<RealType>(function,
  243. "probability parameter is 0, but must be > 0!", Policy());
  244. return -result; // -inf
  245. }
  246. if(p == 1)
  247. {
  248. result = policies::raise_overflow_error<RealType>(function,
  249. "probability parameter is 1, but must be < 1!", Policy());
  250. return result; // inf
  251. }
  252. // Calculate Quantile
  253. RealType scale( dist.scale() );
  254. RealType location( dist.location() );
  255. if (p - 0.5 < 0.0)
  256. result = location + scale*log( static_cast<RealType>(p*2) );
  257. else
  258. result = location - scale*log( static_cast<RealType>(-p*2 + 2) );
  259. return result;
  260. } // quantile
  261. template <class RealType, class Policy>
  262. inline RealType cdf(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
  263. {
  264. // Calculate complement of cdf.
  265. BOOST_MATH_STD_USING // for ADL of std functions
  266. RealType scale = c.dist.scale();
  267. RealType location = c.dist.location();
  268. RealType x = c.param;
  269. RealType result = 0;
  270. // Checking function argument.
  271. const char* function = "boost::math::cdf(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
  272. // Check scale and location.
  273. if (false == c.dist.check_parameters(function, &result)) return result;
  274. // Special cdf values.
  275. if((boost::math::isinf)(x))
  276. {
  277. if(x < 0) return 1; // cdf complement -infinity is unity.
  278. return 0; // cdf complement +infinity is zero.
  279. }
  280. if(false == detail::check_x(function, x, &result, Policy()))return result;
  281. // Cdf interval value.
  282. if (-x < -location)
  283. {
  284. result = exp( (-x+location)/scale )/2;
  285. }
  286. else
  287. {
  288. result = 1 - exp( (-location+x)/scale )/2;
  289. }
  290. return result;
  291. } // cdf complement
  292. template <class RealType, class Policy>
  293. inline RealType logcdf(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
  294. {
  295. // Calculate complement of logcdf.
  296. BOOST_MATH_STD_USING // for ADL of std functions
  297. RealType scale = c.dist.scale();
  298. RealType location = c.dist.location();
  299. RealType x = c.param;
  300. RealType result = 0;
  301. // Checking function argument.
  302. const char* function = "boost::math::logcdf(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
  303. // Check scale and location.
  304. if (false == c.dist.check_parameters(function, &result)) return result;
  305. // Special cdf values.
  306. if((boost::math::isinf)(x))
  307. {
  308. if(x < 0)
  309. {
  310. return 1; // cdf complement -infinity is unity.
  311. }
  312. return 0; // cdf complement +infinity is zero.
  313. }
  314. if(false == detail::check_x(function, x, &result, Policy()))return result;
  315. // Cdf interval value.
  316. if (-x < -location)
  317. {
  318. result = (-x+location)/scale - boost::math::constants::ln_two<RealType>();
  319. }
  320. else
  321. {
  322. result = log1p(-exp( (-location+x)/scale )/2, Policy());
  323. }
  324. return result;
  325. } // cdf complement
  326. template <class RealType, class Policy>
  327. inline RealType quantile(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
  328. {
  329. BOOST_MATH_STD_USING // for ADL of std functions.
  330. // Calculate quantile.
  331. RealType scale = c.dist.scale();
  332. RealType location = c.dist.location();
  333. RealType q = c.param;
  334. RealType result = 0;
  335. // Checking function argument.
  336. const char* function = "quantile(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
  337. if (false == c.dist.check_parameters(function, &result)) return result;
  338. // Extreme values.
  339. if(q == 0)
  340. {
  341. return std::numeric_limits<RealType>::infinity();
  342. }
  343. if(q == 1)
  344. {
  345. return -std::numeric_limits<RealType>::infinity();
  346. }
  347. if(false == detail::check_probability(function, q, &result, Policy())) return result;
  348. if (0.5 - q < 0.0)
  349. result = location + scale*log( static_cast<RealType>(-q*2 + 2) );
  350. else
  351. result = location - scale*log( static_cast<RealType>(q*2) );
  352. return result;
  353. } // quantile
  354. template <class RealType, class Policy>
  355. inline RealType mean(const laplace_distribution<RealType, Policy>& dist)
  356. {
  357. return dist.location();
  358. }
  359. template <class RealType, class Policy>
  360. inline RealType standard_deviation(const laplace_distribution<RealType, Policy>& dist)
  361. {
  362. return constants::root_two<RealType>() * dist.scale();
  363. }
  364. template <class RealType, class Policy>
  365. inline RealType mode(const laplace_distribution<RealType, Policy>& dist)
  366. {
  367. return dist.location();
  368. }
  369. template <class RealType, class Policy>
  370. inline RealType median(const laplace_distribution<RealType, Policy>& dist)
  371. {
  372. return dist.location();
  373. }
  374. template <class RealType, class Policy>
  375. inline RealType skewness(const laplace_distribution<RealType, Policy>& /*dist*/)
  376. {
  377. return 0;
  378. }
  379. template <class RealType, class Policy>
  380. inline RealType kurtosis(const laplace_distribution<RealType, Policy>& /*dist*/)
  381. {
  382. return 6;
  383. }
  384. template <class RealType, class Policy>
  385. inline RealType kurtosis_excess(const laplace_distribution<RealType, Policy>& /*dist*/)
  386. {
  387. return 3;
  388. }
  389. template <class RealType, class Policy>
  390. inline RealType entropy(const laplace_distribution<RealType, Policy> & dist)
  391. {
  392. using std::log;
  393. return log(2*dist.scale()*constants::e<RealType>());
  394. }
  395. #ifdef _MSC_VER
  396. # pragma warning(pop)
  397. #endif
  398. } // namespace math
  399. } // namespace boost
  400. // This include must be at the end, *after* the accessors
  401. // for this distribution have been defined, in order to
  402. // keep compilers that support two-phase lookup happy.
  403. #include <boost/math/distributions/detail/derived_accessors.hpp>
  404. #endif // BOOST_STATS_LAPLACE_HPP