]>
Commit | Line | Data |
---|---|---|
7c673cae FG |
1 | // Copyright John Maddock 2006. |
2 | // Use, modification and distribution are subject to the | |
3 | // Boost Software License, Version 1.0. (See accompanying file | |
4 | // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) | |
5 | ||
6 | #ifndef BOOST_STATS_GAMMA_HPP | |
7 | #define BOOST_STATS_GAMMA_HPP | |
8 | ||
9 | // http://www.itl.nist.gov/div898/handbook/eda/section3/eda366b.htm | |
10 | // http://mathworld.wolfram.com/GammaDistribution.html | |
11 | // http://en.wikipedia.org/wiki/Gamma_distribution | |
12 | ||
13 | #include <boost/math/distributions/fwd.hpp> | |
14 | #include <boost/math/special_functions/gamma.hpp> | |
15 | #include <boost/math/distributions/detail/common_error_handling.hpp> | |
16 | #include <boost/math/distributions/complement.hpp> | |
17 | ||
18 | #include <utility> | |
19 | ||
20 | namespace boost{ namespace math | |
21 | { | |
22 | namespace detail | |
23 | { | |
24 | ||
25 | template <class RealType, class Policy> | |
26 | inline bool check_gamma_shape( | |
27 | const char* function, | |
28 | RealType shape, | |
29 | RealType* result, const Policy& pol) | |
30 | { | |
31 | if((shape <= 0) || !(boost::math::isfinite)(shape)) | |
32 | { | |
33 | *result = policies::raise_domain_error<RealType>( | |
34 | function, | |
35 | "Shape parameter is %1%, but must be > 0 !", shape, pol); | |
36 | return false; | |
37 | } | |
38 | return true; | |
39 | } | |
40 | ||
41 | template <class RealType, class Policy> | |
42 | inline bool check_gamma_x( | |
43 | const char* function, | |
44 | RealType const& x, | |
45 | RealType* result, const Policy& pol) | |
46 | { | |
47 | if((x < 0) || !(boost::math::isfinite)(x)) | |
48 | { | |
49 | *result = policies::raise_domain_error<RealType>( | |
50 | function, | |
51 | "Random variate is %1% but must be >= 0 !", x, pol); | |
52 | return false; | |
53 | } | |
54 | return true; | |
55 | } | |
56 | ||
57 | template <class RealType, class Policy> | |
58 | inline bool check_gamma( | |
59 | const char* function, | |
60 | RealType scale, | |
61 | RealType shape, | |
62 | RealType* result, const Policy& pol) | |
63 | { | |
64 | return check_scale(function, scale, result, pol) && check_gamma_shape(function, shape, result, pol); | |
65 | } | |
66 | ||
67 | } // namespace detail | |
68 | ||
69 | template <class RealType = double, class Policy = policies::policy<> > | |
70 | class gamma_distribution | |
71 | { | |
72 | public: | |
73 | typedef RealType value_type; | |
74 | typedef Policy policy_type; | |
75 | ||
76 | gamma_distribution(RealType l_shape, RealType l_scale = 1) | |
77 | : m_shape(l_shape), m_scale(l_scale) | |
78 | { | |
79 | RealType result; | |
80 | detail::check_gamma("boost::math::gamma_distribution<%1%>::gamma_distribution", l_scale, l_shape, &result, Policy()); | |
81 | } | |
82 | ||
83 | RealType shape()const | |
84 | { | |
85 | return m_shape; | |
86 | } | |
87 | ||
88 | RealType scale()const | |
89 | { | |
90 | return m_scale; | |
91 | } | |
92 | private: | |
93 | // | |
94 | // Data members: | |
95 | // | |
96 | RealType m_shape; // distribution shape | |
97 | RealType m_scale; // distribution scale | |
98 | }; | |
99 | ||
100 | // NO typedef because of clash with name of gamma function. | |
101 | ||
102 | template <class RealType, class Policy> | |
103 | inline const std::pair<RealType, RealType> range(const gamma_distribution<RealType, Policy>& /* dist */) | |
104 | { // Range of permissible values for random variable x. | |
105 | using boost::math::tools::max_value; | |
106 | return std::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>()); | |
107 | } | |
108 | ||
109 | template <class RealType, class Policy> | |
110 | inline const std::pair<RealType, RealType> support(const gamma_distribution<RealType, Policy>& /* dist */) | |
111 | { // Range of supported values for random variable x. | |
112 | // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero. | |
113 | using boost::math::tools::max_value; | |
114 | using boost::math::tools::min_value; | |
115 | return std::pair<RealType, RealType>(min_value<RealType>(), max_value<RealType>()); | |
116 | } | |
117 | ||
118 | template <class RealType, class Policy> | |
119 | inline RealType pdf(const gamma_distribution<RealType, Policy>& dist, const RealType& x) | |
120 | { | |
121 | BOOST_MATH_STD_USING // for ADL of std functions | |
122 | ||
123 | static const char* function = "boost::math::pdf(const gamma_distribution<%1%>&, %1%)"; | |
124 | ||
125 | RealType shape = dist.shape(); | |
126 | RealType scale = dist.scale(); | |
127 | ||
128 | RealType result = 0; | |
129 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
130 | return result; | |
131 | if(false == detail::check_gamma_x(function, x, &result, Policy())) | |
132 | return result; | |
133 | ||
134 | if(x == 0) | |
135 | { | |
136 | return 0; | |
137 | } | |
138 | result = gamma_p_derivative(shape, x / scale, Policy()) / scale; | |
139 | return result; | |
140 | ||
141 | ||
142 | template <class RealType, class Policy> | |
143 | inline RealType cdf(const gamma_distribution<RealType, Policy>& dist, const RealType& x) | |
144 | { | |
145 | BOOST_MATH_STD_USING // for ADL of std functions | |
146 | ||
147 | static const char* function = "boost::math::cdf(const gamma_distribution<%1%>&, %1%)"; | |
148 | ||
149 | RealType shape = dist.shape(); | |
150 | RealType scale = dist.scale(); | |
151 | ||
152 | RealType result = 0; | |
153 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
154 | return result; | |
155 | if(false == detail::check_gamma_x(function, x, &result, Policy())) | |
156 | return result; | |
157 | ||
158 | result = boost::math::gamma_p(shape, x / scale, Policy()); | |
159 | return result; | |
160 | } // cdf | |
161 | ||
162 | template <class RealType, class Policy> | |
163 | inline RealType quantile(const gamma_distribution<RealType, Policy>& dist, const RealType& p) | |
164 | { | |
165 | BOOST_MATH_STD_USING // for ADL of std functions | |
166 | ||
167 | static const char* function = "boost::math::quantile(const gamma_distribution<%1%>&, %1%)"; | |
168 | ||
169 | RealType shape = dist.shape(); | |
170 | RealType scale = dist.scale(); | |
171 | ||
172 | RealType result = 0; | |
173 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
174 | return result; | |
175 | if(false == detail::check_probability(function, p, &result, Policy())) | |
176 | return result; | |
177 | ||
178 | if(p == 1) | |
179 | return policies::raise_overflow_error<RealType>(function, 0, Policy()); | |
180 | ||
181 | result = gamma_p_inv(shape, p, Policy()) * scale; | |
182 | ||
183 | return result; | |
184 | } | |
185 | ||
186 | template <class RealType, class Policy> | |
187 | inline RealType cdf(const complemented2_type<gamma_distribution<RealType, Policy>, RealType>& c) | |
188 | { | |
189 | BOOST_MATH_STD_USING // for ADL of std functions | |
190 | ||
191 | static const char* function = "boost::math::quantile(const gamma_distribution<%1%>&, %1%)"; | |
192 | ||
193 | RealType shape = c.dist.shape(); | |
194 | RealType scale = c.dist.scale(); | |
195 | ||
196 | RealType result = 0; | |
197 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
198 | return result; | |
199 | if(false == detail::check_gamma_x(function, c.param, &result, Policy())) | |
200 | return result; | |
201 | ||
202 | result = gamma_q(shape, c.param / scale, Policy()); | |
203 | ||
204 | return result; | |
205 | } | |
206 | ||
207 | template <class RealType, class Policy> | |
208 | inline RealType quantile(const complemented2_type<gamma_distribution<RealType, Policy>, RealType>& c) | |
209 | { | |
210 | BOOST_MATH_STD_USING // for ADL of std functions | |
211 | ||
212 | static const char* function = "boost::math::quantile(const gamma_distribution<%1%>&, %1%)"; | |
213 | ||
214 | RealType shape = c.dist.shape(); | |
215 | RealType scale = c.dist.scale(); | |
216 | RealType q = c.param; | |
217 | ||
218 | RealType result = 0; | |
219 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
220 | return result; | |
221 | if(false == detail::check_probability(function, q, &result, Policy())) | |
222 | return result; | |
223 | ||
224 | if(q == 0) | |
225 | return policies::raise_overflow_error<RealType>(function, 0, Policy()); | |
226 | ||
227 | result = gamma_q_inv(shape, q, Policy()) * scale; | |
228 | ||
229 | return result; | |
230 | } | |
231 | ||
232 | template <class RealType, class Policy> | |
233 | inline RealType mean(const gamma_distribution<RealType, Policy>& dist) | |
234 | { | |
235 | BOOST_MATH_STD_USING // for ADL of std functions | |
236 | ||
237 | static const char* function = "boost::math::mean(const gamma_distribution<%1%>&)"; | |
238 | ||
239 | RealType shape = dist.shape(); | |
240 | RealType scale = dist.scale(); | |
241 | ||
242 | RealType result = 0; | |
243 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
244 | return result; | |
245 | ||
246 | result = shape * scale; | |
247 | return result; | |
248 | } | |
249 | ||
250 | template <class RealType, class Policy> | |
251 | inline RealType variance(const gamma_distribution<RealType, Policy>& dist) | |
252 | { | |
253 | BOOST_MATH_STD_USING // for ADL of std functions | |
254 | ||
255 | static const char* function = "boost::math::variance(const gamma_distribution<%1%>&)"; | |
256 | ||
257 | RealType shape = dist.shape(); | |
258 | RealType scale = dist.scale(); | |
259 | ||
260 | RealType result = 0; | |
261 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
262 | return result; | |
263 | ||
264 | result = shape * scale * scale; | |
265 | return result; | |
266 | } | |
267 | ||
268 | template <class RealType, class Policy> | |
269 | inline RealType mode(const gamma_distribution<RealType, Policy>& dist) | |
270 | { | |
271 | BOOST_MATH_STD_USING // for ADL of std functions | |
272 | ||
273 | static const char* function = "boost::math::mode(const gamma_distribution<%1%>&)"; | |
274 | ||
275 | RealType shape = dist.shape(); | |
276 | RealType scale = dist.scale(); | |
277 | ||
278 | RealType result = 0; | |
279 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
280 | return result; | |
281 | ||
282 | if(shape < 1) | |
283 | return policies::raise_domain_error<RealType>( | |
284 | function, | |
285 | "The mode of the gamma distribution is only defined for values of the shape parameter >= 1, but got %1%.", | |
286 | shape, Policy()); | |
287 | ||
288 | result = (shape - 1) * scale; | |
289 | return result; | |
290 | } | |
291 | ||
292 | //template <class RealType, class Policy> | |
293 | //inline RealType median(const gamma_distribution<RealType, Policy>& dist) | |
294 | //{ // Rely on default definition in derived accessors. | |
295 | //} | |
296 | ||
297 | template <class RealType, class Policy> | |
298 | inline RealType skewness(const gamma_distribution<RealType, Policy>& dist) | |
299 | { | |
300 | BOOST_MATH_STD_USING // for ADL of std functions | |
301 | ||
302 | static const char* function = "boost::math::skewness(const gamma_distribution<%1%>&)"; | |
303 | ||
304 | RealType shape = dist.shape(); | |
305 | RealType scale = dist.scale(); | |
306 | ||
307 | RealType result = 0; | |
308 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
309 | return result; | |
310 | ||
311 | result = 2 / sqrt(shape); | |
312 | return result; | |
313 | } | |
314 | ||
315 | template <class RealType, class Policy> | |
316 | inline RealType kurtosis_excess(const gamma_distribution<RealType, Policy>& dist) | |
317 | { | |
318 | BOOST_MATH_STD_USING // for ADL of std functions | |
319 | ||
320 | static const char* function = "boost::math::kurtosis_excess(const gamma_distribution<%1%>&)"; | |
321 | ||
322 | RealType shape = dist.shape(); | |
323 | RealType scale = dist.scale(); | |
324 | ||
325 | RealType result = 0; | |
326 | if(false == detail::check_gamma(function, scale, shape, &result, Policy())) | |
327 | return result; | |
328 | ||
329 | result = 6 / shape; | |
330 | return result; | |
331 | } | |
332 | ||
333 | template <class RealType, class Policy> | |
334 | inline RealType kurtosis(const gamma_distribution<RealType, Policy>& dist) | |
335 | { | |
336 | return kurtosis_excess(dist) + 3; | |
337 | } | |
338 | ||
339 | } // namespace math | |
340 | } // namespace boost | |
341 | ||
342 | // This include must be at the end, *after* the accessors | |
343 | // for this distribution have been defined, in order to | |
344 | // keep compilers that support two-phase lookup happy. | |
345 | #include <boost/math/distributions/detail/derived_accessors.hpp> | |
346 | ||
347 | #endif // BOOST_STATS_GAMMA_HPP | |
348 | ||
349 |