boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00
Code Issues Projects Releases Wiki Activity

Removed remaining old style distribution code (we can get it back from cvs if necessary).

Browse Source 
Updated error handling.
Added better chi-squared tests.


[SVN r3165]
This commit is contained in:
John Maddock
2006-08-29 16:41:07 +00:00
parent 9a3ade1320
commit 00a347bee0
7 changed files with 54 additions and 1049 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
include/boost/math/distributions/students_t.hpp
View File

@@ -261,7 +261,7 @@ RealType students_t_distribution<RealType>::estimate_degrees_of_freedom(
" either there is no answer to how many degrees of freedom are required"
" or the answer is infinite. Current best guess is %1%", result);
}
return result + 1;
return result;
}
template <class RealType>

249
include/boost/math/special_functions/chisqr.hpp
View File

@@ -1,249 +0,0 @@
// chisqr.hpp chisqr cumulative distribution function.
// Copyright John Maddock & Paul A. Bristow 2006.
// 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)
// The Chisqr Cumulative Distribution Function (CDF) is given
// by an incomplete gamma integral function.
// http://en.wikipedia.org/wiki/Chi-squared_distribution
// The chisqr-distribution is used for statistical testing:
// http://en.wikipedia.org/wiki/Pearson%27s_chi-square_test
// "Pearson's is one of a variety of chi-square tests –
// statistical procedures whose results are evaluated
// by reference to the chi-square distribution.
// It tests a null hypothesis that the relative frequencies
// of occurrence of observed events follow a specified frequency distribution.
// The events must be mutually exclusive.
// One of the simplest examples is the hypothesis that
// an ordinary six-sided die is "fair": all six outcomes occur equally often.
// Chi-square is calculated by finding the difference between
// each observed and theoretical frequency, squaring them,
// dividing each by the theoretical frequency, and taking the sum of the results:
// chi^2 = sum_{i=1 to 6} {(O_i - E_i)^2 / E_i}
// where:
// O_i = an observed frequency
// E_i = an expected (theoretical) frequency, asserted by the null hypothesis.
// For very large values of both degrees_of_freedom1 and degrees_of_freedom2,
// greater than 10^5, a normal approximation is used.
// If only one of degrees_of_freedom1 or degrees_of_freedom2 is
// greater than 10^5 then a chisqr approximation is used,
// see Abramowitz and Stegun (1965).
// TODO?
#ifndef BOOST_MATH_SPECIAL_CHISQR_HPP
#define BOOST_MATH_SPECIAL_CHISQR_HPP
#include <boost/math/special_functions/gamma.hpp> // for gamma.
#include <boost/math/tools/roots.hpp>
#include <boost/math/tools/error_handling.hpp> // for domain_error, logic_error.
#include <boost/math/tools/promotion.hpp> // for promotion
#include <boost/type_traits/is_floating_point.hpp>
#include <boost/type_traits/is_integral.hpp>
#include <boost/type_traits/is_same.hpp>
#include <boost/mpl/if.hpp>
namespace boost
{
namespace math
{
namespace detail
{ // Implementations called by actual functions
// with arguments that, if necessary,
// have been promoted from ArithmeticType to RealType.
template <class RealType>
RealType chisqr_imp(RealType degrees_of_freedom, RealType chisqr)
{
// Implementation of probability of chisqr.
// Returns the area under the left hand tail (from 0 to x)
// of the Chi square probability density function
// with v degrees of freedom.
using boost::math::gamma_Q; // gamma_Q(degrees_of_freedom/2, chisqr/2)
using boost::math::tools::domain_error;
using boost::math::tools::logic_error;
// Degrees of freedom argument may be integral, signed or unsigned, or floating-point, or User Defined real.
if(degrees_of_freedom <= 0)
{ // Degrees of freedom must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "degrees of freedom argument is %1%, but must be > 0 !", degrees_of_freedom);
}
if(chisqr < 0)
{ // chisqr must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "chisqr argument is %1%, but must be > 0 !", chisqr);
}
// Calculate probability of chisqr using the incomplete gamma integral function.
RealType probability = gamma_Q(degrees_of_freedom / 2, chisqr / 2);
// Expect 0 <= probability <= 1.
// Numerical errors might cause probability to be slightly outside the range < 0 or > 1.
// This might cause trouble downstream, so warn, possibly throw exception, but constrain to the limits.
if (probability < 0)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "probability %1% is < 0, so has been constrained to zero !", probability);
return 0; // Constrain to zero, if logic_error does not throw.
}
if(probability > 1)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "probability %1% is > 1, so has been constrained to unity!", probability);
return 1; // Constrain to unity, if logic_error does not throw.
}
return probability;
} // chisqr_imp
template <class RealType>
RealType chisqr_c_imp(RealType degrees_of_freedom, RealType chisqr)
{
// Implementation of probability of chisqr complemented.
// Returns the area under the right hand tail (from x to infinity)
// of the chisqr probability density function
// with v degrees of freedom:
using boost::math::gamma_Q; // gamma_Q(degrees_of_freedom/2, chisqr/2)
using boost::math::tools::domain_error;
using boost::math::tools::logic_error;
// Degrees of freedom argument may be integral, signed or unsigned, or floating-point, or User Defined real.
if(degrees_of_freedom <= 0)
{ // Degrees of freedom must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "degrees of freedom argument is %1%, but must be > 0 !", degrees_of_freedom);
}
if(chisqr < 0)
{ // chisqr must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "chisqr argument is %1%, but must be > 0 !", chisqr);
}
// Calculate probability of chisqr using the incomplete beta function.
RealType probability = gamma_Q(degrees_of_freedom / 2, chisqr / 2);
// Expect 0 <= probability <= 1.
// Numerical errors might cause probability to be slightly outside the range < 0 or > 1.
// This might cause trouble downstream, so warn, possibly throw exception, but constrain to the limits.
if (probability < 0)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "probability %1% is < 0, so has been constrained to zero !", probability);
return 0; // Constrain to zero, if logic_error does not throw.
}
if(probability > 1)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "probability %1% is > 1, so has been constrained to unity!", probability);
return 1; // Constrain to unity, if logic_error does not throw.
}
return probability;
} // chisqr_imp
template <class RealType>
RealType chisqr_inv_imp(RealType degrees_of_freedom, RealType probability)
{
// Implementation of inverse of chisqr distribution.
// Finds the chisqr argument x such that the integral
// from x to infinity of the chisqr density
// is equal to the given cumulative probability y.
using boost::math::gamma_Q; // gamma_Q(degrees_of_freedom/2, chisqr/2)
using boost::math::tools::domain_error;
using boost::math::tools::logic_error;
// Degrees of freedom argument may be integral, signed or unsigned, or floating-point, or User Defined real.
if(degrees_of_freedom <= 0)
{ // Degrees of freedom must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "degrees of freedom argument is %1%, but must be > 0 !", degrees_of_freedom);
}
if((probability < 0) || (probability > 1))
{ // probability must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "probability argument is %1%, but must be >= 0 and =< 1 !", probability);
}
// Calculate chisqr from probability & degrees_of_freedom using the inverse gamma integral function.
RealType chisqr = gamma_Q_inv(degrees_of_freedom / 2, probability) * 2;
// Expect chisqr > 0.
// Numerical errors might cause probability to be slightly outside the range < 0 or > 1.
// This might cause trouble downstream, so warn, possibly throw exception, but constrain to the limits.
if (chisqr < 0)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "chisqr %1% is < 0, so has been constrained to zero !", probability);
return 0; // Constrain to zero, if logic_error does not throw.
}
return probability;
} // chisqr_inv_imp
template <class RealType>
RealType chisqr_df_inv_imp(RealType chisqr, RealType probability)
{
// Implementation of inverse (complemented?) chisqr distribution.
// Finds the degress_fo_freedom argument x such that the
// integral from x to infinity of the chisqr density
// is equal to the given cumulative probability y.
using boost::math::gamma_Q_inv; // gamma_Q(degrees_of_freedom/2, chisqr/2)
using boost::math::tools::domain_error;
using boost::math::tools::logic_error;
// Degrees of freedom argument may be integral, signed or unsigned, or floating-point, or User Defined real.
if(chisqr <= 0)
{ // chisqr must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "chisqr argument is %1%, but must be > 0 !", chisqr);
}
if ((probability < 0) || (probability > 1))
{ // chisqr must be > 0!
return domain_error<RealType>(BOOST_CURRENT_FUNCTION, "probability argument is %1%, but must be >= 0 and =< 1 !", chisqr);
}
// Calculate degrees_of_freedom from probability & chisqr using the inverse gamma integral function??
RealType degrees_of_freedom = gamma_Q_inv(chisqr / 2, probability) * 2;// TODO <<<<<<<<<<<<<<<<<<<<<<<<<<<<<
// Expect degrees_of_freedom >= 0.
if (degrees_of_freedom < 0)
{
logic_error<RealType>(BOOST_CURRENT_FUNCTION, "degrees_of_freedom %1% is < 0, so has been constrained to zero !", probability);
return 0; // Constrain to zero, if logic_error does not throw.
}
return probability;
} // chisqr_inv_imp
} // namespace detail
template <class ArithmeticType, class RealType> // probability chisqr(degrees_of_freedom, chisqr)
inline typename tools::promote_arg2<RealType, ArithmeticType>::type
// return type is the wider of the two (?promoted) floating point types.
chisqr(ArithmeticType degrees_of_freedom, RealType chisqr)
{
typedef typename tools::promote_arg2<RealType, ArithmeticType>::type promote_type; // Arguments type.
return detail::chisqr_imp(static_cast<promote_type>(degrees_of_freedom), static_cast<promote_type>(chisqr));
} // chisqr
template <class ArithmeticType, class RealType> // complement probability chisqr_c(degrees_of_freedom, chisqr)
inline typename tools::promote_arg2<RealType, ArithmeticType>::type
chisqr_c(ArithmeticType degrees_of_freedom, RealType chisqr)
{
typedef typename tools::promote_arg2<RealType, ArithmeticType>::type promote_type; // Arguments type.
return detail::chisqr_c_imp(static_cast<promote_type>(degrees_of_freedom), static_cast<promote_type>(chisqr));
} // chisqr_c
template <class ArithmeticType, class RealType> // chisqr = chisqr_inv(degrees_of_freedom, probability)
inline typename tools::promote_arg2<RealType, ArithmeticType>::type
chisqr_inv(ArithmeticType degrees_of_freedom, RealType probability)
{
typedef typename tools::promote_arg2<RealType, ArithmeticType>::type promote_type; // Arguments type.
return detail::chisqr_inv_imp(static_cast<promote_type>(degrees_of_freedom), static_cast<promote_type>(probability));
} // chisqr_inv
template <class ArithmeticType, class RealType> // degrees_of_freedom = chisqr_inv_df(chisqr, probability)
inline typename tools::promote_arg2<RealType, RealType>::type // <RealType, ArithmeticType> but both RealType.
chisqr_df_inv(RealType chisqr, ArithmeticType probability)
{
typedef typename tools::promote_arg2<RealType, ArithmeticType>::type promote_type; // Arguments type.
return detail::chisqr_df_inv_imp(static_cast<promote_type>(chisqr), static_cast<promote_type>(probability));
} // chisqr_inv_df
} // namespace math
} // namespace boost
#endif // BOOST_MATH_SPECIAL_CHISQR_HPP

6
test/test_chi_squared.cpp
View File

@@ -292,7 +292,7 @@ void test_spots(RealType)
// Basic sanity checks, test data is to three decimal places only
// so set tolerance to 0.001 expressed as a persentage.
RealType tolerance = 0.001 * 100;
RealType tolerance = 0.001f * 100;
cout << "Tolerance = " << tolerance << "%." << endl;
@@ -301,7 +301,7 @@ void test_spots(RealType)
using ::boost::math::cdf;
using ::boost::math::pdf;
for(int i = 0; i < sizeof(upper_critical_values) / sizeof(upper_critical_values[0]); ++i)
for(unsigned i = 0; i < sizeof(upper_critical_values) / sizeof(upper_critical_values[0]); ++i)
{
test_spot(
static_cast<RealType>(upper_critical_values[i][0]), // degrees of freedom
@@ -335,7 +335,7 @@ void test_spots(RealType)
tolerance);
}
for(int i = 0; i < sizeof(lower_critical_values) / sizeof(lower_critical_values[0]); ++i)
for(unsigned i = 0; i < sizeof(lower_critical_values) / sizeof(lower_critical_values[0]); ++i)
{
test_spot(
static_cast<RealType>(lower_critical_values[i][0]), // degrees of freedom

200
test/test_chisqr.cpp
View File

@@ -1,200 +0,0 @@
// test_chisqr.cpp
// Copyright Paul A. Bristow 2006.
// 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)
// Basic sanity test for chisqr Cumulative Distribution Function.
#define BOOST_MATH_THROW_ON_DOMAIN_ERROR
#ifdef _MSC_VER
# pragma warning(disable: 4127) // conditional expression is constant.
# pragma warning(disable: 4100) // unreferenced formal parameter.
# pragma warning(disable: 4512) // assignment operator could not be generated.
# pragma warning(disable: 4996) // 'std::char_traits<char>::copy' was declared deprecated.
#endif
#include <boost/math/special_functions/chisqr.hpp> // for chisqr
using ::boost::math::chisqr;
#include <boost/math/concepts/real_concept.hpp> // for real_concept
using ::boost::math::concepts::real_concept;
#include <boost/test/included/test_exec_monitor.hpp> // for test_main
#include <boost/test/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE
#include <iostream>
using std::cout;
using std::endl;
#include <limits>
using std::numeric_limits;
template <class FPT> // Any floating-point type FPT.
void test_spots(FPT)
{
// Basic sanity checks, tolerance is about numeric_limits<FPT>::digits10 decimal places,
// guaranteed for type FPT, eg 6 for float, 15 for double,
// expressed as a percentage (so -2) for BOOST_CHECK_CLOSE,
int decdigits = numeric_limits<FPT>::digits10;
decdigits -= 0; // Perhaps allow some decimal digits margin of numerical error.
FPT tolerance = static_cast<FPT>(std::pow(10., -(decdigits-2))); // 1e-6 (-2 so as %)
tolerance *= 1; // Allow some bit(s) small margin (2 means + or - 1 bit) of numerical error.
// Typically 2e-13% = 2e-15 as fraction for double.
// Sources of spot test values:
// http://www.vias.org/tmdatanaleng/cc_distri_calculator.html give some useful info,
// BUT the calculator link is broken.
// http://www.vias.org/simulations/simusoft_distcalc.html
// Distcalc version 1.2 Copyright 2002 H Lohninger, TU Wein
// H.Lohninger: Teach/Me Data Analysis, Springer-Verlag, Berlin-New York-Tokyo, 1999. ISBN 3-540-14743-8
// The Windows calculator is available zipped distcalc.exe for download at:
// http://www.vias.org/simulations/simu_stat.html
// Also usd a Java Chi Square calculator at:
// http://www.stat.sc.edu/~west/applets/chisqdemo.html
// but seems even less accurate.
// Many be some combinations for which the result is 'exact', or at least is to 40 decimal digits.
// 40 decimal digits includes 128-bit significand User Defined Floating-Point types,
// but these are as yet undiscovered for chi-sqr.
// Best source of accurate values is
// Mathworld online calculator (40 decimal digits precision, suitable for up to 128-bit significands)
// http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp?name=GammaRegularized
// GammaRegularized is same as gamma incomplete, gamma or gamma_Q(a, x) or Q(a, z).
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(5.), // degrees_of_freedom - as floating-point.
static_cast<FPT>(11.0705)), // chisqr
static_cast<FPT>(0.04999995542804363455469285561840275602943), // probability.
// Source Mathworld online calculator (40 decimal digits precision)
// Q(5.0/2, 11.0705/2) = 0.049999955428043634554692855618402756029434129418384
tolerance);
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(4.), // degrees_of_freedom - as floating-point.
static_cast<FPT>(4.)), // chisqr
static_cast<FPT>(0.4060058497098380756819984849174532102229), // probability.
// Q(4/2, 4/2) = 0.049999955428043634554692855618402756029434129418384
tolerance);
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(1), // degrees_of_freedom - as floating-point.
static_cast<FPT>(10.)), // chisqr
static_cast<FPT>(0.001565402258002549677499803978385902310435), // probability.
tolerance);
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(10.), // degrees_of_freedom - as floating-point.
static_cast<FPT>(1.)), // chisqr
static_cast<FPT>(0.9998278843700441592218882959620240287207), // probability.
tolerance);
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(100.), // degrees_of_freedom - as floating-point.
static_cast<FPT>(5.)), // chisqr
static_cast<FPT>(0.999999999999999999814527311613020069945), // probability.
tolerance);
/*
// Try bad degrees_of_freedom argument for chisqr.
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(1), // degrees_of_freedom - as floating-point.
static_cast<FPT>(-2.)), // chisqr
static_cast<FPT>(0.999999999999999999), // probability.
tolerance);
// std::domain_error: Error in function float __cdecl boost::math::detail::chisqr_imp<float>(float,float): degrees of freedom argument is -1, but must be > 0 !
// Try bad chisqr argument for chisqr.
BOOST_CHECK_CLOSE(
chisqr(
static_cast<FPT>(1), // degrees_of_freedom - as floating-point.
static_cast<FPT>(-2.)), // chisqr
static_cast<FPT>(0.999999999999999999), // probability.
tolerance);
// std::domain_error: Error in function float __cdecl boost::math::detail::chisqr_imp<float>(float,float): chisqr argument is -2, but must be > 0 !
*/
} // template <class FPT>void test_spots(FPT)
int test_main(int, char* [])
{
// Basic sanity-check spot values.
#ifdef BOOST_MATH_THROW_ON_DOMAIN_ERROR
cout << "BOOST_MATH_THROW_ON_DOMAIN_ERROR" << " is defined to throw on domain error." << endl;
#else
cout << "BOOST_MATH_THROW_ON_DOMAIN_ERROR" << " is NOT defined, so NO throw on domain error." << endl;
#endif
// (Parameter value, arbitrarily zero, only communicates the floating point type).
test_spots(0.0F); // Test float.
test_spots(0.0); // Test double.
test_spots(0.0L); // Test long double.
test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
return 0;
} // int test_main(int, char* [])
/*
Output:
Compiling...
test_chisqr.cpp
Linking...
Embedding manifest...
Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe"
Running 1 test case...
BOOST_MATH_THROW_ON_DOMAIN_ERROR is NOT defined so NO throw on domain error.
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(116): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(-1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 0.0001%
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(124): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 0.0001%
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(116): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(-1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 1e-013%
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(124): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 1e-013%
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(116): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(-1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 1e-013%
../../../../../../boost-sandbox/libs/math_functions/test/test_chisqr.cpp(124): error in "test_main_caller( argc, argv )": difference between chisqr( static_cast<FPT>(1), static_cast<FPT>(-2.)){1.#QNAN} and static_cast<FPT>(0.999999999999999999){1} exceeds 1e-013%
unknown location(0): fatal error in "test_main_caller( argc, argv )": std::domain_error: Error in function class boost::math::concepts::real_concept __cdecl boost::math::detail::chisqr_imp<class boost::math::concepts::real_concept>(class boost::math::concepts::real_concept,class boost::math::concepts::real_concept): degrees of freedom argument is -1, but must be > 0 !
..\..\..\..\..\..\boost-sandbox\libs\math_functions\test\test_chisqr.cpp(116): last checkpoint
*** 7 failures detected in test suite "Test Program"
Project : error PRJ0019: A tool returned an error code from "Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe""
Compiling...
test_chisqr.cpp
Linking...
Embedding manifest...
Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe"
Running 1 test case...
BOOST_MATH_THROW_ON_DOMAIN_ERROR is defined to throw on domain error.
unknown location(0): fatal error in "test_main_caller( argc, argv )": std::domain_error: Error in function float __cdecl boost::math::detail::chisqr_imp<float>(float,float): degrees of freedom argument is -1, but must be > 0 !
..\..\..\..\..\..\boost-sandbox\libs\math_functions\test\test_chisqr.cpp(116): last checkpoint
*** 1 failure detected in test suite "Test Program"
Project : error PRJ0019: A tool returned an error code from "Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe""
Compiling...
test_chisqr.cpp
Linking...
Embedding manifest...
Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe"
Running 1 test case...
BOOST_MATH_THROW_ON_DOMAIN_ERROR is defined to throw on domain error.
unknown location(0): fatal error in "test_main_caller( argc, argv )": std::domain_error: Error in function float __cdecl boost::math::detail::chisqr_imp<float>(float,float): chisqr argument is -2, but must be > 0 !
..\..\..\..\..\..\boost-sandbox\libs\math_functions\test\test_chisqr.cpp(116): last checkpoint
*** 1 failure detected in test suite "Test Program"
Project : error PRJ0019: A tool returned an error code from "Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_chisqr.exe""
*/

344
test/test_fisher.cpp
View File

@@ -1,344 +0,0 @@
// Copyright Paul A. Bristow 2006.
// 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)
// test_fisher.cpp
// Basic sanity test for Fisher Cumulative Distribution Function.
#ifdef _MSC_VER
# pragma warning(disable: 4127) // conditional expression is constant.
# pragma warning(disable: 4100) // unreferenced formal parameter.
# pragma warning(disable: 4512) // assignment operator could not be generated.
# pragma warning(disable: 4996) // 'std::char_traits<char>::copy' was declared deprecated.
#endif
#define BOOST_MATH_THROW_ON_DOMAIN_ERROR
#include <boost/math/concepts/real_concept.hpp> // for real_concept
#include <boost/math/special_functions/fisher.hpp>
#include <boost/test/included/test_exec_monitor.hpp> // test_main
#include <boost/test/floating_point_comparison.hpp> // BOOST_CHECK_CLOSE
#include <limits>
using std::numeric_limits;
template <class FPT>
void test_spots(FPT)
{
// Basic sanity checks, tolerance is 6 decimal places
// expressed as a percentage (so -2) for BOOST_CHECK_CLOSE,
// One check per domain of the implementation:
int decdigits = numeric_limits<FPT>::digits10; // guaranteed for this type, eg 6 for float, 15 for double.
decdigits -= 0; // Perhaps allow some decimal digits margin of numerical error.
FPT tolerance = static_cast<FPT>(std::pow(10., -(decdigits-2))); // 1e-6 (as %)
tolerance *= 2; // Allow some small margin (2 means 1 least decimal digit?) of numerical error.
// Typically 2e-13% = 2e-15 as fraction for double.
// Whereas tolerance *= 1.5 fails thus
// difference between ::boost::math::fisher( static_cast<FPT>(1.), static_cast<FPT>(2.), static_cast<FPT>(6.5333333333333))
// {0.12500000000000022} and static_cast<FPT>(0.125){0.125} exceeds 1.5e-013%
// FPT tolerance = static_cast<FPT>(std::pow(10., -(4-2))); // 1e-4 (as %) OK with misc test values.
// PT tolerance = static_cast<FPT>(std::pow(10., -(6-2))); // 1e-6 (as %) fails for most values.
// This test only passes at 1e-5 because probability value is less accurate,
// a digit in 6th decimal place, although calculated using values from
// http://www.danielsoper.com/statcalc/calc04.aspx
// BUT I conclude that these values are only accurate to about 1e-4,
// so not very useful for testing!
/*
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(1.), // df2
static_cast<FPT>(161.4476)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(19.0000)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(10.), // df2
static_cast<FPT>(4.9646)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(100.), // df2
static_cast<FPT>(3.9361)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(1000000.), // df2
static_cast<FPT>(3.8415)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(10.), // df1
static_cast<FPT>(1.), // df2
static_cast<FPT>(241.8818)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(8.5263)), // F
static_cast<FPT>(0.1), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(18.5128)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(98.5025)), // F
static_cast<FPT>(0.01), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(998.4979)), // F
static_cast<FPT>(0.001), // probability.
tolerance);
*/
// http://www.vias.org/tmdatanaleng/cc_distri_calculator.html give some useful info,
// BUT the calculator link is broken.
// http://www.vias.org/simulations/simusoft_distcalc.html
// Distcalc version 1.2 Copyright 2002 H Lohninger, TU Wein
// H.Lohninger: Teach/Me Data Analysis, Springer-Verlag, Berlin-New York-Tokyo, 1999. ISBN 3-540-14743-8
// The Windows calculator is available zipped distcalc.exe for download at:
// http://www.vias.org/simulations/simu_stat.html
// This interactive Windows program was used to find some combination for which the
// result appears to be exact. No doubt this can be done analytically too,
// by mathematicians!
// Some combinations for which the result is 'exact', or at least is to 40 decimal digits.
// 40 decimal digits includes 128-bit significand User Defined Floating-Point types.
// These all pass tests at near epsilon accuracy for the floating-point type.
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
// static_cast<FPT>(0.6666666666666666666666666666666666667)), // F
static_cast<FPT>(2.)/static_cast<FPT>(3.) ), // F
static_cast<FPT>(0.5), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(1.6)), // F
static_cast<FPT>(0.333333333333333333333333333333333333), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(1.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(6.5333333333333333333333333333333333)), // F
static_cast<FPT>(0.125), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(1.)), // F
static_cast<FPT>(0.5), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(3.)), // F
static_cast<FPT>(0.25), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(3.)), // F
static_cast<FPT>(0.25), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(7.)), // F
static_cast<FPT>(0.125), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(9.)), // F
static_cast<FPT>(0.1), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(19.)), // F
static_cast<FPT>(0.05), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(29.)), // F
static_cast<FPT>(0.03333333333333333333333333333333333333333), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(99.)), // F
static_cast<FPT>(0.01), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(4.), // df1
static_cast<FPT>(4.), // df2
static_cast<FPT>(9.)), // F
static_cast<FPT>(0.028), // probability.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher(
static_cast<FPT>(8.), // df1
static_cast<FPT>(8.), // df2
static_cast<FPT>(1.)), // F
static_cast<FPT>(0.5), // probability.
tolerance);
// These two don't pass at the 2 eps tolerance:
//BOOST_CHECK_CLOSE(
// ::boost::math::fisher(
// static_cast<FPT>(2.), // df1
// static_cast<FPT>(2.), // df2
// static_cast<FPT>(999.)), // F
// static_cast<FPT>(0.001), // probability.== 0.000999987125
//tolerance);
// And this even with 1 decimal digit tolerance.
// BOOST_CHECK_CLOSE(
// ::boost::math::fisher(
// static_cast<FPT>(2.), // df1
// static_cast<FPT>(2.), // df2
// static_cast<FPT>(9999.)), // F
// static_cast<FPT>(0.0001), // probability. == 0.000100016594 or 9.9999999999988987e-005
//tolerance);
//BOOST_CHECK_CLOSE(
// ::boost::math::fisher(
// static_cast<FPT>(2.), // df1
// static_cast<FPT>(7.), // df2
// static_cast<FPT>(99.)), // F
// static_cast<FPT>(0.00001), // probability.
//tolerance);
//BOOST_CHECK_CLOSE(
// ::boost::math::fisher(
// static_cast<FPT>(2.), // df1
// static_cast<FPT>(3.), // df2
// static_cast<FPT>(999.)), // F
// static_cast<FPT>(0.0000666666666666666666666666), // probability.
//tolerance);
// Inverse tests
BOOST_CHECK_CLOSE(
::boost::math::fisher_c_inv(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(0.03333333333333333333333333333333333333333)), // probability
static_cast<FPT>(29.), // F expected.
tolerance);
BOOST_CHECK_CLOSE(
::boost::math::fisher_inv(
static_cast<FPT>(2.), // df1
static_cast<FPT>(2.), // df2
static_cast<FPT>(1 - 0.03333333333333333333333333333333333333333)), // probability
static_cast<FPT>(29.), // F expected.
tolerance);
// Also note limit cases for F(1, infinity) == normal distribution
// F(1, n2) == Student's t distribution
// F(n1, infinity) == Chisq distribution
// These might allow some further cross checks?
} // template <class FPT>void test_spots(FPT)
int test_main(int, char* [])
{
// Basic sanity-check spot values.
// (Parameter value, arbitrarily zero, only communicates the floating point type).
test_spots(0.0F); // Test float.
test_spots(0.0); // Test double.
test_spots(0.0L); // Test long double.
test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
return 0;
} // int test_main(int, char* [])
/*
Output:
1>------ Rebuild All started: Project: test_fisher, Configuration: Debug Win32 ------
1>Deleting intermediate and output files for project 'test_fisher', configuration 'Debug|Win32'
1>Compiling...
1>test_fisher.cpp
1>Linking...
1>Embedding manifest...
1>Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_fisher.exe"
1>Running 1 test case...
1>*** No errors detected
1>Build Time 0:06
1>Build log was saved at "file://i:\boost-06-05-03-1300\libs\math\test\Math_test\test_fisher\Debug\BuildLog.htm"
1>test_fisher - 0 error(s), 0 warning(s)
========== Rebuild All: 1 succeeded, 0 failed, 0 skipped ==========
*/

252
test/test_negative_binomial.cpp
View File

@@ -1,252 +0,0 @@
// test_negative_binomial.cpp
// Copyright Paul A. Bristow 2006.
// 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)
// Basic sanity test for Negative Binomial Cumulative Distribution Function.
#define BOOST_MATH_THROW_ON_DOMAIN_ERROR
#ifdef _MSC_VER
//# pragma warning(disable: 4127) // conditional expression is constant.
//# pragma warning(disable: 4100) // unreferenced formal parameter.
//# pragma warning(disable: 4512) // assignment operator could not be generated.
//# pragma warning(disable: 4996) // 'std::char_traits<char>::copy' was declared deprecated.
//# pragma warning(disable: 4244) // conversion from 'double' to 'const float', possible loss of data.
#endif
#include <boost/math/special_functions/negative_binomial.hpp> // for negative_binomial
using ::boost::math::negative_binomial;
using ::boost::math::negative_binomial_c;
using ::boost::math::negative_binomial_inv;
#include <boost/math/concepts/real_concept.hpp> // for real_concept
using ::boost::math::concepts::real_concept;
#include <boost/test/included/test_exec_monitor.hpp> // for test_main
#include <boost/test/floating_point_comparison.hpp> // for BOOST_CHECK_CLOSE
#include <iostream>
using std::cout;
using std::endl;
#include <limits>
using std::numeric_limits;
template <class RealType> // RealType is any floating-point type.
void test_spots(RealType)
{
// Basic sanity checks, tolerance is about numeric_limits<RealType>::digits10 decimal places,
// guaranteed for type RealType, eg 6 for float, 15 for double,
// expressed as a percentage (so -2) for BOOST_CHECK_CLOSE,
int decdigits = numeric_limits<RealType>::digits10;
decdigits -= 1; // Perhaps allow some decimal digits margin of numerical error.
RealType tolerance = static_cast<RealType>(std::pow(10., -(decdigits-2))); // 1e-6 (-2 so as %)
tolerance *= 1; // Allow some bit(s) small margin (2 means + or - 1 bit) of numerical error.
// Typically 2e-13% = 2e-15 as fraction for double.
// Sources of spot test values:
// MathCAD
// pnbinom(k, n, p)
// Returns cumulative negative binomial distribution.
// qnbinom(p, n, q)
// Returns the inverse negative binomial distribution function,
// the smallest integer k so that pnbinom(k, n, q) >= p
// MathWorld definition: negaBin[n, p]
// number of failures that occur before achieving n successes, where probability of success in a trial is p.
// Cephes probability that k or fewer failures preceeds the nth success.
// Wikipedia For k + r trials with success probability p,
// gives the probability of k failures and r successes, with sucess on the last trial.
// Or probability of number of failures before the rth success.
// Many be some combinations for which the result is 'exact', or at least is to 40 decimal digits.
// 40 decimal digits includes 128-bit significand User Defined Floating-Point types,
// but these are as yet undiscovered for negative binomial. TODO
// Test negative binomial.
// negative_binomial(k, n, p)
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(0), // k. pnbinom(1,2,0.5) = 0.5
static_cast<RealType>(2), // n
static_cast<RealType>(0.5)), // mean
static_cast<RealType>(0.25), // probability 1/4
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(1), // k
static_cast<RealType>(4), // n
static_cast<RealType>(0.5)), // p
static_cast<RealType>(0.1875), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(1), // k
static_cast<RealType>(2), // n
static_cast<RealType>(0.5)), // pnbinom(1,2,0.5) = 0.5
static_cast<RealType>(0.5), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(2), // k
static_cast<RealType>(2), // n
static_cast<RealType>(0.5)), // pnbinom(1,2,0.5) = 0.5
static_cast<RealType>(0.6875), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(1), // k
static_cast<RealType>(3), // n
static_cast<RealType>(0.5)), // pnbinom(1,3,0.5) = 0.5
static_cast<RealType>(0.3125), // 0.6875
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(2), // k
static_cast<RealType>(2), // n
static_cast<RealType>(0.1)), //
static_cast<RealType>(0.0523), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial_c(
static_cast<RealType>(2), // k
static_cast<RealType>(2), // n
static_cast<RealType>(0.1)), //
static_cast<RealType>(1- 0.0523), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(26), // k
static_cast<RealType>(5), // n
static_cast<RealType>(0.4)), // p
static_cast<RealType>(0.998968624661119),
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(25), // k
static_cast<RealType>(5), // n
static_cast<RealType>(0.4)), // p
static_cast<RealType>(0.998489925933617),
tolerance);
BOOST_CHECK_CLOSE(negative_binomial_c( // complement.
static_cast<RealType>(26), // k
static_cast<RealType>(5), // n
static_cast<RealType>(0.4)), // p
static_cast<RealType>(1 - 0.998968624661119),
tolerance * 100); // fails at tolerance * 10
// 0.0010313753388809799
// 0.0010313753388808430
// Shows a loss of about 3 decimal digits accuracy.
// negative_binomial_inv(k, n, p) tests.
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(5), // k.
static_cast<RealType>(10), // n
static_cast<RealType>(0.5)), // trial probability
static_cast<RealType>(0.15087890625000011), // probability of failure.
tolerance);
BOOST_CHECK_CLOSE(negative_binomial_inv(
static_cast<RealType>(5), // k.
static_cast<RealType>(10), // n
static_cast<RealType>(0.15087890625000011)), // probability of failure = negative_binomial(5, 10, 0.5).
static_cast<RealType>(0.5), // result is probability of success
tolerance);
BOOST_CHECK_CLOSE(negative_binomial(
static_cast<RealType>(5), // k failures.
static_cast<RealType>(10), // n successes.
static_cast<RealType>(0.1)), // trial probability of success.
static_cast<RealType>(1.8662024800000033e-007), //
tolerance);
BOOST_CHECK_CLOSE(negative_binomial_inv(
static_cast<RealType>(5), // k failures.
static_cast<RealType>(10), // n successes.
static_cast<RealType>(1.8662024800000033e-007)), // probability of failure = negative_binomial(5, 10, 0.1)
static_cast<RealType>(0.1), // results is trial probability
tolerance);
} // template <class RealType>void test_spots(RealType)
int test_main(int, char* [])
{
// Basic sanity-check spot values.
#ifdef BOOST_MATH_THROW_ON_DOMAIN_ERROR
cout << "BOOST_MATH_THROW_ON_DOMAIN_ERROR" << " is defined to throw on domain error." << endl;
#else
cout << "BOOST_MATH_THROW_ON_DOMAIN_ERROR" << " is NOT defined, so NO throw on domain error." << endl;
#endif
//for (int i = 0; i < 10; i++)
//{
// cout << i << ' ' << negative_binomial(i, 2, 0.5) << endl;
//}
//cout << endl;
//cout.precision(17);
//for (double p = 0.; p < 1; p+= 0.1)
//{
// double y = negative_binomial(5, 10, p);
// double z = negative_binomial_inv(5, 10, y);
// cout << p << ' '<< y << ' ' << z << ' ' << p - z << endl;
//}
// (Parameter value, arbitrarily zero, only communicates the floating point type).
test_spots(0.0F); // Test float.
test_spots(0.0); // Test double.
test_spots(0.0L); // Test long double.
test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
return 0;
} // int test_main(int, char* [])
/*
Output:
------ Build started: Project: test_negative_binomial, Configuration: Debug Win32 ------
Compiling...
test_negative_binomial.cpp
Linking...
Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_negative_binomial.exe"
Running 1 test case...
BOOST_MATH_THROW_ON_DOMAIN_ERROR is defined to throw on domain error.
0 0.25
1 0.5
2 0.6875
3 0.8125
4 0.890625
5 0.9375
6 0.964844
7 0.980469
8 0.989258
9 0.994141
0 0 -1.#INF 1.#INF
0.10000000000000001 1.8662024800000033e-007 0.10000000000000002 -1.3877787807814457e-017
0.20000000000000001 0.00011322566246400017 0.20000000000000004 -2.7755575615628914e-017
0.30000000000000004 0.0036525210084360034 0.29999999999999999 5.5511151231257827e-017
0.40000000000000002 0.033833302884352018 0.39999999999999997 5.5511151231257827e-017
0.5 0.15087890625000011 0.5 0
0.59999999999999998 0.40321555041484813 0.60000000000000009 -1.1102230246251565e-016
0.69999999999999996 0.72162144020436414 0.70000000000000007 -1.1102230246251565e-016
0.79999999999999993 0.93894857038233592 0.79999999999999982 1.1102230246251565e-016
0.89999999999999991 0.99775032991495194 0.89999999999999947 4.4408920985006262e-016
0.99999999999999989 1 1.#INF -1.#INF
*** No errors detected
Build Time 0:05
Build log was saved at "file://i:\boost-06-05-03-1300\libs\math\test\Math_test\test_negative_binomial\Debug\BuildLog.htm"
test_negative_binomial - 0 error(s), 0 warning(s)
========== Build: 1 succeeded, 0 failed, 0 up-to-date, 0 skipped ==========
*/

50
test/test_students_t.cpp
View File

@@ -322,6 +322,56 @@ void test_spots(RealType T)
coefficient_of_variation(dist),
std::overflow_error);
// Parameter estimation. These results are close to but
// not identical to those reported on the NIST website at
// http://www.itl.nist.gov/div898/handbook/prc/section2/prc222.htm
// the NIST results appear to be calculated using a normal
// approximation, which slightly under-estimates the degrees of
// freedom required, particularly when the result is small.
//
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(0.5),
static_cast<RealType>(0.005),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
99);
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(1.5),
static_cast<RealType>(0.005),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
14);
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(0.5),
static_cast<RealType>(0.025),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
76);
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(1.5),
static_cast<RealType>(0.025),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
11);
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(0.5),
static_cast<RealType>(0.05),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
65);
BOOST_CHECK_EQUAL(
ceil(students_t_distribution<RealType>::estimate_degrees_of_freedom(
static_cast<RealType>(1.5),
static_cast<RealType>(0.05),
static_cast<RealType>(0.01),
static_cast<RealType>(1.0))),
9);
} // template <class RealType>void test_spots(RealType)
int test_main(int, char* [])