1 //===----------------------------------------------------------------------===// 2 // 3 // The LLVM Compiler Infrastructure 4 // 5 // This file is dual licensed under the MIT and the University of Illinois Open 6 // Source Licenses. See LICENSE.TXT for details. 7 // 8 //===----------------------------------------------------------------------===// 9 // 10 // REQUIRES: long_tests 11 12 // <random> 13 14 // template<class RealType = double> 15 // class student_t_distribution 16 17 // template<class _URNG> result_type operator()(_URNG& g); 18 19 #include <random> 20 #include <cassert> 21 #include <vector> 22 #include <numeric> 23 24 template <class T> 25 inline 26 T 27 sqr(T x) 28 { 29 return x * x; 30 } 31 32 int main() 33 { 34 { 35 typedef std::student_t_distribution<> D; 36 typedef std::minstd_rand G; 37 G g; 38 D d(5.5); 39 const int N = 1000000; 40 std::vector<D::result_type> u; 41 for (int i = 0; i < N; ++i) 42 u.push_back(d(g)); 43 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 44 double var = 0; 45 double skew = 0; 46 double kurtosis = 0; 47 for (unsigned i = 0; i < u.size(); ++i) 48 { 49 double dbl = (u[i] - mean); 50 double d2 = sqr(dbl); 51 var += d2; 52 skew += dbl * d2; 53 kurtosis += d2 * d2; 54 } 55 var /= u.size(); 56 double dev = std::sqrt(var); 57 skew /= u.size() * dev * var; 58 kurtosis /= u.size() * var * var; 59 kurtosis -= 3; 60 double x_mean = 0; 61 double x_var = d.n() / (d.n() - 2); 62 double x_skew = 0; 63 double x_kurtosis = 6 / (d.n() - 4); 64 assert(std::abs(mean - x_mean) < 0.01); 65 assert(std::abs((var - x_var) / x_var) < 0.01); 66 assert(std::abs(skew - x_skew) < 0.01); 67 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.2); 68 } 69 { 70 typedef std::student_t_distribution<> D; 71 typedef std::minstd_rand G; 72 G g; 73 D d(10); 74 const int N = 1000000; 75 std::vector<D::result_type> u; 76 for (int i = 0; i < N; ++i) 77 u.push_back(d(g)); 78 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 79 double var = 0; 80 double skew = 0; 81 double kurtosis = 0; 82 for (unsigned i = 0; i < u.size(); ++i) 83 { 84 double dbl = (u[i] - mean); 85 double d2 = sqr(dbl); 86 var += d2; 87 skew += dbl * d2; 88 kurtosis += d2 * d2; 89 } 90 var /= u.size(); 91 double dev = std::sqrt(var); 92 skew /= u.size() * dev * var; 93 kurtosis /= u.size() * var * var; 94 kurtosis -= 3; 95 double x_mean = 0; 96 double x_var = d.n() / (d.n() - 2); 97 double x_skew = 0; 98 double x_kurtosis = 6 / (d.n() - 4); 99 assert(std::abs(mean - x_mean) < 0.01); 100 assert(std::abs((var - x_var) / x_var) < 0.01); 101 assert(std::abs(skew - x_skew) < 0.01); 102 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.04); 103 } 104 { 105 typedef std::student_t_distribution<> D; 106 typedef std::minstd_rand G; 107 G g; 108 D d(100); 109 const int N = 1000000; 110 std::vector<D::result_type> u; 111 for (int i = 0; i < N; ++i) 112 u.push_back(d(g)); 113 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 114 double var = 0; 115 double skew = 0; 116 double kurtosis = 0; 117 for (unsigned i = 0; i < u.size(); ++i) 118 { 119 double dbl = (u[i] - mean); 120 double d2 = sqr(dbl); 121 var += d2; 122 skew += dbl * d2; 123 kurtosis += d2 * d2; 124 } 125 var /= u.size(); 126 double dev = std::sqrt(var); 127 skew /= u.size() * dev * var; 128 kurtosis /= u.size() * var * var; 129 kurtosis -= 3; 130 double x_mean = 0; 131 double x_var = d.n() / (d.n() - 2); 132 double x_skew = 0; 133 double x_kurtosis = 6 / (d.n() - 4); 134 assert(std::abs(mean - x_mean) < 0.01); 135 assert(std::abs((var - x_var) / x_var) < 0.01); 136 assert(std::abs(skew - x_skew) < 0.01); 137 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02); 138 } 139 } 140
