students_t_example2.cpp 4.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126
  1. // students_t_example2.cpp
  2. // Copyright Paul A. Bristow 2006.
  3. // Use, modification and distribution are subject to the
  4. // Boost Software License, Version 1.0.
  5. // (See accompanying file LICENSE_1_0.txt
  6. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  7. // Example 2 of using Student's t
  8. // A general guide to Student's t is at
  9. // http://en.wikipedia.org/wiki/Student's_t-test
  10. // (and many other elementary and advanced statistics texts).
  11. // It says:
  12. // The t statistic was invented by William Sealy Gosset
  13. // for cheaply monitoring the quality of beer brews.
  14. // "Student" was his pen name.
  15. // Gosset was statistician for Guinness brewery in Dublin, Ireland,
  16. // hired due to Claude Guinness's innovative policy of recruiting the
  17. // best graduates from Oxford and Cambridge for applying biochemistry
  18. // and statistics to Guinness's industrial processes.
  19. // Gosset published the t test in Biometrika in 1908,
  20. // but was forced to use a pen name by his employer who regarded the fact
  21. // that they were using statistics as a trade secret.
  22. // In fact, Gosset's identity was unknown not only to fellow statisticians
  23. // but to his employer - the company insisted on the pseudonym
  24. // so that it could turn a blind eye to the breach of its rules.
  25. // The Students't distribution function is described at
  26. // http://en.wikipedia.org/wiki/Student%27s_t_distribution
  27. #include <boost/math/distributions/students_t.hpp>
  28. using boost::math::students_t; // Probability of students_t(df, t).
  29. #include <iostream>
  30. using std::cout;
  31. using std::endl;
  32. #include <iomanip>
  33. using std::setprecision;
  34. using std::setw;
  35. #include <cmath>
  36. using std::sqrt;
  37. // This example of a one-sided test is from:
  38. //
  39. // from Statistics for Analytical Chemistry, 3rd ed. (1994), pp 59-60
  40. // J. C. Miller and J. N. Miller, Ellis Horwood ISBN 0 13 0309907.
  41. // An acid-base titrimetric method has a significant indicator error and
  42. // thus tends to give results with a positive systematic error (+bias).
  43. // To test this an exactly 0.1 M solution of acid is used to titrate
  44. // 25.00 ml of exactly 0.1 M solution of alkali,
  45. // with the following results (ml):
  46. double reference = 25.00; // 'True' result.
  47. const int values = 6; // titrations.
  48. double data [values] = {25.06, 25.18, 24.87, 25.51, 25.34, 25.41};
  49. int main()
  50. {
  51. cout << "Example2 using Student's t function. ";
  52. #if defined(__FILE__) && defined(__TIMESTAMP__) && defined(_MSC_FULL_VER)
  53. cout << " " << __FILE__ << ' ' << __TIMESTAMP__ << ' '<< _MSC_FULL_VER;
  54. #endif
  55. cout << endl;
  56. double sum = 0.;
  57. for (int value = 0; value < values; value++)
  58. { // Echo data and calculate mean.
  59. sum += data[value];
  60. cout << setw(4) << value << ' ' << setw(14) << data[value] << endl;
  61. }
  62. double mean = sum /static_cast<double>(values);
  63. cout << "Mean = " << mean << endl; // 25.2283
  64. double sd = 0.;
  65. for (int value = 0; value < values; value++)
  66. { // Calculate standard deviation.
  67. sd +=(data[value] - mean) * (data[value] - mean);
  68. }
  69. int degrees_of_freedom = values - 1; // Use the n-1 formula.
  70. sd /= degrees_of_freedom; // == variance.
  71. sd= sqrt(sd);
  72. cout << "Standard deviation = " << sd<< endl; // = 0.238279
  73. double t = (mean - reference) * sqrt(static_cast<double>(values))/ sd; //
  74. cout << "Student's t = " << t << ", with " << degrees_of_freedom << " degrees of freedom." << endl; // = 2.34725
  75. cout << "Probability of positive bias is " << cdf(students_t(degrees_of_freedom), t) << "."<< endl; // = 0.967108.
  76. // A 1-sided test because only testing for a positive bias.
  77. // If > 0.95 then greater than 1 in 20 conventional (arbitrary) requirement.
  78. return 0;
  79. } // int main()
  80. /*
  81. Output is:
  82. ------ Build started: Project: students_t_example2, Configuration: Debug Win32 ------
  83. Compiling...
  84. students_t_example2.cpp
  85. Linking...
  86. Autorun "i:\boost-06-05-03-1300\libs\math\test\Math_test\debug\students_t_example2.exe"
  87. Example2 using Student's t function. ..\..\..\..\..\..\boost-sandbox\libs\math_functions\example\students_t_example2.cpp Sat Aug 12 16:55:59 2006 140050727
  88. 0 25.06
  89. 1 25.18
  90. 2 24.87
  91. 3 25.51
  92. 4 25.34
  93. 5 25.41
  94. Mean = 25.2283
  95. Standard deviation = 0.238279
  96. Student's t = 2.34725, with 5 degrees of freedom.
  97. Probability of positive bias is 0.967108.
  98. Build Time 0:03
  99. Build log was saved at "file://i:\boost-06-05-03-1300\libs\math\test\Math_test\students_t_example2\Debug\BuildLog.htm"
  100. students_t_example2 - 0 error(s), 0 warning(s)
  101. ========== Build: 1 succeeded, 0 failed, 0 up-to-date, 0 skipped ==========
  102. */