1 // Copyright John Maddock 2006.
2 // Copyright Paul A. Bristow 2007
3
4 // Use, modification and distribution are subject to the
5 // Boost Software License, Version 1.0.
6 // (See accompanying file LICENSE_1_0.txt
7 // or copy at http://www.boost.org/LICENSE_1_0.txt)
8
9 // test_lognormal.cpp
10
11 #include <boost/math/concepts/real_concept.hpp> // for real_concept
12 #define BOOST_TEST_MAIN
13 #include <boost/test/unit_test.hpp> // Boost.Test
14 #include <boost/test/floating_point_comparison.hpp>
15
16 #include <boost/math/distributions/lognormal.hpp>
17 using boost::math::lognormal_distribution;
18 #include <boost/math/tools/test.hpp>
19 #include "test_out_of_range.hpp"
20
21 #include <iostream>
22 using std::cout;
23 using std::endl;
24 using std::setprecision;
25 #include <limits>
26 using std::numeric_limits;
27 #include <cassert>
28
29 template <class RealType>
check_lognormal(RealType loc,RealType scale,RealType x,RealType p,RealType q,RealType tol)30 void check_lognormal(RealType loc, RealType scale, RealType x, RealType p, RealType q, RealType tol)
31 {
32 BOOST_CHECK_CLOSE(
33 ::boost::math::cdf(
34 lognormal_distribution<RealType>(loc, scale), // distribution.
35 x), // random variable.
36 p, // probability.
37 tol); // %tolerance.
38 BOOST_CHECK_CLOSE(
39 ::boost::math::cdf(
40 complement(
41 lognormal_distribution<RealType>(loc, scale), // distribution.
42 x)), // random variable.
43 q, // probability complement.
44 tol); // %tolerance.
45 BOOST_CHECK_CLOSE(
46 ::boost::math::quantile(
47 lognormal_distribution<RealType>(loc, scale), // distribution.
48 p), // probability.
49 x, // random variable.
50 tol); // %tolerance.
51 BOOST_CHECK_CLOSE(
52 ::boost::math::quantile(
53 complement(
54 lognormal_distribution<RealType>(loc, scale), // distribution.
55 q)), // probability complement.
56 x, // random variable.
57 tol); // %tolerance.
58 }
59
60 template <class RealType>
test_spots(RealType)61 void test_spots(RealType)
62 {
63
64 // Basic sanity checks.
65 RealType tolerance = 5e-3 * 100;
66 // Some tests only pass at 1e-4 because values generated by
67 // http://faculty.vassar.edu/lowry/VassarStats.html
68 // give only 5 or 6 *fixed* places, so small values have fewer digits.
69
70 cout << "Tolerance for type " << typeid(RealType).name() << " is " << tolerance << " %" << endl;
71
72 using std::exp;
73
74 //
75 // These test values were generated for the normal distribution
76 // using the online calculator at http://faculty.vassar.edu/lowry/VassarStats.html
77 // and then exponentiating the random variate x.
78 //
79
80 check_lognormal(
81 static_cast<RealType>(0), // location
82 static_cast<RealType>(5), // scale
83 static_cast<RealType>(1), // x
84 static_cast<RealType>(0.5), // p
85 static_cast<RealType>(0.5), // q
86 tolerance);
87 check_lognormal(
88 static_cast<RealType>(2), // location
89 static_cast<RealType>(2), // scale
90 static_cast<RealType>(exp(1.8)), // x
91 static_cast<RealType>(0.46017), // p
92 static_cast<RealType>(1-0.46017), // q
93 tolerance);
94 check_lognormal(
95 static_cast<RealType>(2), // location
96 static_cast<RealType>(2), // scale
97 static_cast<RealType>(exp(2.2)), // x
98 static_cast<RealType>(1-0.46017), // p
99 static_cast<RealType>(0.46017), // q
100 tolerance);
101 check_lognormal(
102 static_cast<RealType>(2), // location
103 static_cast<RealType>(2), // scale
104 static_cast<RealType>(exp(-1.4)), // x
105 static_cast<RealType>(0.04457), // p
106 static_cast<RealType>(1-0.04457), // q
107 tolerance);
108 check_lognormal(
109 static_cast<RealType>(2), // location
110 static_cast<RealType>(2), // scale
111 static_cast<RealType>(exp(5.4)), // x
112 static_cast<RealType>(1-0.04457), // p
113 static_cast<RealType>(0.04457), // q
114 tolerance);
115
116 check_lognormal(
117 static_cast<RealType>(-3), // location
118 static_cast<RealType>(5), // scale
119 static_cast<RealType>(exp(-5.0)), // x
120 static_cast<RealType>(0.34458), // p
121 static_cast<RealType>(1-0.34458), // q
122 tolerance);
123 check_lognormal(
124 static_cast<RealType>(-3), // location
125 static_cast<RealType>(5), // scale
126 static_cast<RealType>(exp(-1.0)), // x
127 static_cast<RealType>(1-0.34458), // p
128 static_cast<RealType>(0.34458), // q
129 tolerance);
130 check_lognormal(
131 static_cast<RealType>(-3), // location
132 static_cast<RealType>(5), // scale
133 static_cast<RealType>(exp(-9.0)), // x
134 static_cast<RealType>(0.11507), // p
135 static_cast<RealType>(1-0.11507), // q
136 tolerance);
137 check_lognormal(
138 static_cast<RealType>(-3), // location
139 static_cast<RealType>(5), // scale
140 static_cast<RealType>(exp(3.0)), // x
141 static_cast<RealType>(1-0.11507), // p
142 static_cast<RealType>(0.11507), // q
143 tolerance);
144
145 //
146 // Tests for PDF
147 //
148 tolerance = boost::math::tools::epsilon<RealType>() * 5 * 100; // 5 eps as a percentage
149 cout << "Tolerance for type " << typeid(RealType).name() << " is " << tolerance << " %" << endl;
150 BOOST_CHECK_CLOSE(
151 pdf(lognormal_distribution<RealType>(), static_cast<RealType>(1)),
152 static_cast<RealType>(0.3989422804014326779399460599343818684759L), // 1/sqrt(2*pi)
153 tolerance);
154 BOOST_CHECK_CLOSE(
155 pdf(lognormal_distribution<RealType>(3), exp(static_cast<RealType>(3))),
156 static_cast<RealType>(0.3989422804014326779399460599343818684759L) / exp(static_cast<RealType>(3)),
157 tolerance);
158 BOOST_CHECK_CLOSE(
159 pdf(lognormal_distribution<RealType>(3, 5), exp(static_cast<RealType>(3))),
160 static_cast<RealType>(0.3989422804014326779399460599343818684759L / (5 * exp(static_cast<RealType>(3)))),
161 tolerance);
162 //
163 // Spot checks for location = -5, scale = 6,
164 // use relation to normal to test:
165 //
166 for(RealType x = -15; x < 5; x += 0.125)
167 {
168 BOOST_CHECK_CLOSE(
169 pdf(lognormal_distribution<RealType>(-5, 6), exp(x)),
170 pdf(boost::math::normal_distribution<RealType>(-5, 6), x) / exp(x),
171 tolerance);
172 }
173
174 //
175 // These test values were obtained by punching numbers into
176 // a calculator, using the formulas at http://mathworld.wolfram.com/LogNormalDistribution.html
177 //
178 tolerance = (std::max)(
179 boost::math::tools::epsilon<RealType>(),
180 static_cast<RealType>(boost::math::tools::epsilon<double>())) * 5 * 100; // 5 eps as a percentage
181 cout << "Tolerance for type " << typeid(RealType).name() << " is " << tolerance << " %" << endl;
182 lognormal_distribution<RealType> dist(8, 3);
183 RealType x = static_cast<RealType>(0.125);
184
185 BOOST_MATH_STD_USING // ADL of std math lib names.
186
187 // mean:
188 BOOST_CHECK_CLOSE(
189 mean(dist)
190 , static_cast<RealType>(268337.28652087445695647967378715L), tolerance);
191 // variance:
192 BOOST_CHECK_CLOSE(
193 variance(dist)
194 , static_cast<RealType>(583389737628117.49553037857325892L), tolerance);
195 // std deviation:
196 BOOST_CHECK_CLOSE(
197 standard_deviation(dist)
198 , static_cast<RealType>(24153462.228594009489719473727471L), tolerance);
199 // hazard:
200 BOOST_CHECK_CLOSE(
201 hazard(dist, x)
202 , pdf(dist, x) / cdf(complement(dist, x)), tolerance);
203 // cumulative hazard:
204 BOOST_CHECK_CLOSE(
205 chf(dist, x)
206 , -log(cdf(complement(dist, x))), tolerance);
207 // coefficient_of_variation:
208 BOOST_CHECK_CLOSE(
209 coefficient_of_variation(dist)
210 , standard_deviation(dist) / mean(dist), tolerance);
211 // mode:
212 BOOST_CHECK_CLOSE(
213 mode(dist)
214 , static_cast<RealType>(0.36787944117144232159552377016146L), tolerance);
215
216 BOOST_CHECK_CLOSE(
217 median(dist)
218 , static_cast<RealType>(exp(dist.location())), tolerance);
219
220 BOOST_CHECK_CLOSE(
221 median(dist),
222 quantile(dist, static_cast<RealType>(0.5)), tolerance);
223
224 // skewness:
225 BOOST_CHECK_CLOSE(
226 skewness(dist)
227 , static_cast<RealType>(729551.38304660255658441529235697L), tolerance);
228 // kertosis:
229 BOOST_CHECK_CLOSE(
230 kurtosis(dist)
231 , static_cast<RealType>(4312295840576303.2363383232038251L), tolerance);
232 // kertosis excess:
233 BOOST_CHECK_CLOSE(
234 kurtosis_excess(dist)
235 , static_cast<RealType>(4312295840576300.2363383232038251L), tolerance);
236
237 BOOST_CHECK_CLOSE(
238 range(dist).first
239 , static_cast<RealType>(0), tolerance);
240
241 //
242 // Special cases:
243 //
244 BOOST_CHECK(pdf(dist, 0) == 0);
245 BOOST_CHECK(cdf(dist, 0) == 0);
246 BOOST_CHECK(cdf(complement(dist, 0)) == 1);
247 BOOST_CHECK(quantile(dist, 0) == 0);
248 BOOST_CHECK(quantile(complement(dist, 1)) == 0);
249
250 //
251 // Error checks:
252 //
253 BOOST_CHECK_THROW(lognormal_distribution<RealType>(0, 0), std::domain_error);
254 BOOST_CHECK_THROW(lognormal_distribution<RealType>(2, -1), std::domain_error);
255 BOOST_CHECK_THROW(pdf(dist, -1), std::domain_error);
256 BOOST_CHECK_THROW(cdf(dist, -1), std::domain_error);
257 BOOST_CHECK_THROW(cdf(complement(dist, -1)), std::domain_error);
258 BOOST_CHECK_THROW(quantile(dist, 1), std::overflow_error);
259 BOOST_CHECK_THROW(quantile(complement(dist, 0)), std::overflow_error);
260 check_out_of_range<lognormal_distribution<RealType> >(1, 2);
261
262 } // template <class RealType>void test_spots(RealType)
263
BOOST_AUTO_TEST_CASE(test_main)264 BOOST_AUTO_TEST_CASE( test_main )
265 {
266
267 // Check that can generate lognormal distribution using the two convenience methods:
268 boost::math::lognormal myf1(1., 2); // Using typedef
269 lognormal_distribution<> myf2(1., 2); // Using default RealType double.
270
271 // Test range and support using double only,
272 // because it supports numeric_limits max for a pseudo-infinity.
273 BOOST_CHECK_EQUAL(range(myf2).first, 0); // range 0 to +infinity
274 BOOST_CHECK_EQUAL(range(myf2).second, (std::numeric_limits<double>::max)());
275 BOOST_CHECK_EQUAL(support(myf2).first, 0); // support 0 to + infinity.
276 BOOST_CHECK_EQUAL(support(myf2).second, (std::numeric_limits<double>::max)());
277
278 // Basic sanity-check spot values.
279 // (Parameter value, arbitrarily zero, only communicates the floating point type).
280 test_spots(0.0F); // Test float. OK at decdigits = 0 tolerance = 0.0001 %
281 test_spots(0.0); // Test double. OK at decdigits 7, tolerance = 1e07 %
282 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
283 test_spots(0.0L); // Test long double.
284 #if !BOOST_WORKAROUND(__BORLANDC__, BOOST_TESTED_AT(0x0582))
285 test_spots(boost::math::concepts::real_concept(0.)); // Test real concept.
286 #endif
287 #else
288 std::cout << "<note>The long double tests have been disabled on this platform "
289 "either because the long double overloads of the usual math functions are "
290 "not available at all, or because they are too inaccurate for these tests "
291 "to pass.</note>" << std::cout;
292 #endif
293
294
295 } // BOOST_AUTO_TEST_CASE( test_main )
296
297 /*
298 Running 1 test case...
299 Tolerance for type float is 0.5 %
300 Tolerance for type float is 5.96046e-005 %
301 Tolerance for type float is 5.96046e-005 %
302 Tolerance for type double is 0.5 %
303 Tolerance for type double is 1.11022e-013 %
304 Tolerance for type double is 1.11022e-013 %
305 Tolerance for type long double is 0.5 %
306 Tolerance for type long double is 1.11022e-013 %
307 Tolerance for type long double is 1.11022e-013 %
308 Tolerance for type class boost::math::concepts::real_concept is 0.5 %
309 Tolerance for type class boost::math::concepts::real_concept is 1.11022e-013 %
310 Tolerance for type class boost::math::concepts::real_concept is 1.11022e-013 %
311 *** No errors detected
312 */
313
314
315