libs/math/example/students_t_single_sample.cpp
// Copyright John Maddock 2006 // Copyright Paul A. Bristow 2007, 2010 // 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) #ifdef _MSC_VER # pragma warning(disable: 4512) // assignment operator could not be generated. # pragma warning(disable: 4510) // default constructor could not be generated. # pragma warning(disable: 4610) // can never be instantiated - user defined constructor required. #endif #include <boost/math/distributions/students_t.hpp> // avoid "using namespace std;" and "using namespace boost::math;" // to avoid potential ambiguity with names in std random. #include <iostream> using std::cout; using std::endl; using std::left; using std::fixed; using std::right; using std::scientific; #include <iomanip> using std::setw; using std::setprecision; void confidence_limits_on_mean(double Sm, double Sd, unsigned Sn) { // // Sm = Sample Mean. // Sd = Sample Standard Deviation. // Sn = Sample Size. // // Calculate confidence intervals for the mean. // For example if we set the confidence limit to // 0.95, we know that if we repeat the sampling // 100 times, then we expect that the true mean // will be between out limits on 95 occations. // Note: this is not the same as saying a 95% // confidence interval means that there is a 95% // probability that the interval contains the true mean. // The interval computed from a given sample either // contains the true mean or it does not. // See http://www.itl.nist.gov/div898/handbook/eda/section3/eda352.htm using boost::math::students_t; // Print out general info: cout << "__________________________________\n" "2-Sided Confidence Limits For Mean\n" "__________________________________\n\n"; cout << setprecision(7); cout << setw(40) << left << "Number of Observations" << "= " << Sn << "\n"; cout << setw(40) << left << "Mean" << "= " << Sm << "\n"; cout << setw(40) << left << "Standard Deviation" << "= " << Sd << "\n"; // // Define a table of significance/risk levels: // double alpha[] = { 0.5, 0.25, 0.1, 0.05, 0.01, 0.001, 0.0001, 0.00001 }; // // Start by declaring the distribution we'll need: // students_t dist(Sn - 1); // // Print table header: // cout << "\n\n" "_______________________________________________________________\n" "Confidence T Interval Lower Upper\n" " Value (%) Value Width Limit Limit\n" "_______________________________________________________________\n"; // // Now print out the data for the table rows. // for(unsigned i = 0; i < sizeof(alpha)/sizeof(alpha[0]); ++i) { // Confidence value: cout << fixed << setprecision(3) << setw(10) << right << 100 * (1-alpha[i]); // calculate T: double T = quantile(complement(dist, alpha[i] / 2)); // Print T: cout << fixed << setprecision(3) << setw(10) << right << T; // Calculate width of interval (one sided): double w = T * Sd / sqrt(double(Sn)); // Print width: if(w < 0.01) cout << scientific << setprecision(3) << setw(17) << right << w; else cout << fixed << setprecision(3) << setw(17) << right << w; // Print Limits: cout << fixed << setprecision(5) << setw(15) << right << Sm - w; cout << fixed << setprecision(5) << setw(15) << right << Sm + w << endl; } cout << endl; } // void confidence_limits_on_mean void single_sample_t_test(double M, double Sm, double Sd, unsigned Sn, double alpha) { // // M = true mean. // Sm = Sample Mean. // Sd = Sample Standard Deviation. // Sn = Sample Size. // alpha = Significance Level. // // A Students t test applied to a single set of data. // We are testing the null hypothesis that the true // mean of the sample is M, and that any variation is down // to chance. We can also test the alternative hypothesis // that any difference is not down to chance. // See http://www.itl.nist.gov/div898/handbook/eda/section3/eda352.htm using boost::math::students_t; // Print header: cout << "__________________________________\n" "Student t test for a single sample\n" "__________________________________\n\n"; cout << setprecision(5); cout << setw(55) << left << "Number of Observations" << "= " << Sn << "\n"; cout << setw(55) << left << "Sample Mean" << "= " << Sm << "\n"; cout << setw(55) << left << "Sample Standard Deviation" << "= " << Sd << "\n"; cout << setw(55) << left << "Expected True Mean" << "= " << M << "\n\n"; // // Now we can calculate and output some stats: // // Difference in means: double diff = Sm - M; cout << setw(55) << left << "Sample Mean - Expected Test Mean" << "= " << diff << "\n"; // Degrees of freedom: unsigned v = Sn - 1; cout << setw(55) << left << "Degrees of Freedom" << "= " << v << "\n"; // t-statistic: double t_stat = diff * sqrt(double(Sn)) / Sd; cout << setw(55) << left << "T Statistic" << "= " << t_stat << "\n"; // // Finally define our distribution, and get the probability: // students_t dist(v); double q = cdf(complement(dist, fabs(t_stat))); cout << setw(55) << left << "Probability that difference is due to chance" << "= " << setprecision(3) << scientific << 2 * q << "\n\n"; // // Finally print out results of alternative hypothesis: // cout << setw(55) << left << "Results for Alternative Hypothesis and alpha" << "= " << setprecision(4) << fixed << alpha << "\n\n"; cout << "Alternative Hypothesis Conclusion\n"; cout << "Mean != " << setprecision(3) << fixed << M << " "; if(q < alpha / 2) cout << "NOT REJECTED\n"; else cout << "REJECTED\n"; cout << "Mean < " << setprecision(3) << fixed << M << " "; if(cdf(dist, t_stat) < alpha) cout << "NOT REJECTED\n"; else cout << "REJECTED\n"; cout << "Mean > " << setprecision(3) << fixed << M << " "; if(cdf(complement(dist, t_stat)) < alpha) cout << "NOT REJECTED\n"; else cout << "REJECTED\n"; cout << endl << endl; } // void single_sample_t_test( void single_sample_find_df(double M, double Sm, double Sd) { // // M = true mean. // Sm = Sample Mean. // Sd = Sample Standard Deviation. // using boost::math::students_t; // Print out general info: cout << "_____________________________________________________________\n" "Estimated sample sizes required for various confidence levels\n" "_____________________________________________________________\n\n"; cout << setprecision(5); cout << setw(40) << left << "True Mean" << "= " << M << "\n"; cout << setw(40) << left << "Sample Mean" << "= " << Sm << "\n"; cout << setw(40) << left << "Sample Standard Deviation" << "= " << Sd << "\n"; // // Define a table of significance intervals: // double alpha[] = { 0.5, 0.25, 0.1, 0.05, 0.01, 0.001, 0.0001, 0.00001 }; // // Print table header: // cout << "\n\n" "_______________________________________________________________\n" "Confidence Estimated Estimated\n" " Value (%) Sample Size Sample Size\n" " (one sided test) (two sided test)\n" "_______________________________________________________________\n"; // // Now print out the data for the table rows. // for(unsigned i = 1; i < sizeof(alpha)/sizeof(alpha[0]); ++i) { // Confidence value: cout << fixed << setprecision(3) << setw(10) << right << 100 * (1-alpha[i]); // calculate df for single sided test: double df = students_t::find_degrees_of_freedom( fabs(M - Sm), alpha[i], alpha[i], Sd); // convert to sample size, always one more than the degrees of freedom: double size = ceil(df) + 1; // Print size: cout << fixed << setprecision(0) << setw(16) << right << size; // calculate df for two sided test: df = students_t::find_degrees_of_freedom( fabs(M - Sm), alpha[i]/2, alpha[i], Sd); // convert to sample size: size = ceil(df) + 1; // Print size: cout << fixed << setprecision(0) << setw(16) << right << size << endl; } cout << endl; } // void single_sample_find_df int main() { // // Run tests for Heat Flow Meter data // see http://www.itl.nist.gov/div898/handbook/eda/section4/eda428.htm // The data was collected while calibrating a heat flow meter // against a known value. // confidence_limits_on_mean(9.261460, 0.2278881e-01, 195); single_sample_t_test(5, 9.261460, 0.2278881e-01, 195, 0.05); single_sample_find_df(5, 9.261460, 0.2278881e-01); // // Data for this example from: // P.K.Hou, O. W. Lau & M.C. Wong, Analyst (1983) vol. 108, p 64. // from Statistics for Analytical Chemistry, 3rd ed. (1994), pp 54-55 // J. C. Miller and J. N. Miller, Ellis Horwood ISBN 0 13 0309907 // // Determination of mercury by cold-vapour atomic absorption, // the following values were obtained fusing a trusted // Standard Reference Material containing 38.9% mercury, // which we assume is correct or 'true'. // confidence_limits_on_mean(37.8, 0.964365, 3); // 95% test: single_sample_t_test(38.9, 37.8, 0.964365, 3, 0.05); // 90% test: single_sample_t_test(38.9, 37.8, 0.964365, 3, 0.1); // parameter estimate: single_sample_find_df(38.9, 37.8, 0.964365); return 0; } // int main() /* Output: ------ Rebuild All started: Project: students_t_single_sample, Configuration: Release Win32 ------ students_t_single_sample.cpp Generating code Finished generating code students_t_single_sample.vcxproj -> J:\Cpp\MathToolkit\test\Math_test\Release\students_t_single_sample.exe __________________________________ 2-Sided Confidence Limits For Mean __________________________________ Number of Observations = 195 Mean = 9.26146 Standard Deviation = 0.02278881 _______________________________________________________________ Confidence T Interval Lower Upper Value (%) Value Width Limit Limit _______________________________________________________________ 50.000 0.676 1.103e-003 9.26036 9.26256 75.000 1.154 1.883e-003 9.25958 9.26334 90.000 1.653 2.697e-003 9.25876 9.26416 95.000 1.972 3.219e-003 9.25824 9.26468 99.000 2.601 4.245e-003 9.25721 9.26571 99.900 3.341 5.453e-003 9.25601 9.26691 99.990 3.973 6.484e-003 9.25498 9.26794 99.999 4.537 7.404e-003 9.25406 9.26886 __________________________________ Student t test for a single sample __________________________________ Number of Observations = 195 Sample Mean = 9.26146 Sample Standard Deviation = 0.02279 Expected True Mean = 5.00000 Sample Mean - Expected Test Mean = 4.26146 Degrees of Freedom = 194 T Statistic = 2611.28380 Probability that difference is due to chance = 0.000e+000 Results for Alternative Hypothesis and alpha = 0.0500 Alternative Hypothesis Conclusion Mean != 5.000 NOT REJECTED Mean < 5.000 REJECTED Mean > 5.000 NOT REJECTED _____________________________________________________________ Estimated sample sizes required for various confidence levels _____________________________________________________________ True Mean = 5.00000 Sample Mean = 9.26146 Sample Standard Deviation = 0.02279 _______________________________________________________________ Confidence Estimated Estimated Value (%) Sample Size Sample Size (one sided test) (two sided test) _______________________________________________________________ 75.000 2 2 90.000 2 2 95.000 2 2 99.000 2 2 99.900 3 3 99.990 3 3 99.999 4 4 __________________________________ 2-Sided Confidence Limits For Mean __________________________________ Number of Observations = 3 Mean = 37.8000000 Standard Deviation = 0.9643650 _______________________________________________________________ Confidence T Interval Lower Upper Value (%) Value Width Limit Limit _______________________________________________________________ 50.000 0.816 0.455 37.34539 38.25461 75.000 1.604 0.893 36.90717 38.69283 90.000 2.920 1.626 36.17422 39.42578 95.000 4.303 2.396 35.40438 40.19562 99.000 9.925 5.526 32.27408 43.32592 99.900 31.599 17.594 20.20639 55.39361 99.990 99.992 55.673 -17.87346 93.47346 99.999 316.225 176.067 -138.26683 213.86683 __________________________________ Student t test for a single sample __________________________________ Number of Observations = 3 Sample Mean = 37.80000 Sample Standard Deviation = 0.96437 Expected True Mean = 38.90000 Sample Mean - Expected Test Mean = -1.10000 Degrees of Freedom = 2 T Statistic = -1.97566 Probability that difference is due to chance = 1.869e-001 Results for Alternative Hypothesis and alpha = 0.0500 Alternative Hypothesis Conclusion Mean != 38.900 REJECTED Mean < 38.900 REJECTED Mean > 38.900 REJECTED __________________________________ Student t test for a single sample __________________________________ Number of Observations = 3 Sample Mean = 37.80000 Sample Standard Deviation = 0.96437 Expected True Mean = 38.90000 Sample Mean - Expected Test Mean = -1.10000 Degrees of Freedom = 2 T Statistic = -1.97566 Probability that difference is due to chance = 1.869e-001 Results for Alternative Hypothesis and alpha = 0.1000 Alternative Hypothesis Conclusion Mean != 38.900 REJECTED Mean < 38.900 NOT REJECTED Mean > 38.900 REJECTED _____________________________________________________________ Estimated sample sizes required for various confidence levels _____________________________________________________________ True Mean = 38.90000 Sample Mean = 37.80000 Sample Standard Deviation = 0.96437 _______________________________________________________________ Confidence Estimated Estimated Value (%) Sample Size Sample Size (one sided test) (two sided test) _______________________________________________________________ 75.000 3 4 90.000 7 9 95.000 11 13 99.000 20 22 99.900 35 37 99.990 50 53 99.999 66 68 */