123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802 |
- // 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/tools/test.hpp>
- #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.
- #define BOOST_TEST_MAIN
- #include <boost/test/unit_test.hpp> // for test_main
- #include <boost/test/tools/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE_FRACTION
- #include "test_out_of_range.hpp"
- #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_MATH_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_MATH_CHECK_THROW(
- pdf( // Negative success_fraction!
- geometric_distribution<RealType>(static_cast<RealType>(-0.25)),
- static_cast<RealType>(0)), std::domain_error);
- BOOST_MATH_CHECK_THROW(
- pdf( // Success_fraction > 1!
- geometric_distribution<RealType>(static_cast<RealType>(1.25)),
- static_cast<RealType>(0)),
- std::domain_error);
- BOOST_MATH_CHECK_THROW(
- pdf( // Negative k argument !
- geometric_distribution<RealType>(static_cast<RealType>(0.25)),
- static_cast<RealType>(-1)),
- std::domain_error);
- //BOOST_MATH_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_MATH_CHECK_THROW(
- cdf( // Negative k argument !
- geometric_distribution<RealType>(static_cast<RealType>(0.25)),
- static_cast<RealType>(-1)),
- std::domain_error);
- BOOST_MATH_CHECK_THROW(
- cdf( // Negative success_fraction!
- geometric_distribution<RealType>(static_cast<RealType>(-0.25)),
- static_cast<RealType>(0)), std::domain_error);
- BOOST_MATH_CHECK_THROW(
- cdf( // Success_fraction > 1!
- geometric_distribution<RealType>(static_cast<RealType>(1.25)),
- static_cast<RealType>(0)), std::domain_error);
- BOOST_MATH_CHECK_THROW(
- quantile( // Negative success_fraction!
- geometric_distribution<RealType>(static_cast<RealType>(-0.25)),
- static_cast<RealType>(0)), std::domain_error);
- BOOST_MATH_CHECK_THROW(
- quantile( // Success_fraction > 1!
- geometric_distribution<RealType>(static_cast<RealType>(1.25)),
- static_cast<RealType>(0)), std::domain_error);
- check_out_of_range<geometric_distribution<RealType> >(0.5);
- // 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.
- BOOST_AUTO_TEST_CASE( test_main )
- {
- // 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::endl;
- #endif
-
- } // BOOST_AUTO_TEST_CASE( test_main )
- /*
- */
|