| // test_geometric.cpp |
| |
| // Copyright Paul A. Bristow 2010. |
| // Copyright John Maddock 2010. |
| |
| // Use, modification and distribution are subject to the |
| // Boost Software License, Version 1.0. |
| // (See accompanying file LICENSE_1_0.txt |
| // or copy at http://www.boost.org/LICENSE_1_0.txt) |
| |
| // Tests for Geometric Distribution. |
| |
| // Note that these defines must be placed BEFORE #includes. |
| #define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error |
| // because several tests overflow & underflow by design. |
| #define BOOST_MATH_DISCRETE_QUANTILE_POLICY real |
| |
| #ifdef _MSC_VER |
| # pragma warning(disable: 4127) // conditional expression is constant. |
| #endif |
| |
| #if !defined(TEST_FLOAT) && !defined(TEST_DOUBLE) && !defined(TEST_LDOUBLE) && !defined(TEST_REAL_CONCEPT) |
| # define TEST_FLOAT |
| # define TEST_DOUBLE |
| # define TEST_LDOUBLE |
| # define TEST_REAL_CONCEPT |
| #endif |
| |
| #include <boost/math/concepts/real_concept.hpp> // for real_concept |
| using ::boost::math::concepts::real_concept; |
| |
| #include <boost/math/distributions/geometric.hpp> // for geometric_distribution |
| using boost::math::geometric_distribution; |
| using boost::math::geometric; // using typedef for geometric_distribution<double> |
| |
| #include <boost/math/distributions/negative_binomial.hpp> // for some comparisons. |
| |
| #include <boost/test/test_exec_monitor.hpp> // for test_main |
| #include <boost/test/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE_FRACTION |
| |
| #include <iostream> |
| using std::cout; |
| using std::endl; |
| using std::setprecision; |
| using std::showpoint; |
| #include <limits> |
| using std::numeric_limits; |
| |
| template <class RealType> |
| void test_spot( // Test a single spot value against 'known good' values. |
| RealType k, // Number of failures. |
| RealType p, // Probability of success_fraction. |
| RealType P, // CDF probability. |
| RealType Q, // Complement of CDF. |
| RealType tol) // Test tolerance. |
| { |
| boost::math::geometric_distribution<RealType> g(p); |
| BOOST_CHECK_EQUAL(p, g.success_fraction()); |
| BOOST_CHECK_CLOSE_FRACTION(cdf(g, k), P, tol); |
| |
| if((P < 0.99) && (Q < 0.99)) |
| { |
| // We can only check this if P is not too close to 1, |
| // so that we can guarantee that Q is free of error: |
| // |
| BOOST_CHECK_CLOSE_FRACTION( |
| cdf(complement(g, k)), Q, tol); |
| if(k != 0) |
| { |
| BOOST_CHECK_CLOSE_FRACTION( |
| quantile(g, P), k, tol); |
| } |
| else |
| { |
| // Just check quantile is very small: |
| if((std::numeric_limits<RealType>::max_exponent <= std::numeric_limits<double>::max_exponent) |
| && (boost::is_floating_point<RealType>::value)) |
| { |
| // Limit where this is checked: if exponent range is very large we may |
| // run out of iterations in our root finding algorithm. |
| BOOST_CHECK(quantile(g, P) < boost::math::tools::epsilon<RealType>() * 10); |
| } |
| } |
| if(k != 0) |
| { |
| BOOST_CHECK_CLOSE_FRACTION( |
| quantile(complement(g, Q)), k, tol); |
| } |
| else |
| { |
| // Just check quantile is very small: |
| if((std::numeric_limits<RealType>::max_exponent <= std::numeric_limits<double>::max_exponent) |
| && (boost::is_floating_point<RealType>::value)) |
| { |
| // Limit where this is checked: if exponent range is very large we may |
| // run out of iterations in our root finding algorithm. |
| BOOST_CHECK(quantile(complement(g, Q)) < boost::math::tools::epsilon<RealType>() * 10); |
| } |
| } |
| } // if((P < 0.99) && (Q < 0.99)) |
| |
| // Parameter estimation test: estimate success ratio: |
| BOOST_CHECK_CLOSE_FRACTION( |
| geometric_distribution<RealType>::find_lower_bound_on_p( |
| 1+k, P), |
| p, 0.02); // Wide tolerance needed for some tests. |
| // Note we bump up the sample size here, purely for the sake of the test, |
| // internally the function has to adjust the sample size so that we get |
| // the right upper bound, our test undoes this, so we can verify the result. |
| BOOST_CHECK_CLOSE_FRACTION( |
| geometric_distribution<RealType>::find_upper_bound_on_p( |
| 1+k+1, Q), |
| p, 0.02); |
| |
| if(Q < P) |
| { |
| // |
| // We check two things here, that the upper and lower bounds |
| // are the right way around, and that they do actually bracket |
| // the naive estimate of p = successes / (sample size) |
| // |
| BOOST_CHECK( |
| geometric_distribution<RealType>::find_lower_bound_on_p( |
| 1+k, Q) |
| <= |
| geometric_distribution<RealType>::find_upper_bound_on_p( |
| 1+k, Q) |
| ); |
| BOOST_CHECK( |
| geometric_distribution<RealType>::find_lower_bound_on_p( |
| 1+k, Q) |
| <= |
| 1 / (1+k) |
| ); |
| BOOST_CHECK( |
| 1 / (1+k) |
| <= |
| geometric_distribution<RealType>::find_upper_bound_on_p( |
| 1+k, Q) |
| ); |
| } |
| else |
| { |
| // As above but when P is small. |
| BOOST_CHECK( |
| geometric_distribution<RealType>::find_lower_bound_on_p( |
| 1+k, P) |
| <= |
| geometric_distribution<RealType>::find_upper_bound_on_p( |
| 1+k, P) |
| ); |
| BOOST_CHECK( |
| geometric_distribution<RealType>::find_lower_bound_on_p( |
| 1+k, P) |
| <= |
| 1 / (1+k) |
| ); |
| BOOST_CHECK( |
| 1 / (1+k) |
| <= |
| geometric_distribution<RealType>::find_upper_bound_on_p( |
| 1+k, P) |
| ); |
| } |
| |
| // Estimate sample size: |
| BOOST_CHECK_CLOSE_FRACTION( |
| geometric_distribution<RealType>::find_minimum_number_of_trials( |
| k, p, P), |
| 1+k, 0.02); // Can differ 50 to 51 for small p |
| BOOST_CHECK_CLOSE_FRACTION( |
| geometric_distribution<RealType>::find_maximum_number_of_trials( |
| k, p, Q), |
| 1+k, 0.02); |
| |
| } // test_spot |
| |
| template <class RealType> // Any floating-point type RealType. |
| void test_spots(RealType) |
| { |
| // Basic sanity checks. |
| // Most test data is to double precision (17 decimal digits) only, |
| |
| cout << "Floating point Type is " << typeid(RealType).name() << endl; |
| |
| // so set tolerance to 1000 eps expressed as a fraction, |
| // or 1000 eps of type double expressed as a fraction, |
| // whichever is the larger. |
| |
| RealType tolerance = (std::max) |
| (boost::math::tools::epsilon<RealType>(), |
| static_cast<RealType>(std::numeric_limits<double>::epsilon())); |
| tolerance *= 10; // 10 eps |
| |
| cout << "Tolerance = " << tolerance << "." << endl; |
| |
| RealType tol1eps = boost::math::tools::epsilon<RealType>(); // Very tight, suit exact values. |
| //RealType tol2eps = boost::math::tools::epsilon<RealType>() * 2; // Tight, values. |
| RealType tol5eps = boost::math::tools::epsilon<RealType>() * 5; // Wider 5 epsilon. |
| cout << "Tolerance 5 eps = " << tol5eps << "." << endl; |
| |
| |
| // Sources of spot test values are mainly R. |
| |
| using boost::math::geometric_distribution; |
| using boost::math::geometric; |
| using boost::math::cdf; |
| using boost::math::pdf; |
| using boost::math::quantile; |
| using boost::math::complement; |
| |
| BOOST_MATH_STD_USING // for std math functions |
| |
| // Test geometric using cdf spot values R |
| // These test quantiles and complements as well. |
| |
| test_spot( // |
| static_cast<RealType>(2), // Number of failures, k |
| static_cast<RealType>(0.5), // Probability of success as fraction, p |
| static_cast<RealType>(0.875L), // Probability of result (CDF), P |
| static_cast<RealType>(0.125L), // complement CCDF Q = 1 - P |
| tolerance); |
| |
| test_spot( // |
| static_cast<RealType>(0), // Number of failures, k |
| static_cast<RealType>(0.25), // Probability of success as fraction, p |
| static_cast<RealType>(0.25), // Probability of result (CDF), P |
| static_cast<RealType>(0.75), // Q = 1 - P |
| tolerance); |
| |
| test_spot( |
| // R formatC(pgeom(10,0.25), digits=17) [1] "0.95776486396789551" |
| // formatC(pgeom(10,0.25, FALSE), digits=17) [1] "0.042235136032104499" |
| |
| static_cast<RealType>(10), // Number of failures, k |
| static_cast<RealType>(0.25), // Probability of success, p |
| static_cast<RealType>(0.95776486396789551L), // Probability of result (CDF), P |
| static_cast<RealType>(0.042235136032104499L), // Q = 1 - P |
| tolerance); |
| |
| test_spot( // |
| // > R formatC(pgeom(50,0.25, TRUE), digits=17) [1] "0.99999957525875771" |
| // > R formatC(pgeom(50,0.25, FALSE), digits=17) [1] "4.2474124232020353e-07" |
| static_cast<RealType>(50), // Number of failures, k |
| static_cast<RealType>(0.25), // Probability of success, p |
| static_cast<RealType>(0.99999957525875771), // Probability of result (CDF), P |
| static_cast<RealType>(4.2474124232020353e-07), // Q = 1 - P |
| tolerance); |
| /* |
| // This causes failures in find_upper_bound_on_p p is small branch. |
| test_spot( // formatC(pgeom(50,0.01, TRUE), digits=17)[1] "0.40104399353383874" |
| // > formatC(pgeom(50,0.01, FALSE), digits=17) [1] "0.59895600646616121" |
| static_cast<RealType>(50), // Number of failures, k |
| static_cast<RealType>(0.01), // Probability of success, p |
| static_cast<RealType>(0.40104399353383874), // Probability of result (CDF), P |
| static_cast<RealType>(0.59895600646616121), // Q = 1 - P |
| tolerance); |
| */ |
| |
| test_spot( // > formatC(pgeom(50,0.99, TRUE), digits=17) [1] " 1" |
| // formatC(pgeom(50,0.99, FALSE), digits=17) [1] "1.0000000000000364e-102" |
| static_cast<RealType>(50), // Number of failures, k |
| static_cast<RealType>(0.99), // Probability of success, p |
| static_cast<RealType>(1), // Probability of result (CDF), P |
| static_cast<RealType>(1.0000000000000364e-102), // Q = 1 - P |
| tolerance); |
| |
| test_spot( // > formatC(pgeom(1,0.99, TRUE), digits=17) [1] "0.99990000000000001" |
| // > formatC(pgeom(1,0.99, FALSE), digits=17) [1] "0.00010000000000000009" |
| static_cast<RealType>(1), // Number of failures, k |
| static_cast<RealType>(0.99), // Probability of success, p |
| static_cast<RealType>(0.9999), // Probability of result (CDF), P |
| static_cast<RealType>(0.0001), // Q = 1 - P |
| tolerance); |
| |
| if(std::numeric_limits<RealType>::is_specialized) |
| { // An extreme value test that is more accurate than using negative binomial. |
| // Since geometric only uses exp and log functions. |
| test_spot( // > formatC(pgeom(10000, 0.001, TRUE), digits=17) [1] "0.99995487182736897" |
| // > formatC(pgeom(10000,0.001, FALSE), digits=17) [1] "4.5128172631071587e-05" |
| static_cast<RealType>(10000L), // Number of failures, k |
| static_cast<RealType>(0.001L), // Probability of success, p |
| static_cast<RealType>(0.99995487182736897L), // Probability of result (CDF), P |
| static_cast<RealType>(4.5128172631071587e-05L), // Q = 1 - P |
| tolerance); // |
| } // numeric_limit is specialized |
| // End of single spot tests using RealType |
| |
| // Tests on PDF: |
| |
| BOOST_CHECK_CLOSE_FRACTION( //> formatC(dgeom(0,0.5), digits=17)[1] " 0.5" |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(0.0) ), // Number of failures, k is very small but not integral, |
| static_cast<RealType>(0.5), // nearly success probability. |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( //> formatC(dgeom(0,0.5), digits=17)[1] " 0.5" |
| // R treates geom as a discrete distribution. |
| // > formatC(dgeom(1.999999,0.5, FALSE), digits=17) [1] " 0" |
| // Warning message: |
| // In dgeom(1.999999, 0.5, FALSE) : non-integer x = 1.999999 |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(0.0001L) ), // Number of failures, k is very small but not integral, |
| static_cast<RealType>(0.4999653438420768L), // nearly success probability. |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( // > formatC(pgeom(0.0001,0.5, TRUE), digits=17)[1] " 0.5" |
| // > formatC(pgeom(0.0001,0.5, FALSE), digits=17) [1] " 0.5" |
| // R treates geom as a discrete distribution. |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(0.0001L) ), // Number of failures, k is very small but not integral, |
| static_cast<RealType>(0.4999653438420768L), // nearly success probability. |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( // formatC(dgeom(1,0.01), digits=17)[1] "0.0099000000000000008" |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.01L)), |
| static_cast<RealType>(1) ), // Number of failures, k |
| static_cast<RealType>(0.0099000000000000008), // |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( //> formatC(dgeom(1,0.99), digits=17)[1] "0.0099000000000000043" |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.99L)), |
| static_cast<RealType>(1) ), // Number of failures, k |
| static_cast<RealType>(0.00990000000000000043L), // |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( //> > formatC(dgeom(0,0.99), digits=17)[1] "0.98999999999999999" |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.99L)), |
| static_cast<RealType>(0) ), // Number of failures, k |
| static_cast<RealType>(0.98999999999999999L), // |
| tolerance); |
| |
| // p near unity. |
| BOOST_CHECK_CLOSE_FRACTION( // > formatC(dgeom(100,0.99), digits=17)[1] "9.9000000000003448e-201" |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.99L)), |
| static_cast<RealType>(100) ), // Number of failures, k |
| static_cast<RealType>(9.9000000000003448e-201L), // |
| 100 * tolerance); // Note difference |
| |
| // p nearer unity. |
| BOOST_CHECK_CLOSE_FRACTION( // |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.9999)), |
| static_cast<RealType>(10) ), // Number of failures, k |
| // static_cast<double>(9.9989999999889024e-41), // Boost.Math |
| // static_cast<float>(1.00156406e-040) |
| static_cast<RealType>(9.999e-41), // exact from 100 digit calculator. |
| 2e3 * tolerance); // Note bigger tolerance needed. |
| |
| // Moshier Cephes 100 digits calculator says 9.999e-41 |
| //0.9999*pow(1-0.9999,10) |
| // 9.9990000000000000000000000000000000000000000000000000000000000000000000E-41 |
| // 9.998999999988988e-041 |
| // > formatC(dgeom(10, 0.9999), digits=17) [1] "9.9989999999889024e-41" |
| // p * pow(q, k) 9.9989999999889880e-041 |
| // exp(p * k * log1p(-p)) 9.9989999999889024e-041 |
| |
| |
| |
| // 0.9999999999 * pow(1-0.9999999999,10)= 9.9999999990E-101 |
| // > formatC(dgeom(10,0.9999999999), digits=17) [1] "1.0000008273040127e-100" |
| BOOST_CHECK_CLOSE_FRACTION( // |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.9999999999L)), |
| static_cast<RealType>(10) ), // |
| static_cast<RealType>(9.9999999990E-101L), // 1.0000008273040179e-100 |
| 1e9 * tolerance); // Note big tolerance needed. |
| // 1.0000008273040179e-100 Boost.Math |
| // 1.0000008273040127e-100 R |
| // 0.9999999990000004e-100 100 digit calculator 'exact' |
| |
| BOOST_CHECK_CLOSE_FRACTION( // |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.00000000001L)), |
| static_cast<RealType>(10) ), // |
| static_cast<RealType>(9.999999999e-12L), // get 9.9999999989999994e-012 |
| 1 * tolerance); // Note small tolerance needed. |
| |
| |
| BOOST_CHECK_CLOSE_FRACTION( // |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.00000000001L)), |
| static_cast<RealType>(1000) ), // |
| static_cast<RealType>(9.9999999e-12L), // get 9.9999998999999913e-012 |
| tolerance); // Note small tolerance needed. |
| |
| |
| /////////////////////////////////////////////////// |
| BOOST_CHECK_CLOSE_FRACTION( // |
| // > formatC(dgeom(0.0001,0.5, FALSE), digits=17) [1] " 0.5" |
| // R treates geom as a discrete distribution. |
| // But Boost.Math is continuous, so if you want R behaviour, |
| // make number of failures, k into an integer with the floor function. |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(floor(0.0001L)) ), // Number of failures, k is very small but MADE integral, |
| static_cast<RealType>(0.5), // nearly success probability. |
| tolerance); |
| |
| // R switches over at about 1e7 from k = 0, returning 0.5, to k = 1, returning 0.25. |
| // Boost.Math does not do this, even for 0.9999999999999999 |
| // > formatC(pgeom(0.999999,0.5, FALSE), digits=17) [1] " 0.5" |
| // > formatC(pgeom(0.9999999,0.5, FALSE), digits=17) [1] " 0.25" |
| |
| BOOST_CHECK_CLOSE_FRACTION( // > formatC(pgeom(0.0001,0.5, TRUE), digits=17)[1] " 0.5" |
| // > formatC(pgeom(0.0001,0.5, FALSE), digits=17) [1] " 0.5" |
| // R treates geom as a discrete distribution. |
| // But Boost.Math is continuous, so if you want R behaviour, |
| // make number of failures, k into an integer with the floor function. |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(floor(0.9999999999999999L)) ), // Number of failures, k is very small but MADE integral, |
| static_cast<RealType>(0.5), // nearly success probability. |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( // > formatC(pgeom(0.0001,0.5, TRUE), digits=17)[1] " 0.5" |
| // > formatC(pgeom(0.0001,0.5, FALSE), digits=17) [1] " 0.5" |
| // R treates geom as a discrete distribution. |
| // But Boost.Math is continuous, so if you want R behaviour, |
| // make number of failures, k into an integer with the floor function. |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(floor(1. - tolerance)) ), |
| // Number of failures, k is very small but MADE integral, |
| // Need to use tolerance here, |
| // as epsilon is ill-defined for Real concept: |
| // numeric_limits<RealType>::epsilon() 0 |
| static_cast<RealType>(0.5), // nearly success probability. |
| tolerance * 10); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.0001L)), |
| static_cast<RealType>(2)), // k = 2. |
| static_cast<RealType>(9.99800010e-5L), // 'exact ' |
| tolerance); |
| |
| //> formatC(dgeom(2, 0.9999), digits=17) [1] "9.9989999999977806e-09" |
| BOOST_CHECK_CLOSE_FRACTION( |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.9999L)), |
| static_cast<RealType>(2)), // k = 0 |
| static_cast<RealType>(9.999e-9L), // 'exact' |
| 1000*tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.9999L)), |
| static_cast<RealType>(3)), // k = 3 |
| static_cast<RealType>(9.999e-13L), // get |
| 1000*tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| pdf(geometric_distribution<RealType>(static_cast<RealType>(0.9999L)), |
| static_cast<RealType>(5)), // k = 5 |
| static_cast<RealType>(9.999e-21L), // 9.9989999999944947e-021 |
| 1000*tolerance); |
| |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| pdf(geometric_distribution<RealType>( static_cast<RealType>(0.0001L)), |
| static_cast<RealType>(3)), // k = 0. |
| static_cast<RealType>(9.99700029999e-5L), // |
| tolerance); |
| // Tests on cdf: |
| // MathCAD pgeom k, r, p) == failures, successes, probability. |
| |
| BOOST_CHECK_CLOSE_FRACTION(cdf( |
| geometric_distribution<RealType>(static_cast<RealType>(0.5)), // prob 0.5 |
| static_cast<RealType>(0) ), // k = 0 |
| static_cast<RealType>(0.5), // probability =p |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION(cdf(complement( |
| geometric_distribution<RealType>(static_cast<RealType>(0.5)), // |
| static_cast<RealType>(0) )), // k = 0 |
| static_cast<RealType>(0.5), // probability = |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION(cdf( |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), // prob 0.5 |
| static_cast<RealType>(1) ), // k = 0 |
| static_cast<RealType>(0.4375L), // probability =p |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION(cdf(complement( |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), // |
| static_cast<RealType>(1) )), // k = 0 |
| static_cast<RealType>(1-0.4375L), // probability = |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION(cdf(complement( |
| geometric_distribution<RealType>(static_cast<RealType>(0.5)), // |
| static_cast<RealType>(1) )), // k = 0 |
| static_cast<RealType>(0.25), // probability = exact 0.25 |
| tolerance); |
| |
| BOOST_CHECK_CLOSE_FRACTION( // |
| cdf(geometric_distribution<RealType>(static_cast<RealType>(0.5)), |
| static_cast<RealType>(4)), // k =4. |
| static_cast<RealType>(0.96875L), // exact |
| tolerance); |
| |
| |
| // Tests of other functions, mean and other moments ... |
| |
| geometric_distribution<RealType> dist(static_cast<RealType>(0.25)); |
| // mean: |
| BOOST_CHECK_CLOSE_FRACTION( |
| mean(dist), static_cast<RealType>((1 - 0.25) /0.25), tol5eps); |
| BOOST_CHECK_CLOSE_FRACTION( |
| mode(dist), static_cast<RealType>(0), tol1eps); |
| // variance: |
| BOOST_CHECK_CLOSE_FRACTION( |
| variance(dist), static_cast<RealType>((1 - 0.25) / (0.25 * 0.25)), tol5eps); |
| |
| // std deviation: |
| // sqrt(0.75/0.125) |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| standard_deviation(dist), // |
| static_cast<RealType>(sqrt((1.0L - 0.25L) / (0.25L * 0.25L))), // using 100 digit calc |
| tol5eps); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| skewness(dist), // |
| static_cast<RealType>((2-0.25L) /sqrt(0.75L)), |
| // using calculator |
| tol5eps); |
| BOOST_CHECK_CLOSE_FRACTION( |
| kurtosis_excess(dist), // |
| static_cast<RealType>(6 + 0.0625L/0.75L), // |
| tol5eps); |
| // 6.083333333333333 6.166666666666667 |
| BOOST_CHECK_CLOSE_FRACTION( |
| kurtosis(dist), // true |
| static_cast<RealType>(9 + 0.0625L/0.75L), // |
| tol5eps); |
| // hazard: |
| RealType x = static_cast<RealType>(0.125); |
| BOOST_CHECK_CLOSE_FRACTION( |
| hazard(dist, x) |
| , pdf(dist, x) / cdf(complement(dist, x)), tol5eps); |
| // cumulative hazard: |
| BOOST_CHECK_CLOSE_FRACTION( |
| chf(dist, x), -log(cdf(complement(dist, x))), tol5eps); |
| // coefficient_of_variation: |
| BOOST_CHECK_CLOSE_FRACTION( |
| coefficient_of_variation(dist) |
| , standard_deviation(dist) / mean(dist), tol5eps); |
| |
| // Special cases for PDF: |
| BOOST_CHECK_EQUAL( |
| pdf( |
| geometric_distribution<RealType>(static_cast<RealType>(0)), // |
| static_cast<RealType>(0)), |
| static_cast<RealType>(0) ); |
| |
| BOOST_CHECK_EQUAL( |
| pdf( |
| geometric_distribution<RealType>(static_cast<RealType>(0)), |
| static_cast<RealType>(0.0001)), |
| static_cast<RealType>(0) ); |
| |
| BOOST_CHECK_EQUAL( |
| pdf( |
| geometric_distribution<RealType>(static_cast<RealType>(1)), |
| static_cast<RealType>(0.001)), |
| static_cast<RealType>(0) ); |
| |
| BOOST_CHECK_EQUAL( |
| pdf( |
| geometric_distribution<RealType>(static_cast<RealType>(1)), |
| static_cast<RealType>(8)), |
| static_cast<RealType>(0) ); |
| |
| BOOST_CHECK_SMALL( |
| pdf( |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0))- |
| static_cast<RealType>(0.25), |
| 2 * boost::math::tools::epsilon<RealType>() ); // Expect exact, but not quite. |
| // numeric_limits<RealType>::epsilon()); // Not suitable for real concept! |
| |
| // Quantile boundary cases checks: |
| BOOST_CHECK_EQUAL( |
| quantile( // zero P < cdf(0) so should be exactly zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0)), |
| static_cast<RealType>(0)); |
| |
| BOOST_CHECK_EQUAL( |
| quantile( // min P < cdf(0) so should be exactly zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(boost::math::tools::min_value<RealType>())), |
| static_cast<RealType>(0)); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| quantile( // Small P < cdf(0) so should be near zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(boost::math::tools::epsilon<RealType>())), // |
| static_cast<RealType>(0), |
| tol5eps); |
| |
| BOOST_CHECK_CLOSE_FRACTION( |
| quantile( // Small P < cdf(0) so should be exactly zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0.0001)), |
| static_cast<RealType>(0), |
| tolerance); |
| |
| //BOOST_CHECK( // Fails with overflow for real_concept |
| //quantile( // Small P near 1 so k failures should be big. |
| //geometric_distribution<RealType>(static_cast<RealType>(8), static_cast<RealType>(0.25)), |
| //static_cast<RealType>(1 - boost::math::tools::epsilon<RealType>())) <= |
| //static_cast<RealType>(189.56999032670058) // 106.462769 for float |
| //); |
| |
| if(std::numeric_limits<RealType>::has_infinity) |
| { // BOOST_CHECK tests for infinity using std::numeric_limits<>::infinity() |
| // Note that infinity is not implemented for real_concept, so these tests |
| // are only done for types, like built-in float, double.. that have infinity. |
| // Note that these assume that BOOST_MATH_OVERFLOW_ERROR_POLICY is NOT throw_on_error. |
| // #define BOOST_MATH_THROW_ON_OVERFLOW_POLICY == throw_on_error would throw here. |
| // #define BOOST_MAT_DOMAIN_ERROR_POLICY IS defined throw_on_error, |
| // so the throw path of error handling is tested below with BOOST_CHECK_THROW tests. |
| |
| BOOST_CHECK( |
| quantile( // At P == 1 so k failures should be infinite. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(1)) == |
| //static_cast<RealType>(boost::math::tools::infinity<RealType>()) |
| static_cast<RealType>(std::numeric_limits<RealType>::infinity()) ); |
| |
| BOOST_CHECK_EQUAL( |
| quantile( // At 1 == P so should be infinite. |
| geometric_distribution<RealType>( static_cast<RealType>(0.25)), |
| static_cast<RealType>(1)), // |
| std::numeric_limits<RealType>::infinity() ); |
| |
| BOOST_CHECK_EQUAL( |
| quantile(complement( // Q zero 1 so P == 1 < cdf(0) so should be exactly infinity. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0))), |
| std::numeric_limits<RealType>::infinity() ); |
| } // test for infinity using std::numeric_limits<>::infinity() |
| else |
| { // real_concept case, so check it throws rather than returning infinity. |
| BOOST_CHECK_EQUAL( |
| quantile( // At P == 1 so k failures should be infinite. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(1)), |
| boost::math::tools::max_value<RealType>() ); |
| |
| BOOST_CHECK_EQUAL( |
| quantile(complement( // Q zero 1 so P == 1 < cdf(0) so should be exactly infinity. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0))), |
| boost::math::tools::max_value<RealType>()); |
| } // has infinity |
| |
| BOOST_CHECK( // Should work for built-in and real_concept. |
| quantile(complement( // Q near to 1 so P nearly 1, so should be large > 300. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(boost::math::tools::min_value<RealType>()))) |
| >= static_cast<RealType>(300) ); |
| |
| BOOST_CHECK_EQUAL( |
| quantile( // P == 0 < cdf(0) so should be zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(0)), |
| static_cast<RealType>(0)); |
| |
| // Quantile Complement boundary cases: |
| |
| BOOST_CHECK_EQUAL( |
| quantile(complement( // Q = 1 so P = 0 < cdf(0) so should be exactly zero. |
| geometric_distribution<RealType>( static_cast<RealType>(0.25)), |
| static_cast<RealType>(1))), |
| static_cast<RealType>(0) |
| ); |
| |
| BOOST_CHECK_EQUAL( |
| quantile(complement( // Q very near 1 so P == epsilon < cdf(0) so should be exactly zero. |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(1 - boost::math::tools::epsilon<RealType>()))), |
| static_cast<RealType>(0) |
| ); |
| |
| // Check that duff arguments throw domain_error: |
| |
| BOOST_CHECK_THROW( |
| pdf( // Negative success_fraction! |
| geometric_distribution<RealType>(static_cast<RealType>(-0.25)), |
| static_cast<RealType>(0)), std::domain_error); |
| BOOST_CHECK_THROW( |
| pdf( // Success_fraction > 1! |
| geometric_distribution<RealType>(static_cast<RealType>(1.25)), |
| static_cast<RealType>(0)), |
| std::domain_error); |
| BOOST_CHECK_THROW( |
| pdf( // Negative k argument ! |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(-1)), |
| std::domain_error); |
| //BOOST_CHECK_THROW( |
| //pdf( // check limit on k (failures) |
| //geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| //std::numeric_limits<RealType>infinity()), |
| //std::domain_error); |
| BOOST_CHECK_THROW( |
| cdf( // Negative k argument ! |
| geometric_distribution<RealType>(static_cast<RealType>(0.25)), |
| static_cast<RealType>(-1)), |
| std::domain_error); |
| BOOST_CHECK_THROW( |
| cdf( // Negative success_fraction! |
| geometric_distribution<RealType>(static_cast<RealType>(-0.25)), |
| static_cast<RealType>(0)), std::domain_error); |
| BOOST_CHECK_THROW( |
| cdf( // Success_fraction > 1! |
| geometric_distribution<RealType>(static_cast<RealType>(1.25)), |
| static_cast<RealType>(0)), std::domain_error); |
| BOOST_CHECK_THROW( |
| quantile( // Negative success_fraction! |
| geometric_distribution<RealType>(static_cast<RealType>(-0.25)), |
| static_cast<RealType>(0)), std::domain_error); |
| BOOST_CHECK_THROW( |
| quantile( // Success_fraction > 1! |
| geometric_distribution<RealType>(static_cast<RealType>(1.25)), |
| static_cast<RealType>(0)), std::domain_error); |
| // End of check throwing 'duff' out-of-domain values. |
| |
| { // Compare geometric and negative binomial functions. |
| using boost::math::negative_binomial_distribution; |
| using boost::math::geometric_distribution; |
| |
| RealType k = static_cast<RealType>(2.L); |
| RealType alpha = static_cast<RealType>(0.05L); |
| RealType p = static_cast<RealType>(0.5L); |
| |
| BOOST_CHECK_CLOSE_FRACTION( // Successes parameter in negative binomial is 1 for geometric. |
| geometric_distribution<RealType>::find_lower_bound_on_p(k, alpha), |
| negative_binomial_distribution<RealType>::find_lower_bound_on_p(k, static_cast<RealType>(1), alpha), |
| tolerance); |
| BOOST_CHECK_CLOSE_FRACTION( // Successes parameter in negative binomial is 1 for geometric. |
| geometric_distribution<RealType>::find_upper_bound_on_p(k, alpha), |
| negative_binomial_distribution<RealType>::find_upper_bound_on_p(k, static_cast<RealType>(1), alpha), |
| tolerance); |
| BOOST_CHECK_CLOSE_FRACTION( // Should be identical - successes parameter is not used. |
| geometric_distribution<RealType>::find_maximum_number_of_trials(k, p, alpha), |
| negative_binomial_distribution<RealType>::find_maximum_number_of_trials(k, p, alpha), |
| tolerance); |
| } |
| //geometric::find_upper_bound_on_p(k, alpha); |
| return; |
| } // template <class RealType> void test_spots(RealType) // Any floating-point type RealType. |
| |
| int test_main(int, char* []) |
| { |
| // Check that can generate geometric distribution using the two convenience methods: |
| using namespace boost::math; |
| geometric g05d(0.5); // Using typedef - default type is double. |
| geometric_distribution<> g05dd(0.5); // Using default RealType double. |
| |
| // Basic sanity-check spot values. |
| |
| // Test some simple double only examples. |
| geometric_distribution<double> mydist(0.25); |
| // success fraction == 0.25 == 25% or 1 in 4 successes. |
| // Note: double values (matching the distribution definition) avoid the need for any casting. |
| |
| // Check accessor functions return exact values for double at least. |
| BOOST_CHECK_EQUAL(mydist.success_fraction(), static_cast<double>(1./4.)); |
| |
| //cout << numeric_limits<RealType>::epsilon() << endl; |
| |
| // (Parameter value, arbitrarily zero, only communicates the floating point type). |
| #ifdef TEST_FLOAT |
| test_spots(0.0F); // Test float. |
| #endif |
| #ifdef TEST_DOUBLE |
| test_spots(0.0); // Test double. |
| #endif |
| #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS |
| #ifdef TEST_LDOUBLE |
| test_spots(0.0L); // Test long double. |
| #endif |
| #if !BOOST_WORKAROUND(__BORLANDC__, BOOST_TESTED_AT(0x582)) |
| #ifdef TEST_REAL_CONCEPT |
| test_spots(boost::math::concepts::real_concept(0.)); // Test real concept. |
| #endif |
| #endif |
| #else |
| std::cout << "<note>The long double tests have been disabled on this platform " |
| "either because the long double overloads of the usual math functions are " |
| "not available at all, or because they are too inaccurate for these tests " |
| "to pass.</note>" << std::cout; |
| #endif |
| |
| return 0; |
| } // int test_main(int, char* []) |
| |
| /* |
| |
| |
| |
| */ |
