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 // <random> 11 12 // template<class RealType = double> 13 // class lognormal_distribution 14 15 // template<class _URNG> result_type operator()(_URNG& g); 16 17 #include <random> 18 #include <cassert> 19 #include <vector> 20 #include <numeric> 21 22 template <class T> 23 inline 24 T 25 sqr(T x) 26 { 27 return x * x; 28 } 29 30 int main() 31 { 32 { 33 typedef std::lognormal_distribution<> D; 34 typedef D::param_type P; 35 typedef std::mt19937 G; 36 G g; 37 D d(-1./8192, 0.015625); 38 const int N = 1000000; 39 std::vector<D::result_type> u; 40 for (int i = 0; i < N; ++i) 41 { 42 D::result_type v = d(g); 43 assert(v > 0); 44 u.push_back(v); 45 } 46 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 47 double var = 0; 48 double skew = 0; 49 double kurtosis = 0; 50 for (int i = 0; i < u.size(); ++i) 51 { 52 double d = (u[i] - mean); 53 double d2 = sqr(d); 54 var += d2; 55 skew += d * d2; 56 kurtosis += d2 * d2; 57 } 58 var /= u.size(); 59 double dev = std::sqrt(var); 60 skew /= u.size() * dev * var; 61 kurtosis /= u.size() * var * var; 62 kurtosis -= 3; 63 double x_mean = std::exp(d.m() + sqr(d.s())/2); 64 double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 65 double x_skew = (std::exp(sqr(d.s())) + 2) * 66 std::sqrt((std::exp(sqr(d.s())) - 1)); 67 double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 68 3*std::exp(2*sqr(d.s())) - 6; 69 assert(std::abs((mean - x_mean) / x_mean) < 0.01); 70 assert(std::abs((var - x_var) / x_var) < 0.01); 71 assert(std::abs((skew - x_skew) / x_skew) < 0.05); 72 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.25); 73 } 74 { 75 typedef std::lognormal_distribution<> D; 76 typedef D::param_type P; 77 typedef std::mt19937 G; 78 G g; 79 D d(-1./32, 0.25); 80 const int N = 1000000; 81 std::vector<D::result_type> u; 82 for (int i = 0; i < N; ++i) 83 { 84 D::result_type v = d(g); 85 assert(v > 0); 86 u.push_back(v); 87 } 88 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 89 double var = 0; 90 double skew = 0; 91 double kurtosis = 0; 92 for (int i = 0; i < u.size(); ++i) 93 { 94 double d = (u[i] - mean); 95 double d2 = sqr(d); 96 var += d2; 97 skew += d * d2; 98 kurtosis += d2 * d2; 99 } 100 var /= u.size(); 101 double dev = std::sqrt(var); 102 skew /= u.size() * dev * var; 103 kurtosis /= u.size() * var * var; 104 kurtosis -= 3; 105 double x_mean = std::exp(d.m() + sqr(d.s())/2); 106 double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 107 double x_skew = (std::exp(sqr(d.s())) + 2) * 108 std::sqrt((std::exp(sqr(d.s())) - 1)); 109 double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 110 3*std::exp(2*sqr(d.s())) - 6; 111 assert(std::abs((mean - x_mean) / x_mean) < 0.01); 112 assert(std::abs((var - x_var) / x_var) < 0.01); 113 assert(std::abs((skew - x_skew) / x_skew) < 0.01); 114 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03); 115 } 116 { 117 typedef std::lognormal_distribution<> D; 118 typedef D::param_type P; 119 typedef std::mt19937 G; 120 G g; 121 D d(-1./8, 0.5); 122 const int N = 1000000; 123 std::vector<D::result_type> u; 124 for (int i = 0; i < N; ++i) 125 { 126 D::result_type v = d(g); 127 assert(v > 0); 128 u.push_back(v); 129 } 130 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 131 double var = 0; 132 double skew = 0; 133 double kurtosis = 0; 134 for (int i = 0; i < u.size(); ++i) 135 { 136 double d = (u[i] - mean); 137 double d2 = sqr(d); 138 var += d2; 139 skew += d * d2; 140 kurtosis += d2 * d2; 141 } 142 var /= u.size(); 143 double dev = std::sqrt(var); 144 skew /= u.size() * dev * var; 145 kurtosis /= u.size() * var * var; 146 kurtosis -= 3; 147 double x_mean = std::exp(d.m() + sqr(d.s())/2); 148 double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 149 double x_skew = (std::exp(sqr(d.s())) + 2) * 150 std::sqrt((std::exp(sqr(d.s())) - 1)); 151 double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 152 3*std::exp(2*sqr(d.s())) - 6; 153 assert(std::abs((mean - x_mean) / x_mean) < 0.01); 154 assert(std::abs((var - x_var) / x_var) < 0.01); 155 assert(std::abs((skew - x_skew) / x_skew) < 0.02); 156 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.05); 157 } 158 { 159 typedef std::lognormal_distribution<> D; 160 typedef D::param_type P; 161 typedef std::mt19937 G; 162 G g; 163 D d; 164 const int N = 1000000; 165 std::vector<D::result_type> u; 166 for (int i = 0; i < N; ++i) 167 { 168 D::result_type v = d(g); 169 assert(v > 0); 170 u.push_back(v); 171 } 172 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 173 double var = 0; 174 double skew = 0; 175 double kurtosis = 0; 176 for (int i = 0; i < u.size(); ++i) 177 { 178 double d = (u[i] - mean); 179 double d2 = sqr(d); 180 var += d2; 181 skew += d * d2; 182 kurtosis += d2 * d2; 183 } 184 var /= u.size(); 185 double dev = std::sqrt(var); 186 skew /= u.size() * dev * var; 187 kurtosis /= u.size() * var * var; 188 kurtosis -= 3; 189 double x_mean = std::exp(d.m() + sqr(d.s())/2); 190 double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 191 double x_skew = (std::exp(sqr(d.s())) + 2) * 192 std::sqrt((std::exp(sqr(d.s())) - 1)); 193 double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 194 3*std::exp(2*sqr(d.s())) - 6; 195 assert(std::abs((mean - x_mean) / x_mean) < 0.01); 196 assert(std::abs((var - x_var) / x_var) < 0.02); 197 assert(std::abs((skew - x_skew) / x_skew) < 0.08); 198 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.4); 199 } 200 { 201 typedef std::lognormal_distribution<> D; 202 typedef D::param_type P; 203 typedef std::mt19937 G; 204 G g; 205 D d(-0.78125, 1.25); 206 const int N = 1000000; 207 std::vector<D::result_type> u; 208 for (int i = 0; i < N; ++i) 209 { 210 D::result_type v = d(g); 211 assert(v > 0); 212 u.push_back(v); 213 } 214 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 215 double var = 0; 216 double skew = 0; 217 double kurtosis = 0; 218 for (int i = 0; i < u.size(); ++i) 219 { 220 double d = (u[i] - mean); 221 double d2 = sqr(d); 222 var += d2; 223 skew += d * d2; 224 kurtosis += d2 * d2; 225 } 226 var /= u.size(); 227 double dev = std::sqrt(var); 228 skew /= u.size() * dev * var; 229 kurtosis /= u.size() * var * var; 230 kurtosis -= 3; 231 double x_mean = std::exp(d.m() + sqr(d.s())/2); 232 double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 233 double x_skew = (std::exp(sqr(d.s())) + 2) * 234 std::sqrt((std::exp(sqr(d.s())) - 1)); 235 double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 236 3*std::exp(2*sqr(d.s())) - 6; 237 assert(std::abs((mean - x_mean) / x_mean) < 0.01); 238 assert(std::abs((var - x_var) / x_var) < 0.04); 239 assert(std::abs((skew - x_skew) / x_skew) < 0.2); 240 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.7); 241 } 242 } 243