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