students_t_two_samples.cpp 9.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245
  1. // Copyright John Maddock 2006.
  2. // Copyright Paul A. Bristow 2007, 2010
  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. #ifdef _MSC_VER
  8. # pragma warning(disable: 4512) // assignment operator could not be generated.
  9. # pragma warning(disable: 4510) // default constructor could not be generated.
  10. # pragma warning(disable: 4610) // can never be instantiated - user defined constructor required.
  11. #endif
  12. #include <iostream>
  13. using std::cout; using std::endl;
  14. using std::left; using std::fixed; using std::right; using std::scientific;
  15. #include <iomanip>
  16. using std::setw;
  17. using std::setprecision;
  18. #include <boost/math/distributions/students_t.hpp>
  19. using boost::math::students_t;
  20. void two_samples_t_test_equal_sd(
  21. double Sm1, // Sm1 = Sample Mean 1.
  22. double Sd1, // Sd1 = Sample Standard Deviation 1.
  23. unsigned Sn1, // Sn1 = Sample Size 1.
  24. double Sm2, // Sm2 = Sample Mean 2.
  25. double Sd2, // Sd2 = Sample Standard Deviation 2.
  26. unsigned Sn2, // Sn2 = Sample Size 2.
  27. double alpha) // alpha = Significance Level.
  28. {
  29. // A Students t test applied to two sets of data.
  30. // We are testing the null hypothesis that the two
  31. // samples have the same mean and that any difference
  32. // if due to chance.
  33. // See http://www.itl.nist.gov/div898/handbook/eda/section3/eda353.htm
  34. //
  35. using namespace std;
  36. // using namespace boost::math;
  37. using boost::math::students_t;
  38. // Print header:
  39. cout <<
  40. "_______________________________________________\n"
  41. "Student t test for two samples (equal variances)\n"
  42. "_______________________________________________\n\n";
  43. cout << setprecision(5);
  44. cout << setw(55) << left << "Number of Observations (Sample 1)" << "= " << Sn1 << "\n";
  45. cout << setw(55) << left << "Sample 1 Mean" << "= " << Sm1 << "\n";
  46. cout << setw(55) << left << "Sample 1 Standard Deviation" << "= " << Sd1 << "\n";
  47. cout << setw(55) << left << "Number of Observations (Sample 2)" << "= " << Sn2 << "\n";
  48. cout << setw(55) << left << "Sample 2 Mean" << "= " << Sm2 << "\n";
  49. cout << setw(55) << left << "Sample 2 Standard Deviation" << "= " << Sd2 << "\n";
  50. //
  51. // Now we can calculate and output some stats:
  52. //
  53. // Degrees of freedom:
  54. double v = Sn1 + Sn2 - 2;
  55. cout << setw(55) << left << "Degrees of Freedom" << "= " << v << "\n";
  56. // Pooled variance and hence standard deviation:
  57. double sp = sqrt(((Sn1-1) * Sd1 * Sd1 + (Sn2-1) * Sd2 * Sd2) / v);
  58. cout << setw(55) << left << "Pooled Standard Deviation" << "= " << sp << "\n";
  59. // t-statistic:
  60. double t_stat = (Sm1 - Sm2) / (sp * sqrt(1.0 / Sn1 + 1.0 / Sn2));
  61. cout << setw(55) << left << "T Statistic" << "= " << t_stat << "\n";
  62. //
  63. // Define our distribution, and get the probability:
  64. //
  65. students_t dist(v);
  66. double q = cdf(complement(dist, fabs(t_stat)));
  67. cout << setw(55) << left << "Probability that difference is due to chance" << "= "
  68. << setprecision(3) << scientific << 2 * q << "\n\n";
  69. //
  70. // Finally print out results of alternative hypothesis:
  71. //
  72. cout << setw(55) << left <<
  73. "Results for Alternative Hypothesis and alpha" << "= "
  74. << setprecision(4) << fixed << alpha << "\n\n";
  75. cout << "Alternative Hypothesis Conclusion\n";
  76. cout << "Sample 1 Mean != Sample 2 Mean " ;
  77. if(q < alpha / 2)
  78. cout << "NOT REJECTED\n";
  79. else
  80. cout << "REJECTED\n";
  81. cout << "Sample 1 Mean < Sample 2 Mean ";
  82. if(cdf(dist, t_stat) < alpha)
  83. cout << "NOT REJECTED\n";
  84. else
  85. cout << "REJECTED\n";
  86. cout << "Sample 1 Mean > Sample 2 Mean ";
  87. if(cdf(complement(dist, t_stat)) < alpha)
  88. cout << "NOT REJECTED\n";
  89. else
  90. cout << "REJECTED\n";
  91. cout << endl << endl;
  92. }
  93. void two_samples_t_test_unequal_sd(
  94. double Sm1, // Sm1 = Sample Mean 1.
  95. double Sd1, // Sd1 = Sample Standard Deviation 1.
  96. unsigned Sn1, // Sn1 = Sample Size 1.
  97. double Sm2, // Sm2 = Sample Mean 2.
  98. double Sd2, // Sd2 = Sample Standard Deviation 2.
  99. unsigned Sn2, // Sn2 = Sample Size 2.
  100. double alpha) // alpha = Significance Level.
  101. {
  102. // A Students t test applied to two sets of data.
  103. // We are testing the null hypothesis that the two
  104. // samples have the same mean and
  105. // that any difference is due to chance.
  106. // See http://www.itl.nist.gov/div898/handbook/eda/section3/eda353.htm
  107. //
  108. using namespace std;
  109. using boost::math::students_t;
  110. // Print header:
  111. cout <<
  112. "_________________________________________________\n"
  113. "Student t test for two samples (unequal variances)\n"
  114. "_________________________________________________\n\n";
  115. cout << setprecision(5);
  116. cout << setw(55) << left << "Number of Observations (Sample 1)" << "= " << Sn1 << "\n";
  117. cout << setw(55) << left << "Sample 1 Mean" << "= " << Sm1 << "\n";
  118. cout << setw(55) << left << "Sample 1 Standard Deviation" << "= " << Sd1 << "\n";
  119. cout << setw(55) << left << "Number of Observations (Sample 2)" << "= " << Sn2 << "\n";
  120. cout << setw(55) << left << "Sample 2 Mean" << "= " << Sm2 << "\n";
  121. cout << setw(55) << left << "Sample 2 Standard Deviation" << "= " << Sd2 << "\n";
  122. //
  123. // Now we can calculate and output some stats:
  124. //
  125. // Degrees of freedom:
  126. double v = Sd1 * Sd1 / Sn1 + Sd2 * Sd2 / Sn2;
  127. v *= v;
  128. double t1 = Sd1 * Sd1 / Sn1;
  129. t1 *= t1;
  130. t1 /= (Sn1 - 1);
  131. double t2 = Sd2 * Sd2 / Sn2;
  132. t2 *= t2;
  133. t2 /= (Sn2 - 1);
  134. v /= (t1 + t2);
  135. cout << setw(55) << left << "Degrees of Freedom" << "= " << v << "\n";
  136. // t-statistic:
  137. double t_stat = (Sm1 - Sm2) / sqrt(Sd1 * Sd1 / Sn1 + Sd2 * Sd2 / Sn2);
  138. cout << setw(55) << left << "T Statistic" << "= " << t_stat << "\n";
  139. //
  140. // Define our distribution, and get the probability:
  141. //
  142. students_t dist(v);
  143. double q = cdf(complement(dist, fabs(t_stat)));
  144. cout << setw(55) << left << "Probability that difference is due to chance" << "= "
  145. << setprecision(3) << scientific << 2 * q << "\n\n";
  146. //
  147. // Finally print out results of alternative hypothesis:
  148. //
  149. cout << setw(55) << left <<
  150. "Results for Alternative Hypothesis and alpha" << "= "
  151. << setprecision(4) << fixed << alpha << "\n\n";
  152. cout << "Alternative Hypothesis Conclusion\n";
  153. cout << "Sample 1 Mean != Sample 2 Mean " ;
  154. if(q < alpha / 2)
  155. cout << "NOT REJECTED\n";
  156. else
  157. cout << "REJECTED\n";
  158. cout << "Sample 1 Mean < Sample 2 Mean ";
  159. if(cdf(dist, t_stat) < alpha)
  160. cout << "NOT REJECTED\n";
  161. else
  162. cout << "REJECTED\n";
  163. cout << "Sample 1 Mean > Sample 2 Mean ";
  164. if(cdf(complement(dist, t_stat)) < alpha)
  165. cout << "NOT REJECTED\n";
  166. else
  167. cout << "REJECTED\n";
  168. cout << endl << endl;
  169. }
  170. int main()
  171. {
  172. //
  173. // Run tests for Car Mileage sample data
  174. // http://www.itl.nist.gov/div898/handbook/eda/section3/eda3531.htm
  175. // from the NIST website http://www.itl.nist.gov. The data compares
  176. // miles per gallon of US cars with miles per gallon of Japanese cars.
  177. //
  178. two_samples_t_test_equal_sd(20.14458, 6.414700, 249, 30.48101, 6.107710, 79, 0.05);
  179. two_samples_t_test_unequal_sd(20.14458, 6.414700, 249, 30.48101, 6.107710, 79, 0.05);
  180. return 0;
  181. } // int main()
  182. /*
  183. Output is:
  184. _______________________________________________
  185. Student t test for two samples (equal variances)
  186. _______________________________________________
  187. Number of Observations (Sample 1) = 249
  188. Sample 1 Mean = 20.145
  189. Sample 1 Standard Deviation = 6.4147
  190. Number of Observations (Sample 2) = 79
  191. Sample 2 Mean = 30.481
  192. Sample 2 Standard Deviation = 6.1077
  193. Degrees of Freedom = 326
  194. Pooled Standard Deviation = 6.3426
  195. T Statistic = -12.621
  196. Probability that difference is due to chance = 5.273e-030
  197. Results for Alternative Hypothesis and alpha = 0.0500
  198. Alternative Hypothesis Conclusion
  199. Sample 1 Mean != Sample 2 Mean NOT REJECTED
  200. Sample 1 Mean < Sample 2 Mean NOT REJECTED
  201. Sample 1 Mean > Sample 2 Mean REJECTED
  202. _________________________________________________
  203. Student t test for two samples (unequal variances)
  204. _________________________________________________
  205. Number of Observations (Sample 1) = 249
  206. Sample 1 Mean = 20.14458
  207. Sample 1 Standard Deviation = 6.41470
  208. Number of Observations (Sample 2) = 79
  209. Sample 2 Mean = 30.48101
  210. Sample 2 Standard Deviation = 6.10771
  211. Degrees of Freedom = 136.87499
  212. T Statistic = -12.94627
  213. Probability that difference is due to chance = 1.571e-025
  214. Results for Alternative Hypothesis and alpha = 0.0500
  215. Alternative Hypothesis Conclusion
  216. Sample 1 Mean != Sample 2 Mean NOT REJECTED
  217. Sample 1 Mean < Sample 2 Mean NOT REJECTED
  218. Sample 1 Mean > Sample 2 Mean REJECTED
  219. */