laplace_example.cpp 6.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169
  1. // laplace_example.cpp
  2. // Copyright Paul A. Bristow 2008, 2010.
  3. // Use, modification and distribution are subject to the
  4. // Boost Software License, Version 1.0.
  5. // (See accompanying file LICENSE_1_0.txt
  6. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  7. // Example of using laplace (& comparing with normal) distribution.
  8. // Note that this file contains Quickbook mark-up as well as code
  9. // and comments, don't change any of the special comment mark-ups!
  10. //[laplace_example1
  11. /*`
  12. First we need some includes to access the laplace & normal distributions
  13. (and some std output of course).
  14. */
  15. #include <boost/math/distributions/laplace.hpp> // for laplace_distribution
  16. using boost::math::laplace; // typedef provides default type is double.
  17. #include <boost/math/distributions/normal.hpp> // for normal_distribution
  18. using boost::math::normal; // typedef provides default type is double.
  19. #include <iostream>
  20. using std::cout; using std::endl; using std::left; using std::showpoint; using std::noshowpoint;
  21. #include <iomanip>
  22. using std::setw; using std::setprecision;
  23. #include <limits>
  24. using std::numeric_limits;
  25. int main()
  26. {
  27. cout << "Example: Laplace distribution." << endl;
  28. try
  29. {
  30. { // Traditional tables and values.
  31. /*`Let's start by printing some traditional tables.
  32. */
  33. double step = 1.; // in z
  34. double range = 4; // min and max z = -range to +range.
  35. //int precision = 17; // traditional tables are only computed to much lower precision.
  36. int precision = 4; // traditional table at much lower precision.
  37. int width = 10; // for use with setw.
  38. // Construct standard laplace & normal distributions l & s
  39. normal s; // (default location or mean = zero, and scale or standard deviation = unity)
  40. cout << "Standard normal distribution, mean or location = "<< s.location()
  41. << ", standard deviation or scale = " << s.scale() << endl;
  42. laplace l; // (default mean = zero, and standard deviation = unity)
  43. cout << "Laplace normal distribution, location = "<< l.location()
  44. << ", scale = " << l.scale() << endl;
  45. /*` First the probability distribution function (pdf).
  46. */
  47. cout << "Probability distribution function values" << endl;
  48. cout << " z PDF normal laplace (difference)" << endl;
  49. cout.precision(5);
  50. for (double z = -range; z < range + step; z += step)
  51. {
  52. cout << left << setprecision(3) << setw(6) << z << " "
  53. << setprecision(precision) << setw(width) << pdf(s, z) << " "
  54. << setprecision(precision) << setw(width) << pdf(l, z)<< " ("
  55. << setprecision(precision) << setw(width) << pdf(l, z) - pdf(s, z) // difference.
  56. << ")" << endl;
  57. }
  58. cout.precision(6); // default
  59. /*`Notice how the laplace is less at z = 1 , but has 'fatter' tails at 2 and 3.
  60. And the area under the normal curve from -[infin] up to z,
  61. the cumulative distribution function (cdf).
  62. */
  63. // For a standard distribution
  64. cout << "Standard location = "<< s.location()
  65. << ", scale = " << s.scale() << endl;
  66. cout << "Integral (area under the curve) from - infinity up to z " << endl;
  67. cout << " z CDF normal laplace (difference)" << endl;
  68. for (double z = -range; z < range + step; z += step)
  69. {
  70. cout << left << setprecision(3) << setw(6) << z << " "
  71. << setprecision(precision) << setw(width) << cdf(s, z) << " "
  72. << setprecision(precision) << setw(width) << cdf(l, z) << " ("
  73. << setprecision(precision) << setw(width) << cdf(l, z) - cdf(s, z) // difference.
  74. << ")" << endl;
  75. }
  76. cout.precision(6); // default
  77. /*`
  78. Pretty-printing a traditional 2-dimensional table is left as an exercise for the student,
  79. but why bother now that the Boost Math Toolkit lets you write
  80. */
  81. double z = 2.;
  82. cout << "Area for gaussian z = " << z << " is " << cdf(s, z) << endl; // to get the area for z.
  83. cout << "Area for laplace z = " << z << " is " << cdf(l, z) << endl; //
  84. /*`
  85. Correspondingly, we can obtain the traditional 'critical' values for significance levels.
  86. For the 95% confidence level, the significance level usually called alpha,
  87. is 0.05 = 1 - 0.95 (for a one-sided test), so we can write
  88. */
  89. cout << "95% of gaussian area has a z below " << quantile(s, 0.95) << endl;
  90. cout << "95% of laplace area has a z below " << quantile(l, 0.95) << endl;
  91. // 95% of area has a z below 1.64485
  92. // 95% of laplace area has a z below 2.30259
  93. /*`and a two-sided test (a comparison between two levels, rather than a one-sided test)
  94. */
  95. cout << "95% of gaussian area has a z between " << quantile(s, 0.975)
  96. << " and " << -quantile(s, 0.975) << endl;
  97. cout << "95% of laplace area has a z between " << quantile(l, 0.975)
  98. << " and " << -quantile(l, 0.975) << endl;
  99. // 95% of area has a z between 1.95996 and -1.95996
  100. // 95% of laplace area has a z between 2.99573 and -2.99573
  101. /*`Notice how much wider z has to be to enclose 95% of the area.
  102. */
  103. }
  104. //] [/[laplace_example1]
  105. }
  106. catch(const std::exception& e)
  107. { // Always useful to include try & catch blocks because default policies
  108. // are to throw exceptions on arguments that cause errors like underflow, overflow.
  109. // Lacking try & catch blocks, the program will abort without a message below,
  110. // which may give some helpful clues as to the cause of the exception.
  111. std::cout <<
  112. "\n""Message from thrown exception was:\n " << e.what() << std::endl;
  113. }
  114. return 0;
  115. } // int main()
  116. /*
  117. Output is:
  118. Example: Laplace distribution.
  119. Standard normal distribution, mean or location = 0, standard deviation or scale = 1
  120. Laplace normal distribution, location = 0, scale = 1
  121. Probability distribution function values
  122. z PDF normal laplace (difference)
  123. -4 0.0001338 0.009158 (0.009024 )
  124. -3 0.004432 0.02489 (0.02046 )
  125. -2 0.05399 0.06767 (0.01368 )
  126. -1 0.242 0.1839 (-0.05803 )
  127. 0 0.3989 0.5 (0.1011 )
  128. 1 0.242 0.1839 (-0.05803 )
  129. 2 0.05399 0.06767 (0.01368 )
  130. 3 0.004432 0.02489 (0.02046 )
  131. 4 0.0001338 0.009158 (0.009024 )
  132. Standard location = 0, scale = 1
  133. Integral (area under the curve) from - infinity up to z
  134. z CDF normal laplace (difference)
  135. -4 3.167e-005 0.009158 (0.009126 )
  136. -3 0.00135 0.02489 (0.02354 )
  137. -2 0.02275 0.06767 (0.04492 )
  138. -1 0.1587 0.1839 (0.02528 )
  139. 0 0.5 0.5 (0 )
  140. 1 0.8413 0.8161 (-0.02528 )
  141. 2 0.9772 0.9323 (-0.04492 )
  142. 3 0.9987 0.9751 (-0.02354 )
  143. 4 1 0.9908 (-0.009126 )
  144. Area for gaussian z = 2 is 0.97725
  145. Area for laplace z = 2 is 0.932332
  146. 95% of gaussian area has a z below 1.64485
  147. 95% of laplace area has a z below 2.30259
  148. 95% of gaussian area has a z between 1.95996 and -1.95996
  149. 95% of laplace area has a z between 2.99573 and -2.99573
  150. */