1 // Copyright John Maddock 2007.
2 // Copyright Paul A. Bristow 2007.
4 // Use, modification and distribution are subject to the
5 // Boost Software License, Version 1.0. (See accompanying file
6 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
8 #ifndef BOOST_STATS_FIND_SCALE_HPP
9 #define BOOST_STATS_FIND_SCALE_HPP
11 #include <boost/math/distributions/fwd.hpp> // for all distribution signatures.
12 #include <boost/math/distributions/complement.hpp>
13 #include <boost/math/policies/policy.hpp>
14 // using boost::math::policies::policy;
15 #include <boost/math/tools/traits.hpp>
16 #include <boost/static_assert.hpp>
17 #include <boost/math/special_functions/fpclassify.hpp>
18 #include <boost/math/policies/error_handling.hpp>
19 // using boost::math::complement; // will be needed by users who want complement,
20 // but NOT placed here to avoid putting it in global scope.
26 // Function to find location of random variable z
27 // to give probability p (given scale)
28 // Applies to normal, lognormal, extreme value, Cauchy, (and symmetrical triangular),
29 // distributions that have scale.
30 // BOOST_STATIC_ASSERTs, see below, are used to enforce this.
32 template <class Dist, class Policy>
34 typename Dist::value_type find_scale( // For example, normal mean.
35 typename Dist::value_type z, // location of random variable z to give probability, P(X > z) == p.
36 // For example, a nominal minimum acceptable weight z, so that p * 100 % are > z
37 typename Dist::value_type p, // probability value desired at x, say 0.95 for 95% > z.
38 typename Dist::value_type location, // location parameter, for example, normal distribution mean.
42 #if !defined(BOOST_NO_SFINAE) && !BOOST_WORKAROUND(__SUNPRO_CC, BOOST_TESTED_AT(0x590))
43 BOOST_STATIC_ASSERT(::boost::math::tools::is_distribution<Dist>::value);
44 BOOST_STATIC_ASSERT(::boost::math::tools::is_scaled_distribution<Dist>::value);
46 static const char* function = "boost::math::find_scale<Dist, Policy>(%1%, %1%, %1%, Policy)";
48 if(!(boost::math::isfinite)(p) || (p < 0) || (p > 1))
50 return policies::raise_domain_error<typename Dist::value_type>(
51 function, "Probability parameter was %1%, but must be >= 0 and <= 1!", p, pol);
53 if(!(boost::math::isfinite)(z))
55 return policies::raise_domain_error<typename Dist::value_type>(
56 function, "find_scale z parameter was %1%, but must be finite!", z, pol);
58 if(!(boost::math::isfinite)(location))
60 return policies::raise_domain_error<typename Dist::value_type>(
61 function, "find_scale location parameter was %1%, but must be finite!", location, pol);
64 //cout << "z " << z << ", p " << p << ", quantile(Dist(), p) "
65 //<< quantile(Dist(), p) << ", z - mean " << z - location
66 //<<", sd " << (z - location) / quantile(Dist(), p) << endl;
68 //quantile(N01, 0.001) -3.09023
69 //quantile(N01, 0.01) -2.32635
70 //quantile(N01, 0.05) -1.64485
71 //quantile(N01, 0.333333) -0.430728
72 //quantile(N01, 0.5) 0
73 //quantile(N01, 0.666667) 0.430728
74 //quantile(N01, 0.9) 1.28155
75 //quantile(N01, 0.95) 1.64485
76 //quantile(N01, 0.99) 2.32635
77 //quantile(N01, 0.999) 3.09023
79 typename Dist::value_type result =
80 (z - location) // difference between desired x and current location.
81 / quantile(Dist(), p); // standard distribution.
84 { // If policy isn't to throw, return the scale <= 0.
85 policies::raise_evaluation_error<typename Dist::value_type>(function,
86 "Computed scale (%1%) is <= 0!" " Was the complement intended?",
90 } // template <class Dist, class Policy> find_scale
93 inline // with default policy.
94 typename Dist::value_type find_scale( // For example, normal mean.
95 typename Dist::value_type z, // location of random variable z to give probability, P(X > z) == p.
96 // For example, a nominal minimum acceptable z, so that p * 100 % are > z
97 typename Dist::value_type p, // probability value desired at x, say 0.95 for 95% > z.
98 typename Dist::value_type location) // location parameter, for example, mean.
99 { // Forward to find_scale using the default policy.
100 return (find_scale<Dist>(z, p, location, policies::policy<>()));
103 template <class Dist, class Real1, class Real2, class Real3, class Policy>
104 inline typename Dist::value_type find_scale(
105 complemented4_type<Real1, Real2, Real3, Policy> const& c)
107 //cout << "cparam1 q " << c.param1 // q
108 // << ", c.dist z " << c.dist // z
109 // << ", c.param2 l " << c.param2 // l
110 // << ", quantile (Dist(), c.param1 = q) "
111 // << quantile(Dist(), c.param1) //q
114 #if !defined(BOOST_NO_SFINAE) && !BOOST_WORKAROUND(__SUNPRO_CC, BOOST_TESTED_AT(0x590))
115 BOOST_STATIC_ASSERT(::boost::math::tools::is_distribution<Dist>::value);
116 BOOST_STATIC_ASSERT(::boost::math::tools::is_scaled_distribution<Dist>::value);
118 static const char* function = "boost::math::find_scale<Dist, Policy>(complement(%1%, %1%, %1%, Policy))";
120 // Checks on arguments, as not complemented version,
122 typename Dist::value_type q = c.param1;
123 if(!(boost::math::isfinite)(q) || (q < 0) || (q > 1))
125 return policies::raise_domain_error<typename Dist::value_type>(
126 function, "Probability parameter was %1%, but must be >= 0 and <= 1!", q, c.param3);
128 typename Dist::value_type z = c.dist;
129 if(!(boost::math::isfinite)(z))
131 return policies::raise_domain_error<typename Dist::value_type>(
132 function, "find_scale z parameter was %1%, but must be finite!", z, c.param3);
134 typename Dist::value_type location = c.param2;
135 if(!(boost::math::isfinite)(location))
137 return policies::raise_domain_error<typename Dist::value_type>(
138 function, "find_scale location parameter was %1%, but must be finite!", location, c.param3);
141 typename Dist::value_type result =
142 (c.dist - c.param2) // difference between desired x and current location.
143 / quantile(complement(Dist(), c.param1));
144 // ( z - location) / (quantile(complement(Dist(), q))
146 { // If policy isn't to throw, return the scale <= 0.
147 policies::raise_evaluation_error<typename Dist::value_type>(function,
148 "Computed scale (%1%) is <= 0!" " Was the complement intended?",
152 } // template <class Dist, class Policy, class Real1, class Real2, class Real3> typename Dist::value_type find_scale
154 // So the user can start from the complement q = (1 - p) of the probability p,
155 // for example, s = find_scale<normal>(complement(z, q, l));
157 template <class Dist, class Real1, class Real2, class Real3>
158 inline typename Dist::value_type find_scale(
159 complemented3_type<Real1, Real2, Real3> const& c)
161 //cout << "cparam1 q " << c.param1 // q
162 // << ", c.dist z " << c.dist // z
163 // << ", c.param2 l " << c.param2 // l
164 // << ", quantile (Dist(), c.param1 = q) "
165 // << quantile(Dist(), c.param1) //q
168 #if !defined(BOOST_NO_SFINAE) && !BOOST_WORKAROUND(__SUNPRO_CC, BOOST_TESTED_AT(0x590))
169 BOOST_STATIC_ASSERT(::boost::math::tools::is_distribution<Dist>::value);
170 BOOST_STATIC_ASSERT(::boost::math::tools::is_scaled_distribution<Dist>::value);
172 static const char* function = "boost::math::find_scale<Dist, Policy>(complement(%1%, %1%, %1%, Policy))";
174 // Checks on arguments, as not complemented version,
175 // default policy policies::policy<>().
176 typename Dist::value_type q = c.param1;
177 if(!(boost::math::isfinite)(q) || (q < 0) || (q > 1))
179 return policies::raise_domain_error<typename Dist::value_type>(
180 function, "Probability parameter was %1%, but must be >= 0 and <= 1!", q, policies::policy<>());
182 typename Dist::value_type z = c.dist;
183 if(!(boost::math::isfinite)(z))
185 return policies::raise_domain_error<typename Dist::value_type>(
186 function, "find_scale z parameter was %1%, but must be finite!", z, policies::policy<>());
188 typename Dist::value_type location = c.param2;
189 if(!(boost::math::isfinite)(location))
191 return policies::raise_domain_error<typename Dist::value_type>(
192 function, "find_scale location parameter was %1%, but must be finite!", location, policies::policy<>());
195 typename Dist::value_type result =
196 (z - location) // difference between desired x and current location.
197 / quantile(complement(Dist(), q));
198 // ( z - location) / (quantile(complement(Dist(), q))
200 { // If policy isn't to throw, return the scale <= 0.
201 policies::raise_evaluation_error<typename Dist::value_type>(function,
202 "Computed scale (%1%) is <= 0!" " Was the complement intended?",
203 result, policies::policy<>()); // This is only the default policy - also Want a version with Policy here.
206 } // template <class Dist, class Real1, class Real2, class Real3> typename Dist::value_type find_scale
211 #endif // BOOST_STATS_FIND_SCALE_HPP