3 // Copyright John Maddock 2008.
5 // Use, modification and distribution are subject to the
6 // Boost Software License, Version 1.0.
7 // (See accompanying file LICENSE_1_0.txt
8 // or copy at http://www.boost.org/LICENSE_1_0.txt)
13 #pragma warning (disable:4127 4512)
16 #if !defined(TEST_FLOAT) && !defined(TEST_DOUBLE) && !defined(TEST_LDOUBLE) && !defined(TEST_REAL_CONCEPT)
20 # define TEST_REAL_CONCEPT
23 #include <boost/math/tools/test.hpp>
24 #include <boost/math/concepts/real_concept.hpp> // for real_concept
25 #include <boost/math/distributions/non_central_f.hpp> // for chi_squared_distribution
26 #define BOOST_TEST_MAIN
27 #include <boost/test/unit_test.hpp> // for test_main
28 #include <boost/test/results_collector.hpp>
29 #include <boost/test/unit_test.hpp>
30 #include <boost/test/tools/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE
31 #include "test_out_of_range.hpp"
33 #include "functor.hpp"
34 #include "handle_test_result.hpp"
41 using std::numeric_limits
;
43 #define BOOST_CHECK_CLOSE_EX(a, b, prec, i) \
45 unsigned int failures = boost::unit_test::results_collector.results( boost::unit_test::framework::current_test_case().p_id ).p_assertions_failed;\
46 BOOST_CHECK_CLOSE(a, b, prec); \
47 if(failures != boost::unit_test::results_collector.results( boost::unit_test::framework::current_test_case().p_id ).p_assertions_failed)\
49 std::cerr << "Failure was at row " << i << std::endl;\
50 std::cerr << std::setprecision(35); \
51 std::cerr << "{ " << data[i][0] << " , " << data[i][1] << " , " << data[i][2];\
52 std::cerr << " , " << data[i][3] << " , " << data[i][4] << " } " << std::endl;\
56 #define BOOST_CHECK_EX(a, i) \
58 unsigned int failures = boost::unit_test::results_collector.results( boost::unit_test::framework::current_test_case().p_id ).p_assertions_failed;\
60 if(failures != boost::unit_test::results_collector.results( boost::unit_test::framework::current_test_case().p_id ).p_assertions_failed)\
62 std::cerr << "Failure was at row " << i << std::endl;\
63 std::cerr << std::setprecision(35); \
64 std::cerr << "{ " << data[i][0] << " , " << data[i][1] << " , " << data[i][2];\
65 std::cerr << " , " << data[i][3] << " , " << data[i][4] << " } " << std::endl;\
69 void expected_results()
72 // Printing out the compiler/stdlib/platform names,
73 // we do this to make it easier to mark up expected error rates.
75 std::cout
<< "Tests run with " << BOOST_COMPILER
<< ", "
76 << BOOST_STDLIB
<< ", " << BOOST_PLATFORM
<< std::endl
;
79 template <class RealType
>
80 RealType
naive_pdf(RealType v1
, RealType v2
, RealType l
, RealType f
)
86 RealType term
= -l
/2 + log(l
/2) * i
;
87 term
+= boost::math::lgamma(v2
/2 + v1
/2+i
) - (boost::math::lgamma(v2
/2) + boost::math::lgamma(v1
/2+i
));
88 term
-= boost::math::lgamma(RealType(i
+1));
89 term
+= log(v1
/v2
) * (v1
/2+i
) + log(v2
/ (v2
+ v1
* f
)) * ((v1
+ v2
) / 2 + i
);
90 term
+= log(f
) * (v1
/2 - 1 + i
);
93 if((term
/sum
< boost::math::tools::epsilon
<RealType
>()) || (term
== 0))
99 template <class RealType
>
103 RealType ncp
, // non-centrality param
104 RealType x
, // F statistic
106 RealType Q
, // Complement of CDF
108 RealType tol
) // Test tolerance
110 boost::math::non_central_f_distribution
<RealType
> dist(a
, b
, ncp
);
112 cdf(dist
, x
), P
, tol
);
114 pdf(dist
, x
), D
, tol
);
115 if(boost::math::tools::digits
<RealType
>() > 50)
118 // The "naive" pdf calculation fails at float precision.
121 pdf(dist
, x
), naive_pdf(a
, b
, ncp
, x
), tol
);
124 if((P
< 0.99) && (Q
< 0.99))
127 // We can only check this if P is not too close to 1,
128 // so that we can guarantee Q is reasonably free of error:
131 cdf(complement(dist
, x
)), Q
, tol
);
133 quantile(dist
, P
), x
, tol
* 10);
135 quantile(complement(dist
, Q
)), x
, tol
* 10);
137 if(boost::math::tools::digits
<RealType
>() > 50)
140 // Sanity check mode:
142 RealType m
= mode(dist
);
143 RealType p
= pdf(dist
, m
);
144 BOOST_CHECK(pdf(dist
, m
* (1 + sqrt(tol
) * 10)) <= p
);
145 BOOST_CHECK(pdf(dist
, m
* (1 - sqrt(tol
)) * 10) <= p
);
149 template <class RealType
> // Any floating-point type RealType.
150 void test_spots(RealType
)
152 RealType tolerance
= boost::math::tools::epsilon
<RealType
>() * 10000;
154 cout
<< "Tolerance = " << (tolerance
/ 100) << "%." << endl
;
157 // Spot tests from Mathematica computed values:
160 RealType(5), // alpha
162 RealType(1), // non-centrality param
163 RealType(1.5), // F statistic
164 RealType(0.49845842011686358665786775091245664L), // CDF
165 RealType(1 - 0.49845842011686358665786775091245664L), // Complement of CDF
166 RealType(0.20251311620629730205859816288225385L), // PDF
167 RealType(tolerance
));
170 RealType(2), // alpha
172 RealType(1), // non-centrality param
173 RealType(2), // F statistic
174 RealType(0.64938711196845800322066756609406894L), // CDF
175 RealType(1 - 0.64938711196845800322066756609406894L), // Complement of CDF
176 RealType(0.15512617916132011524583796078456003L), // PDF
177 RealType(tolerance
));
179 RealType(100), // alpha
181 RealType(15), // non-centrality param
182 RealType(105), // F statistic
183 RealType(0.99996207325249555786258005958906310L), // CDF
184 RealType(0.000037926747504442137419940410936905407L), // Complement of CDF
185 RealType(8.9562292619539161551049126260104435e-7), // PDF
186 RealType(tolerance
* 10));
188 RealType(100), // alpha
190 RealType(15), // non-centrality param
191 RealType(1.5), // F statistic
192 RealType(0.57592315596686179870591317303126895L), // CDF
193 RealType(1 - 0.57592315596686179870591317303126895L), // Complement of CDF
194 RealType(0.36743745541686900593212039288061162L), // PDF
195 RealType(tolerance
* 5));
197 RealType(5), // alpha
198 RealType(100), // beta
199 RealType(102), // non-centrality param
200 RealType(25), // F statistic
201 RealType(0.74993383259829917593356265102363267L), // CDF
202 RealType(1 - 0.74993383259829917593356265102363267L), // Complement of CDF
203 RealType(0.054467600423154020554785779421659007L), // PDF
204 RealType(tolerance
* 5));
206 RealType(85), // alpha
207 RealType(100), // beta
208 RealType(0.5), // non-centrality param
209 RealType(1.25), // F statistic
210 RealType(0.85228624977948142884820398473385455L), // CDF
211 RealType(0.14771375022051857115179601526614545L), // Complement of CDF
212 RealType(0.88510283331631848299643323511414868L), // PDF
213 RealType(tolerance
* 5));
216 // Spot tests use values computed by the R statistical
217 // package and the pbeta and dbeta functions:
219 tolerance
= (std::max
)(
220 boost::math::tools::epsilon
<RealType
>() * 100,
221 (RealType
)1e-6) * 100;
224 RealType(85), // alpha
225 RealType(100), // beta
226 RealType(245), // non-centrality param
227 RealType(3.5), // F statistic
228 RealType(0.2697244), // CDF
229 RealType(1 - 0.2697244), // Complement of CDF
230 RealType(0.54352369104452836465948073900030320L), // PDF
231 RealType(tolerance
));
236 // 5 eps expressed as a persentage, otherwise the limit of the test data:
238 RealType tol2
= (std::max
)(boost::math::tools::epsilon
<RealType
>() * 500, RealType(1e-25));
241 boost::math::non_central_f_distribution
<RealType
> dist(20, 15, 30);
245 , static_cast<RealType
>(2.8846153846153846153846153846154L), tol2
);
249 , static_cast<RealType
>(2.1422807961269499731038192576654L), tol2
);
252 standard_deviation(dist
)
253 , sqrt(variance(dist
)), tol2
);
257 , pdf(dist
, x
) / cdf(complement(dist
, x
)), tol2
);
258 // cumulative hazard:
261 , -log(cdf(complement(dist
, x
))), tol2
);
262 // coefficient_of_variation:
264 coefficient_of_variation(dist
)
265 , standard_deviation(dist
) / mean(dist
), tol2
);
270 static_cast<RealType
>(0.5)), static_cast<RealType
>(tol2
));
274 , static_cast<RealType
>(2.070019130232759428074835788815387743293972985542L), sqrt(tolerance
));
278 , static_cast<RealType
>(2.1011821125804540669752380228510691320707051692719L), tol2
);
282 , 3 + kurtosis_excess(dist
), tol2
);
285 kurtosis_excess(dist
)
286 , static_cast<RealType
>(13.225781681053154767604638331440974359675882226873L), tol2
);
288 // Error handling checks:
289 check_out_of_range
<boost::math::non_central_f_distribution
<RealType
> >(1, 1, 1);
290 BOOST_MATH_CHECK_THROW(pdf(boost::math::non_central_f_distribution
<RealType
>(0, 1, 1), 0), std::domain_error
);
291 BOOST_MATH_CHECK_THROW(pdf(boost::math::non_central_f_distribution
<RealType
>(-1, 1, 1), 0), std::domain_error
);
292 BOOST_MATH_CHECK_THROW(pdf(boost::math::non_central_f_distribution
<RealType
>(1, 0, 1), 0), std::domain_error
);
293 BOOST_MATH_CHECK_THROW(pdf(boost::math::non_central_f_distribution
<RealType
>(1, -1, 1), 0), std::domain_error
);
294 BOOST_MATH_CHECK_THROW(quantile(boost::math::non_central_f_distribution
<RealType
>(1, 1, 1), -1), std::domain_error
);
295 BOOST_MATH_CHECK_THROW(quantile(boost::math::non_central_f_distribution
<RealType
>(1, 1, 1), 2), std::domain_error
);
296 } // template <class RealType>void test_spots(RealType)
298 BOOST_AUTO_TEST_CASE( test_main
)
300 BOOST_MATH_CONTROL_FP
;
301 // Basic sanity-check spot values.
303 // (Parameter value, arbitrarily zero, only communicates the floating point type).
304 test_spots(0.0F
); // Test float.
305 test_spots(0.0); // Test double.
306 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
307 test_spots(0.0L); // Test long double.
308 #if !BOOST_WORKAROUND(__BORLANDC__, BOOST_TESTED_AT(0x582))
309 test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
314 } // BOOST_AUTO_TEST_CASE( test_main )