]>
git.proxmox.com Git - ceph.git/blob - ceph/src/boost/libs/math/test/test_beta_dist.cpp
249c646ef415436375a914c6ab8432a61758ad68
3 // Copyright John Maddock 2006.
4 // Copyright Paul A. Bristow 2007, 2009, 2010, 2012.
6 // Use, modification and distribution are subject to the
7 // Boost Software License, Version 1.0.
8 // (See accompanying file LICENSE_1_0.txt
9 // or copy at http://www.boost.org/LICENSE_1_0.txt)
11 // Basic sanity tests for the beta Distribution.
13 // http://members.aol.com/iandjmsmith/BETAEX.HTM beta distribution calculator
14 // Appreas to be a 64-bit calculator showing 17 decimal digit (last is noisy).
15 // Similar to mathCAD?
17 // http://www.nuhertz.com/statmat/distributions.html#Beta
18 // Pretty graphs and explanations for most distributions.
20 // http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp
21 // provided 40 decimal digits accuracy incomplete beta aka beta regularized == cdf
23 // http://www.ausvet.com.au/pprev/content.php?page=PPscript
24 // mode 0.75 5/95% 0.9 alpha 7.39 beta 3.13
25 // http://www.epi.ucdavis.edu/diagnostictests/betabuster.html
26 // Beta Buster also calculates alpha and beta from mode & percentile estimates.
27 // This is NOT (yet) implemented.
30 # pragma warning(disable: 4127) // conditional expression is constant.
31 # pragma warning (disable : 4996) // POSIX name for this item is deprecated.
32 # pragma warning (disable : 4224) // nonstandard extension used : formal parameter 'arg' was previously defined as a type.
35 #include <boost/math/concepts/real_concept.hpp> // for real_concept
36 using ::boost::math::concepts::real_concept
;
37 #include <boost/math/tools/test.hpp>
39 #include <boost/math/distributions/beta.hpp> // for beta_distribution
40 using boost::math::beta_distribution
;
41 using boost::math::beta
;
43 #define BOOST_TEST_MAIN
44 #include <boost/test/unit_test.hpp> // for test_main
45 #include <boost/test/tools/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE_FRACTION
47 #include "test_out_of_range.hpp"
53 using std::numeric_limits
;
55 template <class RealType
>
57 RealType a
, // alpha a
59 RealType x
, // Probability
60 RealType P
, // CDF of beta(a, b)
61 RealType Q
, // Complement of CDF
62 RealType tol
) // Test tolerance.
64 boost::math::beta_distribution
<RealType
> abeta(a
, b
);
65 BOOST_CHECK_CLOSE_FRACTION(cdf(abeta
, x
), P
, tol
);
66 if((P
< 0.99) && (Q
< 0.99))
67 { // We can only check this if P is not too close to 1,
68 // so that we can guarantee that Q is free of error,
69 // (and similarly for Q)
70 BOOST_CHECK_CLOSE_FRACTION(
71 cdf(complement(abeta
, x
)), Q
, tol
);
74 BOOST_CHECK_CLOSE_FRACTION(
75 quantile(abeta
, P
), x
, tol
);
79 // Just check quantile is very small:
80 if((std::numeric_limits
<RealType
>::max_exponent
<= std::numeric_limits
<double>::max_exponent
)
81 && (boost::is_floating_point
<RealType
>::value
))
83 // Limit where this is checked: if exponent range is very large we may
84 // run out of iterations in our root finding algorithm.
85 BOOST_CHECK(quantile(abeta
, P
) < boost::math::tools::epsilon
<RealType
>() * 10);
90 BOOST_CHECK_CLOSE_FRACTION(quantile(complement(abeta
, Q
)), x
, tol
);
93 { // Just check quantile is very small:
94 if((std::numeric_limits
<RealType
>::max_exponent
<= std::numeric_limits
<double>::max_exponent
) && (boost::is_floating_point
<RealType
>::value
))
95 { // Limit where this is checked: if exponent range is very large we may
96 // run out of iterations in our root finding algorithm.
97 BOOST_CHECK(quantile(complement(abeta
, Q
)) < boost::math::tools::epsilon
<RealType
>() * 10);
100 // Estimate alpha & beta from mean and variance:
102 BOOST_CHECK_CLOSE_FRACTION(
103 beta_distribution
<RealType
>::find_alpha(mean(abeta
), variance(abeta
)),
105 BOOST_CHECK_CLOSE_FRACTION(
106 beta_distribution
<RealType
>::find_beta(mean(abeta
), variance(abeta
)),
109 // Estimate sample alpha and beta from others:
110 BOOST_CHECK_CLOSE_FRACTION(
111 beta_distribution
<RealType
>::find_alpha(abeta
.beta(), x
, P
),
113 BOOST_CHECK_CLOSE_FRACTION(
114 beta_distribution
<RealType
>::find_beta(abeta
.alpha(), x
, P
),
116 } // if((P < 0.99) && (Q < 0.99)
118 } // template <class RealType> void test_spot
120 template <class RealType
> // Any floating-point type RealType.
121 void test_spots(RealType
)
123 // Basic sanity checks with 'known good' values.
124 // MathCAD test data is to double precision only,
125 // so set tolerance to 100 eps expressed as a fraction, or
126 // 100 eps of type double expressed as a fraction,
127 // whichever is the larger.
129 RealType tolerance
= (std::max
)
130 (boost::math::tools::epsilon
<RealType
>(),
131 static_cast<RealType
>(std::numeric_limits
<double>::epsilon())); // 0 if real_concept.
133 cout
<< "Boost::math::tools::epsilon = " << boost::math::tools::epsilon
<RealType
>() <<endl
;
134 cout
<< "std::numeric_limits::epsilon = " << std::numeric_limits
<RealType
>::epsilon() <<endl
;
135 cout
<< "epsilon = " << tolerance
;
137 tolerance
*= 100000; // Note: NO * 100 because is fraction, NOT %.
138 cout
<< ", Tolerance = " << tolerance
* 100 << "%." << endl
;
140 // RealType teneps = boost::math::tools::epsilon<RealType>() * 10;
142 // Sources of spot test values:
144 // MathCAD defines dbeta(x, s1, s2) pdf, s1 == alpha, s2 = beta, x = x in Wolfram
145 // pbeta(x, s1, s2) cdf and qbeta(x, s1, s2) inverse of cdf
146 // returns pr(X ,= x) when random variable X
147 // has the beta distribution with parameters s1)alpha) and s2(beta).
148 // s1 > 0 and s2 >0 and 0 < x < 1 (but allows x == 0! and x == 1!)
150 // dbeta(0.5,1,1) = 1
152 using boost::math::beta_distribution
;
153 using ::boost::math::cdf
;
154 using ::boost::math::pdf
;
156 // Tests that should throw:
157 BOOST_MATH_CHECK_THROW(mode(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1))), std::domain_error
);
158 // mode is undefined, and throws domain_error!
160 // BOOST_MATH_CHECK_THROW(median(beta_distribution<RealType>(static_cast<RealType>(1), static_cast<RealType>(1))), std::domain_error);
161 // median is undefined, and throws domain_error!
162 // But now median IS provided via derived accessor as quantile(half).
165 BOOST_MATH_CHECK_THROW( // For various bad arguments.
167 beta_distribution
<RealType
>(static_cast<RealType
>(-1), static_cast<RealType
>(1)), // bad alpha < 0.
168 static_cast<RealType
>(1)), std::domain_error
);
170 BOOST_MATH_CHECK_THROW(
172 beta_distribution
<RealType
>(static_cast<RealType
>(0), static_cast<RealType
>(1)), // bad alpha == 0.
173 static_cast<RealType
>(1)), std::domain_error
);
175 BOOST_MATH_CHECK_THROW(
177 beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(0)), // bad beta == 0.
178 static_cast<RealType
>(1)), std::domain_error
);
180 BOOST_MATH_CHECK_THROW(
182 beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(-1)), // bad beta < 0.
183 static_cast<RealType
>(1)), std::domain_error
);
185 BOOST_MATH_CHECK_THROW(
187 beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)), // bad x < 0.
188 static_cast<RealType
>(-1)), std::domain_error
);
190 BOOST_MATH_CHECK_THROW(
192 beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)), // bad x > 1.
193 static_cast<RealType
>(999)), std::domain_error
);
195 // Some exact pdf values.
197 BOOST_CHECK_EQUAL( // a = b = 1 is uniform distribution.
198 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)),
199 static_cast<RealType
>(1)), // x
200 static_cast<RealType
>(1));
202 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)),
203 static_cast<RealType
>(0)), // x
204 static_cast<RealType
>(1));
205 BOOST_CHECK_CLOSE_FRACTION(
206 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)),
207 static_cast<RealType
>(0.5)), // x
208 static_cast<RealType
>(1),
212 beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)).alpha(),
213 static_cast<RealType
>(1) ); //
216 mean(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1))),
217 static_cast<RealType
>(0.5) ); // Exact one half.
219 BOOST_CHECK_CLOSE_FRACTION(
220 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
221 static_cast<RealType
>(0.5)), // x
222 static_cast<RealType
>(1.5), // Exactly 3/2
225 BOOST_CHECK_CLOSE_FRACTION(
226 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
227 static_cast<RealType
>(0.5)), // x
228 static_cast<RealType
>(1.5), // Exactly 3/2
232 BOOST_CHECK_CLOSE_FRACTION(
233 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
234 static_cast<RealType
>(0.1)), // x
235 static_cast<RealType
>(0.02800000000000000000000000000000000000000L), // Seems exact.
236 // http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp?name=BetaRegularized&ptype=0&z=0.1&a=2&b=2&digits=40
239 BOOST_CHECK_CLOSE_FRACTION(
240 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
241 static_cast<RealType
>(0.0001)), // x
242 static_cast<RealType
>(2.999800000000000000000000000000000000000e-8L),
243 // http://members.aol.com/iandjmsmith/BETAEX.HTM 2.9998000000004
244 // http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp?name=BetaRegularized&ptype=0&z=0.0001&a=2&b=2&digits=40
248 BOOST_CHECK_CLOSE_FRACTION(
249 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
250 static_cast<RealType
>(0.0001)), // x
251 static_cast<RealType
>(0.0005999400000000004L), // http://members.aol.com/iandjmsmith/BETAEX.HTM
252 // Slightly higher tolerance for real concept:
253 (std::numeric_limits
<RealType
>::is_specialized
? 1 : 10) * tolerance
);
256 BOOST_CHECK_CLOSE_FRACTION(
257 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
258 static_cast<RealType
>(0.9999)), // x
259 static_cast<RealType
>(0.999999970002L), // http://members.aol.com/iandjmsmith/BETAEX.HTM
260 // Wolfram 0.9999999700020000000000000000000000000000
263 BOOST_CHECK_CLOSE_FRACTION(
264 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(2)),
265 static_cast<RealType
>(0.9)), // x
266 static_cast<RealType
>(0.9961174629530394895796514664963063381217L),
270 BOOST_CHECK_CLOSE_FRACTION(
271 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
272 static_cast<RealType
>(0.1)), // x
273 static_cast<RealType
>(0.2048327646991334516491978475505189480977L),
277 BOOST_CHECK_CLOSE_FRACTION(
278 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
279 static_cast<RealType
>(0.9)), // x
280 static_cast<RealType
>(0.7951672353008665483508021524494810519023L),
284 BOOST_CHECK_CLOSE_FRACTION(
285 quantile(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
286 static_cast<RealType
>(0.7951672353008665483508021524494810519023L)), // x
287 static_cast<RealType
>(0.9),
291 BOOST_CHECK_CLOSE_FRACTION(
292 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
293 static_cast<RealType
>(0.6)), // x
294 static_cast<RealType
>(0.5640942168489749316118742861695149357858L),
298 BOOST_CHECK_CLOSE_FRACTION(
299 quantile(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
300 static_cast<RealType
>(0.5640942168489749316118742861695149357858L)), // x
301 static_cast<RealType
>(0.6),
306 BOOST_CHECK_CLOSE_FRACTION(
307 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(0.5)),
308 static_cast<RealType
>(0.6)), // x
309 static_cast<RealType
>(0.1778078083562213736802876784474931812329L),
313 BOOST_CHECK_CLOSE_FRACTION(
314 quantile(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(0.5)),
315 static_cast<RealType
>(0.1778078083562213736802876784474931812329L)), // x
316 static_cast<RealType
>(0.6),
320 BOOST_CHECK_CLOSE_FRACTION(
321 cdf(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)),
322 static_cast<RealType
>(0.1)), // x
323 static_cast<RealType
>(0.1), // 0.1000000000000000000000000000000000000000
327 BOOST_CHECK_CLOSE_FRACTION(
328 quantile(beta_distribution
<RealType
>(static_cast<RealType
>(1), static_cast<RealType
>(1)),
329 static_cast<RealType
>(0.1)), // x
330 static_cast<RealType
>(0.1), // 0.1000000000000000000000000000000000000000
334 BOOST_CHECK_CLOSE_FRACTION(
335 cdf(complement(beta_distribution
<RealType
>(static_cast<RealType
>(0.5), static_cast<RealType
>(0.5)),
336 static_cast<RealType
>(0.1))), // complement of x
337 static_cast<RealType
>(0.7951672353008665483508021524494810519023L),
341 BOOST_CHECK_CLOSE_FRACTION(
342 quantile(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
343 static_cast<RealType
>(0.0280000000000000000000000000000000000L)), // x
344 static_cast<RealType
>(0.1),
349 BOOST_CHECK_CLOSE_FRACTION(
350 cdf(complement(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
351 static_cast<RealType
>(0.1))), // x
352 static_cast<RealType
>(0.9720000000000000000000000000000000000000L), // Exact.
356 BOOST_CHECK_CLOSE_FRACTION(
357 pdf(beta_distribution
<RealType
>(static_cast<RealType
>(2), static_cast<RealType
>(2)),
358 static_cast<RealType
>(0.9999)), // x
359 static_cast<RealType
>(0.0005999399999999344L), // http://members.aol.com/iandjmsmith/BETAEX.HTM
360 tolerance
*10); // Note loss of precision calculating 1-p test value.
363 // RealType a, // alpha a
364 // RealType b, // beta b
365 // RealType x, // Probability
366 // RealType P, // CDF of beta(a, b)
367 // RealType Q, // Complement of CDF
368 // RealType tol) // Test tolerance.
370 // These test quantiles and complements, and parameter estimates as well.
371 // Spot values using, for example:
372 // http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp?name=BetaRegularized&ptype=0&z=0.1&a=0.5&b=3&digits=40
375 static_cast<RealType
>(1), // alpha a
376 static_cast<RealType
>(1), // beta b
377 static_cast<RealType
>(0.1), // Probability p
378 static_cast<RealType
>(0.1), // Probability of result (CDF of beta), P
379 static_cast<RealType
>(0.9), // Complement of CDF Q = 1 - P
380 tolerance
); // Test tolerance.
382 static_cast<RealType
>(2), // alpha a
383 static_cast<RealType
>(2), // beta b
384 static_cast<RealType
>(0.1), // Probability p
385 static_cast<RealType
>(0.0280000000000000000000000000000000000L), // Probability of result (CDF of beta), P
386 static_cast<RealType
>(1 - 0.0280000000000000000000000000000000000L), // Complement of CDF Q = 1 - P
387 tolerance
); // Test tolerance.
391 static_cast<RealType
>(2), // alpha a
392 static_cast<RealType
>(2), // beta b
393 static_cast<RealType
>(0.5), // Probability p
394 static_cast<RealType
>(0.5), // Probability of result (CDF of beta), P
395 static_cast<RealType
>(0.5), // Complement of CDF Q = 1 - P
396 tolerance
); // Test tolerance.
399 static_cast<RealType
>(2), // alpha a
400 static_cast<RealType
>(2), // beta b
401 static_cast<RealType
>(0.9), // Probability p
402 static_cast<RealType
>(0.972000000000000), // Probability of result (CDF of beta), P
403 static_cast<RealType
>(1-0.972000000000000), // Complement of CDF Q = 1 - P
404 tolerance
); // Test tolerance.
407 static_cast<RealType
>(2), // alpha a
408 static_cast<RealType
>(2), // beta b
409 static_cast<RealType
>(0.01), // Probability p
410 static_cast<RealType
>(0.0002980000000000000000000000000000000000000L), // Probability of result (CDF of beta), P
411 static_cast<RealType
>(1-0.0002980000000000000000000000000000000000000L), // Complement of CDF Q = 1 - P
412 tolerance
); // Test tolerance.
415 static_cast<RealType
>(2), // alpha a
416 static_cast<RealType
>(2), // beta b
417 static_cast<RealType
>(0.001), // Probability p
418 static_cast<RealType
>(2.998000000000000000000000000000000000000E-6L), // Probability of result (CDF of beta), P
419 static_cast<RealType
>(1-2.998000000000000000000000000000000000000E-6L), // Complement of CDF Q = 1 - P
420 tolerance
); // Test tolerance.
423 static_cast<RealType
>(2), // alpha a
424 static_cast<RealType
>(2), // beta b
425 static_cast<RealType
>(0.0001), // Probability p
426 static_cast<RealType
>(2.999800000000000000000000000000000000000E-8L), // Probability of result (CDF of beta), P
427 static_cast<RealType
>(1-2.999800000000000000000000000000000000000E-8L), // Complement of CDF Q = 1 - P
428 tolerance
); // Test tolerance.
431 static_cast<RealType
>(2), // alpha a
432 static_cast<RealType
>(2), // beta b
433 static_cast<RealType
>(0.99), // Probability p
434 static_cast<RealType
>(0.9997020000000000000000000000000000000000L), // Probability of result (CDF of beta), P
435 static_cast<RealType
>(1-0.9997020000000000000000000000000000000000L), // Complement of CDF Q = 1 - P
436 tolerance
); // Test tolerance.
439 static_cast<RealType
>(0.5), // alpha a
440 static_cast<RealType
>(2), // beta b
441 static_cast<RealType
>(0.5), // Probability p
442 static_cast<RealType
>(0.8838834764831844055010554526310612991060L), // Probability of result (CDF of beta), P
443 static_cast<RealType
>(1-0.8838834764831844055010554526310612991060L), // Complement of CDF Q = 1 - P
444 tolerance
); // Test tolerance.
447 static_cast<RealType
>(0.5), // alpha a
448 static_cast<RealType
>(3.), // beta b
449 static_cast<RealType
>(0.7), // Probability p
450 static_cast<RealType
>(0.9903963064097119299191611355232156905687L), // Probability of result (CDF of beta), P
451 static_cast<RealType
>(1-0.9903963064097119299191611355232156905687L), // Complement of CDF Q = 1 - P
452 tolerance
); // Test tolerance.
455 static_cast<RealType
>(0.5), // alpha a
456 static_cast<RealType
>(3.), // beta b
457 static_cast<RealType
>(0.1), // Probability p
458 static_cast<RealType
>(0.5545844446520295253493059553548880128511L), // Probability of result (CDF of beta), P
459 static_cast<RealType
>(1-0.5545844446520295253493059553548880128511L), // Complement of CDF Q = 1 - P
460 tolerance
); // Test tolerance.
464 // Construction with 'bad' parameters.
465 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(1, -1), std::domain_error
);
466 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(-1, 1), std::domain_error
);
467 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(1, 0), std::domain_error
);
468 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(0, 1), std::domain_error
);
470 beta_distribution
<> dist
;
471 BOOST_MATH_CHECK_THROW(pdf(dist
, -1), std::domain_error
);
472 BOOST_MATH_CHECK_THROW(cdf(dist
, -1), std::domain_error
);
473 BOOST_MATH_CHECK_THROW(cdf(complement(dist
, -1)), std::domain_error
);
474 BOOST_MATH_CHECK_THROW(quantile(dist
, -1), std::domain_error
);
475 BOOST_MATH_CHECK_THROW(quantile(complement(dist
, -1)), std::domain_error
);
476 BOOST_MATH_CHECK_THROW(quantile(dist
, -1), std::domain_error
);
477 BOOST_MATH_CHECK_THROW(quantile(complement(dist
, -1)), std::domain_error
);
479 // No longer allow any parameter to be NaN or inf, so all these tests should throw.
480 if (std::numeric_limits
<RealType
>::has_quiet_NaN
)
482 // Attempt to construct from non-finite should throw.
483 RealType nan
= std::numeric_limits
<RealType
>::quiet_NaN();
484 #ifndef BOOST_NO_EXCEPTIONS
485 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(nan
), std::domain_error
);
486 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(1, nan
), std::domain_error
);
488 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(nan
), std::domain_error
);
489 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(1, nan
), std::domain_error
);
492 // Non-finite parameters should throw.
493 beta_distribution
<RealType
> w(RealType(1));
494 BOOST_MATH_CHECK_THROW(pdf(w
, +nan
), std::domain_error
); // x = NaN
495 BOOST_MATH_CHECK_THROW(cdf(w
, +nan
), std::domain_error
); // x = NaN
496 BOOST_MATH_CHECK_THROW(cdf(complement(w
, +nan
)), std::domain_error
); // x = + nan
497 BOOST_MATH_CHECK_THROW(quantile(w
, +nan
), std::domain_error
); // p = + nan
498 BOOST_MATH_CHECK_THROW(quantile(complement(w
, +nan
)), std::domain_error
); // p = + nan
501 if (std::numeric_limits
<RealType
>::has_infinity
)
503 // Attempt to construct from non-finite should throw.
504 RealType inf
= std::numeric_limits
<RealType
>::infinity();
505 #ifndef BOOST_NO_EXCEPTIONS
506 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(inf
), std::domain_error
);
507 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(1, inf
), std::domain_error
);
509 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(inf
), std::domain_error
);
510 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(1, inf
), std::domain_error
);
513 // Non-finite parameters should throw.
514 beta_distribution
<RealType
> w(RealType(1));
515 #ifndef BOOST_NO_EXCEPTIONS
516 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(inf
), std::domain_error
);
517 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
> w(1, inf
), std::domain_error
);
519 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(inf
), std::domain_error
);
520 BOOST_MATH_CHECK_THROW(beta_distribution
<RealType
>(1, inf
), std::domain_error
);
522 BOOST_MATH_CHECK_THROW(pdf(w
, +inf
), std::domain_error
); // x = inf
523 BOOST_MATH_CHECK_THROW(cdf(w
, +inf
), std::domain_error
); // x = inf
524 BOOST_MATH_CHECK_THROW(cdf(complement(w
, +inf
)), std::domain_error
); // x = + inf
525 BOOST_MATH_CHECK_THROW(quantile(w
, +inf
), std::domain_error
); // p = + inf
526 BOOST_MATH_CHECK_THROW(quantile(complement(w
, +inf
)), std::domain_error
); // p = + inf
529 // Error handling checks:
530 check_out_of_range
<boost::math::beta_distribution
<RealType
> >(1, 1); // (All) valid constructor parameter values.
531 // and range and non-finite.
534 BOOST_MATH_CHECK_THROW(pdf(boost::math::beta_distribution
<RealType
>(0, 1), 0), std::domain_error
);
535 BOOST_MATH_CHECK_THROW(pdf(boost::math::beta_distribution
<RealType
>(-1, 1), 0), std::domain_error
);
536 BOOST_MATH_CHECK_THROW(quantile(boost::math::beta_distribution
<RealType
>(1, 1), -1), std::domain_error
);
537 BOOST_MATH_CHECK_THROW(quantile(boost::math::beta_distribution
<RealType
>(1, 1), 2), std::domain_error
);
540 } // template <class RealType>void test_spots(RealType)
542 BOOST_AUTO_TEST_CASE( test_main
)
544 BOOST_MATH_CONTROL_FP
;
545 // Check that can generate beta distribution using one convenience methods:
546 beta_distribution
<> mybeta11(1., 1.); // Using default RealType double.
548 // boost::math::beta mybeta1(1., 1.); // Using typedef fails.
549 // error C2039: 'beta' : is not a member of 'boost::math'
551 // Basic sanity-check spot values.
553 // Some simple checks using double only.
554 BOOST_CHECK_EQUAL(mybeta11
.alpha(), 1); //
555 BOOST_CHECK_EQUAL(mybeta11
.beta(), 1);
556 BOOST_CHECK_EQUAL(mean(mybeta11
), 0.5); // 1 / (1 + 1) = 1/2 exactly
557 BOOST_MATH_CHECK_THROW(mode(mybeta11
), std::domain_error
);
558 beta_distribution
<> mybeta22(2., 2.); // pdf is dome shape.
559 BOOST_CHECK_EQUAL(mode(mybeta22
), 0.5); // 2-1 / (2+2-2) = 1/2 exactly.
560 beta_distribution
<> mybetaH2(0.5, 2.); //
561 beta_distribution
<> mybetaH3(0.5, 3.); //
563 // Check a few values using double.
564 BOOST_CHECK_EQUAL(pdf(mybeta11
, 1), 1); // is uniform unity over 0 to 1,
565 BOOST_CHECK_EQUAL(pdf(mybeta11
, 0), 1); // including zero and unity.
566 // Although these next three have an exact result, internally they're
567 // *not* treated as special cases, and may be out by a couple of eps:
568 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta11
, 0.5), 1.0, 5*std::numeric_limits
<double>::epsilon());
569 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta11
, 0.0001), 1.0, 5*std::numeric_limits
<double>::epsilon());
570 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta11
, 0.9999), 1.0, 5*std::numeric_limits
<double>::epsilon());
571 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta11
, 0.1), 0.1, 2 * std::numeric_limits
<double>::epsilon());
572 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta11
, 0.5), 0.5, 2 * std::numeric_limits
<double>::epsilon());
573 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta11
, 0.9), 0.9, 2 * std::numeric_limits
<double>::epsilon());
574 BOOST_CHECK_EQUAL(cdf(mybeta11
, 1), 1.); // Exact unity expected.
576 double tol
= std::numeric_limits
<double>::epsilon() * 10;
577 BOOST_CHECK_EQUAL(pdf(mybeta22
, 1), 0); // is dome shape.
578 BOOST_CHECK_EQUAL(pdf(mybeta22
, 0), 0);
579 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta22
, 0.5), 1.5, tol
); // top of dome, expect exactly 3/2.
580 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta22
, 0.0001), 5.9994000000000E-4, tol
);
581 BOOST_CHECK_CLOSE_FRACTION(pdf(mybeta22
, 0.9999), 5.9994000000000E-4, tol
*50);
583 BOOST_CHECK_EQUAL(cdf(mybeta22
, 0.), 0); // cdf is a curved line from 0 to 1.
584 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.1), 0.028000000000000, tol
);
585 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.5), 0.5, tol
);
586 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.9), 0.972000000000000, tol
);
587 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.0001), 2.999800000000000000000000000000000000000E-8, tol
);
588 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.001), 2.998000000000000000000000000000000000000E-6, tol
);
589 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.01), 0.0002980000000000000000000000000000000000000, tol
);
590 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.1), 0.02800000000000000000000000000000000000000, tol
); // exact
591 BOOST_CHECK_CLOSE_FRACTION(cdf(mybeta22
, 0.99), 0.9997020000000000000000000000000000000000, tol
);
593 BOOST_CHECK_EQUAL(cdf(mybeta22
, 1), 1.); // Exact unity expected.
597 BOOST_CHECK_CLOSE_FRACTION(cdf(complement(mybeta22
, 0.9)), 0.028000000000000, tol
);
600 BOOST_CHECK_CLOSE_FRACTION(quantile(mybeta22
, 0.028), 0.1, tol
);
601 BOOST_CHECK_CLOSE_FRACTION(quantile(complement(mybeta22
, 1 - 0.028)), 0.1, tol
);
602 BOOST_CHECK_EQUAL(kurtosis(mybeta11
), 3+ kurtosis_excess(mybeta11
)); // Check kurtosis_excess = kurtosis - 3;
603 BOOST_CHECK_CLOSE_FRACTION(variance(mybeta22
), 0.05, tol
);
604 BOOST_CHECK_CLOSE_FRACTION(mean(mybeta22
), 0.5, tol
);
605 BOOST_CHECK_CLOSE_FRACTION(mode(mybeta22
), 0.5, tol
);
606 BOOST_CHECK_CLOSE_FRACTION(median(mybeta22
), 0.5, sqrt(tol
)); // Theoretical maximum accuracy using Brent is sqrt(epsilon).
608 BOOST_CHECK_CLOSE_FRACTION(skewness(mybeta22
), 0.0, tol
);
609 BOOST_CHECK_CLOSE_FRACTION(kurtosis_excess(mybeta22
), -144.0 / 168, tol
);
610 BOOST_CHECK_CLOSE_FRACTION(skewness(beta_distribution
<>(3, 5)), 0.30983866769659335081434123198259, tol
);
612 BOOST_CHECK_CLOSE_FRACTION(beta_distribution
<double>::find_alpha(mean(mybeta22
), variance(mybeta22
)), mybeta22
.alpha(), tol
); // mean, variance, probability.
613 BOOST_CHECK_CLOSE_FRACTION(beta_distribution
<double>::find_beta(mean(mybeta22
), variance(mybeta22
)), mybeta22
.beta(), tol
);// mean, variance, probability.
615 BOOST_CHECK_CLOSE_FRACTION(mybeta22
.find_alpha(mybeta22
.beta(), 0.8, cdf(mybeta22
, 0.8)), mybeta22
.alpha(), tol
);
616 BOOST_CHECK_CLOSE_FRACTION(mybeta22
.find_beta(mybeta22
.alpha(), 0.8, cdf(mybeta22
, 0.8)), mybeta22
.beta(), tol
);
619 beta_distribution
<real_concept
> rcbeta22(2, 2); // Using RealType real_concept.
620 cout
<< "numeric_limits<real_concept>::is_specialized " << numeric_limits
<real_concept
>::is_specialized
<< endl
;
621 cout
<< "numeric_limits<real_concept>::digits " << numeric_limits
<real_concept
>::digits
<< endl
;
622 cout
<< "numeric_limits<real_concept>::digits10 " << numeric_limits
<real_concept
>::digits10
<< endl
;
623 cout
<< "numeric_limits<real_concept>::epsilon " << numeric_limits
<real_concept
>::epsilon() << endl
;
625 // (Parameter value, arbitrarily zero, only communicates the floating point type).
626 test_spots(0.0F
); // Test float.
627 test_spots(0.0); // Test double.
628 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
629 test_spots(0.0L); // Test long double.
630 #if !BOOST_WORKAROUND(__BORLANDC__, BOOST_TESTED_AT(0x582))
631 test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
634 } // BOOST_AUTO_TEST_CASE( test_main )
640 -Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_beta_dist.exe"
641 Running 1 test case...
642 numeric_limits<real_concept>::is_specialized 0
643 numeric_limits<real_concept>::digits 0
644 numeric_limits<real_concept>::digits10 0
645 numeric_limits<real_concept>::epsilon 0
646 Boost::math::tools::epsilon = 1.19209e-007
647 std::numeric_limits::epsilon = 1.19209e-007
648 epsilon = 1.19209e-007, Tolerance = 0.0119209%.
649 Boost::math::tools::epsilon = 2.22045e-016
650 std::numeric_limits::epsilon = 2.22045e-016
651 epsilon = 2.22045e-016, Tolerance = 2.22045e-011%.
652 Boost::math::tools::epsilon = 2.22045e-016
653 std::numeric_limits::epsilon = 2.22045e-016
654 epsilon = 2.22045e-016, Tolerance = 2.22045e-011%.
655 Boost::math::tools::epsilon = 2.22045e-016
656 std::numeric_limits::epsilon = 0
657 epsilon = 2.22045e-016, Tolerance = 2.22045e-011%.
658 *** No errors detected