libs/math/test/test_students_t.cpp

   1 // Copyright Paul A. Bristow 2006, 2017.
   2 // Copyright John Maddock 2006.
   3
   4 // Use, modification and distribution are subject to the
   5 // Boost Software License, Version 1.0.
   6 // (See accompanying file LICENSE_1_0.txt
   7 // or copy at http://www.boost.org/LICENSE_1_0.txt)
   8
   9 // test_students_t.cpp
  10
  11 // http://en.wikipedia.org/wiki/Student%27s_t_distribution
  12 // http://www.itl.nist.gov/div898/handbook/eda/section3/eda3664.htm
  13
  14 // Basic sanity test for Student's t probability (quantile) (0. < p < 1).
  15 // and Student's t probability Quantile (0. < p < 1).
  16
  17 #ifdef _MSC_VER
  18 #  pragma warning (disable :4127) // conditional expression is constant.
  19 #endif
  20
  21 #define BOOST_TEST_MAIN
  22 #include <boost/test/unit_test.hpp> // Boost.Test
  23 #include <boost/test/tools/floating_point_comparison.hpp>
  24
  25 #include <boost/math/concepts/real_concept.hpp> // for real_concept
  26 #include <boost/math/tools/test.hpp> // for real_concept
  27 #include "test_out_of_range.hpp"
  28 #include <boost/math/distributions/students_t.hpp>
  29     using boost::math::students_t_distribution;
  30
  31 #include <iostream>
  32    using std::cout;
  33    using std::endl;
  34    using std::setprecision;
  35 #include <limits>
  36   using std::numeric_limits;
  37
  38 template <class RealType>
  39 RealType naive_pdf(RealType v, RealType t)
  40 {
  41    // Calculate the pdf of the students t in a deliberately
  42    // naive way, using equation (5) from
  43    // http://mathworld.wolfram.com/Studentst-Distribution.html
  44    // This is equivalent to, but a different method
  45    // to the one in the actual implementation, so can be used as
  46    // a very basic sanity check.  However some published values
  47    // would be nice....
  48
  49    using namespace std;  // for ADL
  50    using boost::math::beta;
  51
  52    //return pow(v / (v + t*t), (1+v) / 2) / (sqrt(v) * beta(v/2, RealType(0.5f)));
  53    RealType result = boost::math::tgamma_ratio((v+1)/2, v/2);
  54    result /= sqrt(v * boost::math::constants::pi<RealType>());
  55    result /= pow(1 + t*t/v, (v+1)/2);
  56    return result;
  57 }
  58
  59 template <class RealType>
  60 void test_spots(RealType)
  61 {
  62   // Basic sanity checks
  63
  64    RealType tolerance = static_cast<RealType>(1e-4); // 1e-6 (as %)
  65    // Some tests only pass at 1e-5 because probability value is less accurate,
  66    // a digit in 6th decimal place, although calculated using
  67    // a t-distribution generator (claimed 6 decimal digits) at
  68   // http://faculty.vassar.edu/lowry/VassarStats.html
  69    // http://faculty.vassar.edu/lowry/tsamp.html
  70    // df = 5, +/-t = 2.0, 1-tailed = 0.050970, 2-tailed = 0.101939
  71
  72    cout << "Tolerance for type " << typeid(RealType).name()  << " is " << tolerance << " %" << endl;
  73
  74    // http://en.wikipedia.org/wiki/Student%27s_t_distribution#Table_of_selected_values
  75   // Using tabulated value of t = 3.182 for 0.975, 3 df, one-sided.
  76
  77    // http://www.mth.kcl.ac.uk/~shaww/web_page/papers/Tdistribution06.pdf refers to:
  78
  79    // A lookup table of quantiles of the RealType distribution
  80   // for 1 to 25 in steps of 0.1 is provided in CSV form at:
  81   // www.mth.kcl.ac.uk/~shaww/web_page/papers/Tsupp/tquantiles.csv
  82    // gives accurate t of -3.1824463052837 and 3 degrees of freedom.
  83    // Values below are from this source, saved as tquantiles.xls.
  84    // DF are across the columns, probabilities down the rows
  85    // and the t- values (quantiles) are shown.
  86    // These values are probably accurate to nearly 64-bit double
  87   // (perhaps 14 decimal digits).
  88
  89    BOOST_CHECK_CLOSE(
  90       ::boost::math::cdf(
  91          students_t_distribution<RealType>(2),       // degrees_of_freedom
  92          static_cast<RealType>(-6.96455673428326)),  // t
  93          static_cast<RealType>(0.01),                // probability.
  94          tolerance); // %
  95
  96    BOOST_CHECK_CLOSE(
  97       ::boost::math::cdf(
  98          students_t_distribution<RealType>(5),       // degrees_of_freedom
  99          static_cast<RealType>(-3.36492999890721)),  // t
 100          static_cast<RealType>(0.01),                // probability.
 101          tolerance);
 102
 103    BOOST_CHECK_CLOSE(
 104       ::boost::math::cdf(
 105          students_t_distribution<RealType>(1),      // degrees_of_freedom
 106          static_cast<RealType>(-31830.988607907)),  // t
 107          static_cast<RealType>(0.00001),            // probability.
 108          tolerance);
 109
 110    BOOST_CHECK_CLOSE(
 111       ::boost::math::cdf(
 112          students_t_distribution<RealType>(25.),    // degrees_of_freedom
 113          static_cast<RealType>(-5.2410429995425)),  // t
 114          static_cast<RealType>(0.00001),            // probability.
 115          tolerance);
 116
 117    BOOST_CHECK_CLOSE(
 118       ::boost::math::cdf(
 119          students_t_distribution<RealType>(1),   // degrees_of_freedom
 120          static_cast<RealType>(-63661.97723)),   // t
 121          static_cast<RealType>(0.000005),        // probability.
 122          tolerance);
 123
 124     BOOST_CHECK_CLOSE(
 125       ::boost::math::cdf(
 126          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 127          static_cast<RealType>(-17.89686614)),   // t
 128          static_cast<RealType>(0.000005),        // probability.
 129          tolerance);
 130
 131     BOOST_CHECK_CLOSE(
 132       ::boost::math::cdf(
 133          students_t_distribution<RealType>(25.),  // degrees_of_freedom
 134          static_cast<RealType>(-5.510848412)),    // t
 135          static_cast<RealType>(0.000005),         // probability.
 136          tolerance);
 137
 138   BOOST_CHECK_CLOSE(
 139       ::boost::math::cdf(
 140          students_t_distribution<RealType>(10.),  // degrees_of_freedom
 141          static_cast<RealType>(-1.812461123)),    // t
 142          static_cast<RealType>(0.05),             // probability.
 143          tolerance);
 144
 145   BOOST_CHECK_CLOSE(
 146       ::boost::math::cdf(
 147          students_t_distribution<RealType>(10),  // degrees_of_freedom
 148          static_cast<RealType>(1.812461123)),    // t
 149          static_cast<RealType>(0.95),            // probability.
 150          tolerance);
 151
 152   BOOST_CHECK_CLOSE(
 153       ::boost::math::cdf(
 154          complement(
 155             students_t_distribution<RealType>(10),  // degrees_of_freedom
 156             static_cast<RealType>(1.812461123))),    // t
 157          static_cast<RealType>(0.05),            // probability.
 158          tolerance);
 159
 160   BOOST_CHECK_CLOSE(
 161       ::boost::math::cdf(
 162          students_t_distribution<RealType>(10),  // degrees_of_freedom
 163          static_cast<RealType>(9.751995491)),    // t
 164          static_cast<RealType>(0.999999),        // probability.
 165          tolerance);
 166
 167   BOOST_CHECK_CLOSE(
 168       ::boost::math::cdf(
 169          students_t_distribution<RealType>(10.),  // degrees_of_freedom - for ALL degrees_of_freedom!
 170          static_cast<RealType>(0.)),              // t
 171          static_cast<RealType>(0.5),              // probability.
 172          tolerance);
 173
 174
 175    // Student's t Inverse function tests.
 176   // Special cases
 177
 178   BOOST_MATH_CHECK_THROW(boost::math::quantile(
 179          students_t_distribution<RealType>(1.),  // degrees_of_freedom (ignored).
 180          static_cast<RealType>(0)), std::overflow_error); // t == -infinity.
 181
 182   BOOST_MATH_CHECK_THROW(boost::math::quantile(
 183          students_t_distribution<RealType>(1.),  // degrees_of_freedom (ignored).
 184          static_cast<RealType>(1)), std::overflow_error); // t == +infinity.
 185
 186   BOOST_CHECK_EQUAL(boost::math::quantile(
 187          students_t_distribution<RealType>(1.),  // degrees_of_freedom (ignored).
 188          static_cast<RealType>(0.5)),  //  probability == half - special case.
 189          static_cast<RealType>(0)); // t == zero.
 190
 191   BOOST_CHECK_EQUAL(boost::math::quantile(
 192          complement(
 193             students_t_distribution<RealType>(1.),  // degrees_of_freedom (ignored).
 194             static_cast<RealType>(0.5))),  //  probability == half - special case.
 195          static_cast<RealType>(0)); // t == zero.
 196
 197   BOOST_CHECK_CLOSE(boost::math::quantile(
 198          students_t_distribution<RealType>(1.),  // degrees_of_freedom (ignored).
 199          static_cast<RealType>(0.5)),  //  probability == half - special case.
 200          static_cast<RealType>(0), // t == zero.
 201          tolerance);
 202
 203    BOOST_CHECK_CLOSE( // Tests of p middling.
 204       ::boost::math::cdf(
 205          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 206          static_cast<RealType>(-0.559429644)),  // t
 207          static_cast<RealType>(0.3), // probability.
 208          tolerance);
 209    BOOST_CHECK_CLOSE(
 210       ::boost::math::quantile(
 211          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 212          static_cast<RealType>(0.3)),  // probability.
 213          static_cast<RealType>(-0.559429644), // t
 214          tolerance);
 215
 216    BOOST_CHECK_CLOSE(
 217       ::boost::math::quantile(
 218          complement(
 219             students_t_distribution<RealType>(5.),  // degrees_of_freedom
 220             static_cast<RealType>(0.7))),  // probability.
 221          static_cast<RealType>(-0.559429644), // t
 222          tolerance);
 223
 224    BOOST_CHECK_CLOSE( // Tests of p high.
 225       ::boost::math::cdf(
 226          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 227          static_cast<RealType>(1.475884049)),  // t
 228          static_cast<RealType>(0.9), // probability.
 229          tolerance);
 230    BOOST_CHECK_CLOSE(
 231       ::boost::math::quantile(
 232          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 233          static_cast<RealType>(0.9)),  // probability.
 234          static_cast<RealType>(1.475884049), // t
 235          tolerance);
 236
 237    BOOST_CHECK_CLOSE( // Tests of p low.
 238       ::boost::math::cdf(
 239          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 240          static_cast<RealType>(-1.475884049)),  // t
 241          static_cast<RealType>(0.1), // probability.
 242          tolerance);
 243    BOOST_CHECK_CLOSE(
 244       ::boost::math::quantile(
 245          students_t_distribution<RealType>(5.),  // degrees_of_freedom
 246          static_cast<RealType>(0.1)),  // probability.
 247          static_cast<RealType>(-1.475884049), // t
 248          tolerance);
 249
 250    BOOST_CHECK_CLOSE(
 251       ::boost::math::cdf(
 252          students_t_distribution<RealType>(2.),  // degrees_of_freedom
 253          static_cast<RealType>(-6.96455673428326)),  // t
 254          static_cast<RealType>(0.01), // probability.
 255          tolerance);
 256
 257    BOOST_CHECK_CLOSE(
 258       ::boost::math::quantile(
 259          students_t_distribution<RealType>(2.),  // degrees_of_freedom
 260          static_cast<RealType>(0.01)),  // probability.
 261          static_cast<RealType>(-6.96455673428326), // t
 262          tolerance);
 263
 264       //
 265       // Some special tests to exercise the double-precision approximations
 266       // to the quantile:
 267       //
 268       // tolerance is 50 eps expressed as a persent:
 269       //
 270       tolerance = boost::math::tools::epsilon<RealType>() * 5000;
 271       BOOST_CHECK_CLOSE(boost::math::quantile(
 272          students_t_distribution<RealType>(2.00390625L),                     // degrees_of_freedom.
 273          static_cast<RealType>(0.5625L)),                                    //  probability.
 274          static_cast<RealType>(0.178133131573788108465134803511798566L),     // t.
 275          tolerance);
 276       BOOST_CHECK_CLOSE(boost::math::quantile(
 277          students_t_distribution<RealType>(1L),                              // degrees_of_freedom.
 278          static_cast<RealType>(0.03125L)),                                   //  probability.
 279          static_cast<RealType>(-10.1531703876088604621071476634194722L),     // t.
 280          tolerance);
 281       BOOST_CHECK_CLOSE(boost::math::quantile(
 282          students_t_distribution<RealType>(1L),                            // degrees_of_freedom.
 283          static_cast<RealType>(0.875L)),                                   //  probability.
 284          static_cast<RealType>(2.41421356237309504880168872421390942L),    // t.
 285          tolerance);
 286       BOOST_CHECK_CLOSE(boost::math::quantile(
 287          students_t_distribution<RealType>(2L),                              // degrees_of_freedom.
 288          static_cast<RealType>(0.03125L)),                                   //  probability.
 289          static_cast<RealType>(-3.81000381000571500952501666878143315L),     // t.
 290          tolerance);
 291       BOOST_CHECK_CLOSE(boost::math::quantile(
 292          students_t_distribution<RealType>(2L),                            // degrees_of_freedom.
 293          static_cast<RealType>(0.875L)),                                   //  probability.
 294          static_cast<RealType>(1.60356745147454630810732088527854144L),    // t.
 295          tolerance);
 296       BOOST_CHECK_CLOSE(boost::math::quantile(
 297          students_t_distribution<RealType>(4L),                              // degrees_of_freedom.
 298          static_cast<RealType>(0.03125L)),                                   //  probability.
 299          static_cast<RealType>(-2.56208431914409044861223047927635034L),     // t.
 300          tolerance);
 301       BOOST_CHECK_CLOSE(boost::math::quantile(
 302          students_t_distribution<RealType>(4L),                            // degrees_of_freedom.
 303          static_cast<RealType>(0.875L)),                                   //  probability.
 304          static_cast<RealType>(1.34439755550909142430681981315923574L),    // t.
 305          tolerance);
 306       BOOST_CHECK_CLOSE(boost::math::quantile(
 307          students_t_distribution<RealType>(6L),                              // degrees_of_freedom.
 308          static_cast<RealType>(0.03125L)),                                   //  probability.
 309          static_cast<RealType>(-2.28348667906973065861212495010082952L),     // t.
 310          tolerance);
 311       BOOST_CHECK_CLOSE(boost::math::quantile(
 312          students_t_distribution<RealType>(6L),                            // degrees_of_freedom.
 313          static_cast<RealType>(0.875L)),                                   //  probability.
 314          static_cast<RealType>(1.27334930914664286821103236660071906L),    // t.
 315          tolerance);
 316       BOOST_CHECK_CLOSE(boost::math::quantile(
 317          students_t_distribution<RealType>(8L),                              // degrees_of_freedom.
 318          static_cast<RealType>(0.03125L)),                                   //  probability.
 319          static_cast<RealType>(-2.16296475406014719458642055768894376L),     // t.
 320          tolerance);
 321       BOOST_CHECK_CLOSE(boost::math::quantile(
 322          students_t_distribution<RealType>(8L),                            // degrees_of_freedom.
 323          static_cast<RealType>(0.875L)),                                   //  probability.
 324          static_cast<RealType>(1.24031826078267310637634677726479038L),    // t.
 325          tolerance);
 326       BOOST_CHECK_CLOSE(boost::math::quantile(
 327          students_t_distribution<RealType>(10L),                             // degrees_of_freedom.
 328          static_cast<RealType>(0.03125L)),                                   //  probability.
 329          static_cast<RealType>(-2.09596136475109350926340169211429572L),     // t.
 330          tolerance);
 331       BOOST_CHECK_CLOSE(boost::math::quantile(
 332          students_t_distribution<RealType>(10L),                         // degrees_of_freedom.
 333          static_cast<RealType>(0.875L)),                                 //  probability.
 334          static_cast<RealType>(1.2212553950039221407185188573696834L),   // t.
 335          tolerance);
 336       BOOST_CHECK_CLOSE(boost::math::quantile(
 337          students_t_distribution<RealType>(2.125L),                          // degrees_of_freedom.
 338          static_cast<RealType>(0.03125L)),                                   //  probability.
 339          static_cast<RealType>(-3.62246031671091980110493455859296532L),     // t.
 340          tolerance);
 341       BOOST_CHECK_CLOSE(boost::math::quantile(
 342          students_t_distribution<RealType>(2.125L),                        // degrees_of_freedom.
 343          static_cast<RealType>(0.875L)),                                   //  probability.
 344          static_cast<RealType>(1.56905270993307293450392958697861969L),    // t.
 345          tolerance);
 346       BOOST_CHECK_CLOSE(boost::math::quantile(
 347          students_t_distribution<RealType>(3L),                              // degrees_of_freedom.
 348          static_cast<RealType>(0.03125L)),                                   //  probability.
 349          static_cast<RealType>(-2.90004411882995814036141778367917946L),     // t.
 350          tolerance);
 351       BOOST_CHECK_CLOSE(boost::math::quantile(
 352          students_t_distribution<RealType>(3L),                            // degrees_of_freedom.
 353          static_cast<RealType>(0.875L)),                                   //  probability.
 354          static_cast<RealType>(1.42262528146180931868169289781115099L),    // t.
 355          tolerance);
 356
 357       if(boost::is_floating_point<RealType>::value)
 358       {
 359          BOOST_CHECK_CLOSE(boost::math::cdf(
 360             students_t_distribution<RealType>(1e30f),
 361                boost::math::quantile(
 362                   students_t_distribution<RealType>(1e30f), static_cast<RealType>(0.25f))),
 363             static_cast<RealType>(0.25f), tolerance);
 364          BOOST_CHECK_CLOSE(boost::math::cdf(
 365             students_t_distribution<RealType>(1e20f),
 366                boost::math::quantile(
 367                   students_t_distribution<RealType>(1e20f), static_cast<RealType>(0.25f))),
 368             static_cast<RealType>(0.25f), tolerance);
 369          BOOST_CHECK_CLOSE(boost::math::cdf(
 370             students_t_distribution<RealType>(static_cast<RealType>(0x7FFFFFFF)),
 371                boost::math::quantile(
 372                   students_t_distribution<RealType>(static_cast<RealType>(0x7FFFFFFF)), static_cast<RealType>(0.25f))),
 373             static_cast<RealType>(0.25f), tolerance);
 374          BOOST_CHECK_CLOSE(boost::math::cdf(
 375             students_t_distribution<RealType>(static_cast<RealType>(0x10000000)),
 376                boost::math::quantile(
 377                   students_t_distribution<RealType>(static_cast<RealType>(0x10000000)), static_cast<RealType>(0.25f))),
 378             static_cast<RealType>(0.25f), tolerance);
 379          BOOST_CHECK_CLOSE(boost::math::cdf(
 380             students_t_distribution<RealType>(static_cast<RealType>(0x0fffffff)),
 381                boost::math::quantile(
 382                   students_t_distribution<RealType>(static_cast<RealType>(0x0fffffff)), static_cast<RealType>(0.25f))),
 383             static_cast<RealType>(0.25f), tolerance);
 384       }
 385
 386   // Student's t pdf tests.
 387   // for PDF checks, use 100 eps tolerance expressed as a percent:
 388    tolerance = boost::math::tools::epsilon<RealType>() * 10000;
 389
 390    for(unsigned i = 1; i < 20; i += 3)
 391    {
 392       for(RealType r = -10; r < 10; r += 0.125)
 393       {
 394          //std::cout << "df=" << i << " t=" << r << std::endl;
 395          BOOST_CHECK_CLOSE(
 396             boost::math::pdf(
 397                students_t_distribution<RealType>(static_cast<RealType>(i)),
 398                r),
 399             naive_pdf<RealType>(static_cast<RealType>(i), r),
 400             tolerance);
 401       }
 402    }
 403
 404     RealType tol2 = boost::math::tools::epsilon<RealType>() * 5;
 405     students_t_distribution<RealType> dist(8);
 406     RealType x = static_cast<RealType>(0.125);
 407     using namespace std; // ADL of std names.
 408     // mean:
 409     BOOST_CHECK_CLOSE(
 410        mean(dist)
 411        , static_cast<RealType>(0), tol2);
 412     // variance:
 413  //   BOOST_CHECK_CLOSE(
 414  //      variance(dist)
 415  //      , static_cast<RealType>(13.0L / 6.0L), tol2);
 416  //// was     , static_cast<RealType>(8.0L / 6.0L), tol2);
 417     // std deviation:
 418     BOOST_CHECK_CLOSE(
 419        standard_deviation(dist)
 420        , static_cast<RealType>(sqrt(8.0L / 6.0L)), tol2);
 421     // hazard:
 422     BOOST_CHECK_CLOSE(
 423        hazard(dist, x)
 424        , pdf(dist, x) / cdf(complement(dist, x)), tol2);
 425     // cumulative hazard:
 426     BOOST_CHECK_CLOSE(
 427        chf(dist, x)
 428        , -log(cdf(complement(dist, x))), tol2);
 429     // coefficient_of_variation:
 430     BOOST_MATH_CHECK_THROW(
 431        coefficient_of_variation(dist),
 432        std::overflow_error);
 433     // mode:
 434     BOOST_CHECK_CLOSE(
 435        mean(dist)
 436        , static_cast<RealType>(0), tol2);
 437     // median:
 438     BOOST_CHECK_CLOSE(
 439        median(dist)
 440        , static_cast<RealType>(0), tol2);
 441     // skewness:
 442     BOOST_CHECK_CLOSE(
 443        skewness(dist)
 444        , static_cast<RealType>(0), tol2);
 445     // kurtosis:
 446     BOOST_CHECK_CLOSE(
 447        kurtosis(dist)
 448        , static_cast<RealType>(4.5), tol2);
 449     // kurtosis excess:
 450     BOOST_CHECK_CLOSE(
 451        kurtosis_excess(dist)
 452        , static_cast<RealType>(1.5), tol2);
 453
 454     // Parameter estimation. These results are close to but
 455     // not identical to those reported on the NIST website at
 456     // http://www.itl.nist.gov/div898/handbook/prc/section2/prc222.htm
 457     // the NIST results appear to be calculated using a normal
 458     // approximation, which slightly under-estimates the degrees of
 459     // freedom required, particularly when the result is small.
 460     //
 461     BOOST_CHECK_EQUAL(
 462        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 463          static_cast<RealType>(0.5),
 464          static_cast<RealType>(0.005),
 465          static_cast<RealType>(0.01),
 466          static_cast<RealType>(1.0))),
 467          99);
 468     BOOST_CHECK_EQUAL(
 469        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 470          static_cast<RealType>(1.5),
 471          static_cast<RealType>(0.005),
 472          static_cast<RealType>(0.01),
 473          static_cast<RealType>(1.0))),
 474          14);
 475     BOOST_CHECK_EQUAL(
 476        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 477          static_cast<RealType>(0.5),
 478          static_cast<RealType>(0.025),
 479          static_cast<RealType>(0.01),
 480          static_cast<RealType>(1.0))),
 481          76);
 482     BOOST_CHECK_EQUAL(
 483        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 484          static_cast<RealType>(1.5),
 485          static_cast<RealType>(0.025),
 486          static_cast<RealType>(0.01),
 487          static_cast<RealType>(1.0))),
 488          11);
 489     BOOST_CHECK_EQUAL(
 490        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 491          static_cast<RealType>(0.5),
 492          static_cast<RealType>(0.05),
 493          static_cast<RealType>(0.01),
 494          static_cast<RealType>(1.0))),
 495          65);
 496     BOOST_CHECK_EQUAL(
 497        ceil(students_t_distribution<RealType>::find_degrees_of_freedom(
 498          static_cast<RealType>(1.5),
 499          static_cast<RealType>(0.05),
 500          static_cast<RealType>(0.01),
 501          static_cast<RealType>(1.0))),
 502          9);
 503
 504     // Test for large degrees of freedom when should be same as normal.
 505     RealType inf = std::numeric_limits<RealType>::infinity();
 506     RealType nan = std::numeric_limits<RealType>::quiet_NaN();
 507
 508     std::string type = typeid(RealType).name();
 509 //    if (type != "class boost::math::concepts::real_concept") fails for gcc
 510     if (typeid(RealType) != typeid(boost::math::concepts::real_concept))
 511     { // Ordinary floats only.
 512       RealType limit = 1/ boost::math::tools::epsilon<RealType>();
 513       // Default policy to get full accuracy.
 514       // std::cout << "Switch over to normal if df > " << limit << std::endl;
 515       // float Switch over to normal if df > 8.38861e+006
 516       // double Switch over to normal if df > 4.5036e+015
 517       // Can't test real_concept - doesn't converge.
 518
 519       boost::math::normal_distribution<RealType> n(0, 1); //
 520       students_t_distribution<RealType> st(boost::math::tools::max_value<RealType>()); // Well over the switchover point,
 521       // PDF
 522       BOOST_CHECK_EQUAL(pdf(st, 0), pdf(n, 0.)); // Should be exactly equal.
 523
 524       students_t_distribution<RealType> st2(limit /5 ); // Just below the switchover point,
 525       BOOST_CHECK_CLOSE_FRACTION(pdf(st2, 0), pdf(n, 0.), tolerance); // Should be very close to normal.
 526       // CDF
 527       BOOST_CHECK_EQUAL(cdf(st, 0), cdf(n, 0.)); // Should be exactly equal.
 528       BOOST_CHECK_CLOSE_FRACTION(cdf(st2, 0), cdf(n, 0.), tolerance); // Should be very close to normal.
 529
 530       // Tests for df = infinity.
 531       students_t_distribution<RealType> infdf(inf);
 532       BOOST_CHECK_EQUAL(infdf.degrees_of_freedom(), inf);
 533       BOOST_CHECK_EQUAL(mean(infdf), 0); // OK.
 534 #ifndef BOOST_NO_EXCEPTIONS
 535       BOOST_MATH_CHECK_THROW(students_t_distribution<RealType> minfdf(-inf), std::domain_error);
 536       BOOST_MATH_CHECK_THROW(students_t_distribution<RealType> minfdf(nan), std::domain_error);
 537       BOOST_MATH_CHECK_THROW(students_t_distribution<RealType> minfdf(-nan), std::domain_error);
 538 #endif
 539       BOOST_CHECK_EQUAL(pdf(infdf, -inf), 0);
 540       BOOST_CHECK_EQUAL(pdf(infdf, +inf), 0);
 541       BOOST_CHECK_EQUAL(cdf(infdf, -inf), 0);
 542       BOOST_CHECK_EQUAL(cdf(infdf, +inf), 1);
 543
 544      // BOOST_CHECK_CLOSE_FRACTION(pdf(infdf, 0), static_cast<RealType>(0.3989422804014326779399460599343818684759L), tolerance);
 545       BOOST_CHECK_CLOSE_FRACTION(pdf(infdf, 0),boost::math::constants::one_div_root_two_pi<RealType>() , tolerance);
 546       BOOST_CHECK_CLOSE_FRACTION(cdf(infdf, 0),boost::math::constants::half<RealType>() , tolerance);
 547
 548     // Checks added for Trac #7717 report by Thomas Mang.
 549
 550     BOOST_MATH_CHECK_THROW(quantile(dist, -1), std::domain_error);
 551     BOOST_MATH_CHECK_THROW(quantile(dist, 2), std::domain_error);
 552     BOOST_MATH_CHECK_THROW(pdf(students_t_distribution<RealType>(0), 0), std::domain_error);
 553     BOOST_MATH_CHECK_THROW(pdf(students_t_distribution<RealType>(-1), 0), std::domain_error);
 554
 555     // Check on df for mean (moment k = 1)
 556     BOOST_MATH_CHECK_THROW(mean(students_t_distribution<RealType>(nan)), std::domain_error);
 557 //    BOOST_MATH_CHECK_THROW(mean(students_t_distribution<RealType>(inf)), std::domain_error); inf is now OK
 558     BOOST_MATH_CHECK_THROW(mean(students_t_distribution<RealType>(-1)), std::domain_error);
 559     BOOST_MATH_CHECK_THROW(mean(students_t_distribution<RealType>(0)), std::domain_error);
 560     BOOST_MATH_CHECK_THROW(mean(students_t_distribution<RealType>(1)), std::domain_error); // df == k
 561     BOOST_CHECK_EQUAL(mean(students_t_distribution<RealType>(2)), 0); // OK.
 562     BOOST_CHECK_EQUAL(mean(students_t_distribution<RealType>(inf)), 0); // OK.
 563
 564     // Check on df for variance (moment 2)
 565     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(nan)), std::domain_error);
 566 //    BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(inf)), std::domain_error); // inf is now OK.
 567     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(-1)), std::domain_error);
 568     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(0)), std::domain_error);
 569     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(1)), std::domain_error);
 570     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(static_cast<RealType>(1.99999L))), std::domain_error);
 571     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(static_cast<RealType>(1.99999L))), std::domain_error);
 572     BOOST_MATH_CHECK_THROW(variance(students_t_distribution<RealType>(2)), std::domain_error); // df ==
 573     BOOST_CHECK_EQUAL(variance(students_t_distribution<RealType>(2.5)), 5); // OK.
 574     BOOST_CHECK_EQUAL(variance(students_t_distribution<RealType>(3)), 3); // OK.
 575     BOOST_CHECK_EQUAL(variance(students_t_distribution<RealType>(inf)), 1); // OK.
 576
 577     // Check on df for skewness (moment 3)
 578     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(nan)), std::domain_error);
 579     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(-1)), std::domain_error);
 580     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(0)), std::domain_error);
 581     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(1)), std::domain_error);
 582     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(1.5L)), std::domain_error);
 583     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(2)), std::domain_error);
 584     BOOST_MATH_CHECK_THROW(skewness(students_t_distribution<RealType>(3)), std::domain_error); // df == k
 585     BOOST_CHECK_EQUAL(skewness(students_t_distribution<RealType>(3.5)), 0); // OK.
 586     BOOST_CHECK_EQUAL(skewness(students_t_distribution<RealType>(4)), 0); // OK.
 587     BOOST_CHECK_EQUAL(skewness(students_t_distribution<RealType>(inf)), 0); // OK.
 588
 589     // Check on df for kurtosis_excess (moment 4)
 590     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(nan)), std::domain_error);
 591     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(-1)), std::domain_error);
 592     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(0)), std::domain_error);
 593     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(1)), std::domain_error);
 594     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(1.5L)), std::domain_error);
 595     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(2)), std::domain_error);
 596     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(static_cast<RealType>(2.1))), std::domain_error);
 597     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(3)), std::domain_error);
 598     BOOST_MATH_CHECK_THROW(kurtosis_excess(students_t_distribution<RealType>(4)), std::domain_error); // df == k
 599     BOOST_CHECK_EQUAL(kurtosis_excess(students_t_distribution<RealType>(5)), 6); // OK.
 600     BOOST_CHECK_EQUAL(kurtosis_excess(students_t_distribution<RealType>(inf)), 0); // OK.
 601
 602     // Check on df for kurtosis (moment 4)
 603     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(nan)), std::domain_error);
 604     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(-1)), std::domain_error);
 605     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(0)), std::domain_error);
 606     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(1)), std::domain_error);
 607     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(2)), std::domain_error);
 608     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(static_cast<RealType>(2.0001L))), std::domain_error);
 609     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(3)), std::domain_error);
 610     BOOST_MATH_CHECK_THROW(kurtosis(students_t_distribution<RealType>(4)), std::domain_error); // df == k
 611     BOOST_CHECK_EQUAL(kurtosis(students_t_distribution<RealType>(5)), 9); // OK.
 612     BOOST_CHECK_EQUAL(kurtosis(students_t_distribution<RealType>(inf)), 3); // OK.
 613
 614    }
 615
 616
 617     // Use a new distribution ignore_error_students_t with a custom policy to ignore all errors,
 618     // and check returned values are as expected.
 619
 620     /*
 621      Sandia-darwin-intel-12.0 - math - test_students_t / intel-darwin-12.0
 622     ../libs/math/test/test_students_t.cpp(544): error: "domain_error" has already been declared in the current scope
 623     using boost::math::policies::domain_error;
 624
 625 ../libs/math/test/test_students_t.cpp(552): error: "pole_error" has already been declared in the current scope
 626     using boost::math::policies::pole_error;
 627
 628     Unclear where previous declaration is.
 629     Does not seem to be in student_t.hpp or any included files???
 630
 631     So to avoid this perceived problem by this compiler,
 632     the ignore policy below uses fully specified names.
 633     */
 634
 635     using boost::math::policies::policy;
 636   // Types of error whose action can be altered by policies:.
 637   //using boost::math::policies::evaluation_error;
 638   //using boost::math::policies::domain_error;
 639   //using boost::math::policies::overflow_error;
 640   //using boost::math::policies::underflow_error;
 641   //using boost::math::policies::domain_error;
 642   //using boost::math::policies::pole_error;
 643
 644   //// Actions on error (in enum error_policy_type):
 645   //using boost::math::policies::errno_on_error;
 646   //using boost::math::policies::ignore_error;
 647   //using boost::math::policies::throw_on_error;
 648   //using boost::math::policies::denorm_error;
 649   //using boost::math::policies::pole_error;
 650   //using boost::math::policies::user_error;
 651
 652   typedef policy<
 653     boost::math::policies::domain_error<boost::math::policies::ignore_error>,
 654     boost::math::policies::overflow_error<boost::math::policies::ignore_error>,
 655     boost::math::policies::underflow_error<boost::math::policies::ignore_error>,
 656     boost::math::policies::denorm_error<boost::math::policies::ignore_error>,
 657     boost::math::policies::pole_error<boost::math::policies::ignore_error>,
 658     boost::math::policies::evaluation_error<boost::math::policies::ignore_error>
 659               > my_ignore_policy;
 660
 661   typedef students_t_distribution<RealType, my_ignore_policy> ignore_error_students_t;
 662
 663
 664
 665   // Only test NaN and infinity if type has these features (realconcept returns zero).
 666   // Integers are always converted to RealType,
 667   // others requires static cast to RealType from long double.
 668
 669   if(std::numeric_limits<RealType>::has_quiet_NaN)
 670   {
 671   // Mean
 672     BOOST_CHECK((boost::math::isnan)(mean(ignore_error_students_t(-1))));
 673     BOOST_CHECK((boost::math::isnan)(mean(ignore_error_students_t(0))));
 674     BOOST_CHECK((boost::math::isnan)(mean(ignore_error_students_t(1))));
 675
 676     // Variance
 677     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(std::numeric_limits<RealType>::quiet_NaN()))));
 678     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(-1))));
 679     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(0))));
 680     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(1))));
 681     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(static_cast<RealType>(1.7L)))));
 682     BOOST_CHECK((boost::math::isnan)(variance(ignore_error_students_t(2))));
 683
 684   // Skewness
 685     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(std::numeric_limits<RealType>::quiet_NaN()))));
 686     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(-1))));
 687     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(0))));
 688     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(1))));
 689     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(2))));
 690     BOOST_CHECK((boost::math::isnan)(skewness(ignore_error_students_t(3))));
 691
 692   // Kurtosis
 693     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(std::numeric_limits<RealType>::quiet_NaN()))));
 694     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(-1))));
 695     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(0))));
 696     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(1))));
 697     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(2))));
 698     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(static_cast<RealType>(2.0001L)))));
 699     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(3))));
 700     BOOST_CHECK((boost::math::isnan)(kurtosis(ignore_error_students_t(4))));
 701
 702     // Kurtosis excess
 703     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(std::numeric_limits<RealType>::quiet_NaN()))));
 704     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(-1))));
 705     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(0))));
 706     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(1))));
 707     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(2))));
 708     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(static_cast<RealType>(2.0001L)))));
 709     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(3))));
 710     BOOST_CHECK((boost::math::isnan)(kurtosis_excess(ignore_error_students_t(4))));
 711   } // has_quiet_NaN
 712
 713   BOOST_CHECK(boost::math::isfinite(mean(ignore_error_students_t(1 + std::numeric_limits<RealType>::epsilon()))));
 714   BOOST_CHECK(boost::math::isfinite(variance(ignore_error_students_t(2 + 2 * std::numeric_limits<RealType>::epsilon()))));
 715   BOOST_CHECK(boost::math::isfinite(variance(ignore_error_students_t(static_cast<RealType>(2.0001L)))));
 716   BOOST_CHECK(boost::math::isfinite(variance(ignore_error_students_t(2 + 2 * std::numeric_limits<RealType>::epsilon()))));
 717   BOOST_CHECK(boost::math::isfinite(skewness(ignore_error_students_t(3 + 3 * std::numeric_limits<RealType>::epsilon()))));
 718   BOOST_CHECK(boost::math::isfinite(kurtosis(ignore_error_students_t(4 + 4 * std::numeric_limits<RealType>::epsilon()))));
 719   BOOST_CHECK(boost::math::isfinite(kurtosis(ignore_error_students_t(static_cast<RealType>(4.0001L)))));
 720
 721   // check_out_of_range<students_t_distribution<RealType> >(1);
 722   // Cannot be used because fails "exception std::domain_error is expected but not raised"
 723   // if df = +infinity is allowed, must use new version that allows skipping infinity tests.
 724   // Infinite == true
 725
 726   check_support<students_t_distribution<RealType> >(students_t_distribution<RealType>(1), true);
 727
 728 } // template <class RealType>void test_spots(RealType)
 729
 730 BOOST_AUTO_TEST_CASE( test_main )
 731 {
 732   // Check that can construct students_t distribution using the two convenience methods:
 733   using namespace boost::math;
 734   students_t myst1(2); // Using typedef
 735   students_t_distribution<> myst2(2); // Using default RealType double.
 736    //students_t_distribution<double> myst3(2); // Using explicit RealType double.
 737
 738    // Basic sanity-check spot values.
 739    // (Parameter value, arbitrarily zero, only communicates the floating point type).
 740   test_spots(0.0F); // Test float. OK at decdigits = 0 tolerance = 0.0001 %
 741   test_spots(0.0); // Test double. OK at decdigits 7, tolerance = 1e07 %
 742 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
 743   test_spots(0.0L); // Test long double.
 744 #if !BOOST_WORKAROUND(__BORLANDC__, BOOST_TESTED_AT(0x582))
 745   test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
 746 #endif
 747 #else
 748    std::cout << "<note>The long double tests have been disabled on this platform "
 749       "either because the long double overloads of the usual math functions are "
 750       "not available at all, or because they are too inaccurate for these tests "
 751       "to pass.</note>" << std::endl;
 752 #endif
 753
 754
 755
 756 } // BOOST_AUTO_TEST_CASE( test_main )
 757
 758 /*
 759
 760 Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\test_students_t.exe"
 761 Running 1 test case...
 762 Tolerance for type float is 0.0001 %
 763 Tolerance for type double is 0.0001 %
 764 Tolerance for type long double is 0.0001 %
 765 Tolerance for type class boost::math::concepts::real_concept is 0.0001 %
 766 *** No errors detected
 767
 768 */
 769
 770