inverse_gamma_distribution_example.cpp 3.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103
  1. // inverse_gamma_distribution_example.cpp
  2. // Copyright Paul A. Bristow 2010.
  3. // Copyright Thomas Mang 2010.
  4. // Use, modification and distribution are subject to the
  5. // Boost Software License, Version 1.0.
  6. // (See accompanying file LICENSE_1_0.txt
  7. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  8. // Example 1 of using inverse gamma
  9. #include <boost/math/distributions/inverse_gamma.hpp>
  10. using boost::math::inverse_gamma_distribution; // inverse_gamma_distribution.
  11. using boost::math::inverse_gamma;
  12. #include <boost/math/special_functions/gamma.hpp>
  13. using boost::math::tgamma; // Used for naive pdf as a comparison.
  14. #include <boost/math/distributions/gamma.hpp>
  15. using boost::math::inverse_gamma_distribution;
  16. #include <iostream>
  17. using std::cout; using std::endl;
  18. #include <iomanip>
  19. using std::setprecision;
  20. #include <cmath>
  21. using std::sqrt;
  22. int main()
  23. {
  24. cout << "Example using Inverse Gamma distribution. " << endl;
  25. // TODO - awaiting a real example using Bayesian statistics.
  26. #ifdef BOOST_NO_CXX11_NUMERIC_LIMITS
  27. int max_digits10 = 2 + (boost::math::policies::digits<double, boost::math::policies::policy<> >() * 30103UL) / 100000UL;
  28. cout << "BOOST_NO_CXX11_NUMERIC_LIMITS is defined" << endl;
  29. #else
  30. int max_digits10 = std::numeric_limits<double>::max_digits10;
  31. #endif
  32. cout << "Show all potentially significant decimal digits std::numeric_limits<double>::max_digits10 = "
  33. << max_digits10 << endl;
  34. cout.precision(max_digits10); //
  35. double shape = 1.;
  36. double scale = 1.;
  37. double x = 0.5;
  38. // Construction using default RealType double, and default shape and scale..
  39. inverse_gamma_distribution<> my_inverse_gamma(shape, scale); // (alpha, beta)
  40. cout << "my_inverse_gamma.shape() = " << my_inverse_gamma.shape()
  41. << ", scale = "<< my_inverse_gamma.scale() << endl;
  42. cout << "x = " << x << ", pdf = " << pdf(my_inverse_gamma, x)
  43. << ", cdf = " << cdf(my_inverse_gamma, x) << endl;
  44. // Construct using typedef and default shape and scale parameters.
  45. inverse_gamma my_ig;
  46. inverse_gamma my_ig23(2, 3);
  47. cout << "my_inverse_gamma.shape() = " << my_ig23.shape()
  48. << ", scale = "<< my_ig23.scale() << endl;
  49. cout << "x = " << x << ", pdf = " << pdf(my_ig23, x)
  50. << ", cdf = " << cdf(my_ig23, x) << endl;
  51. // Example of providing an 'out of domain' or 'bad' parameter,
  52. // here a shape < 1, for which mean is not defined.
  53. // Try block is essential to catch the exception message.
  54. // (Uses the default policy which is to throw on all errors).
  55. try
  56. {
  57. inverse_gamma if051(0.5, 1);
  58. //inverse_gamma if051(0.5, 1);
  59. cout << "mean(if051) = " << mean(if051) << endl;
  60. }
  61. catch(const std::exception& e)
  62. { // Always useful to include try & catch blocks because default policies
  63. // are to throw exceptions on arguments that cause errors like underflow, overflow.
  64. // Lacking try & catch blocks, the program will abort without a message below,
  65. // which may give some helpful clues as to the cause of the exception.
  66. std::cout <<
  67. "\n""Message from thrown exception was:\n " << e.what() << std::endl;
  68. }
  69. return 0;
  70. } // int main()
  71. /*
  72. Output is:
  73. Example using Inverse Gamma distribution.
  74. std::numeric_limits<double>::max_digits10 = 17
  75. my_inverse_gamma.shape() = 1, scale = 1
  76. x = 0.5, pdf = 0.54134113294645081, cdf = 0.1353352832366127
  77. my_inverse_gamma.shape() = 2, scale = 3
  78. x = 0.5, pdf = 0.17847015671997774, cdf = 0.017351265236664509
  79. Message from thrown exception was:
  80. Error in function boost::math::mean(const inverse_gamma_distribution<double>&): Shape parameter is 0.5, but for a defined mean it must be > 1
  81. */