statrs/distribution/erlang.rs
1use crate::distribution::{Continuous, ContinuousCDF, Gamma};
2use crate::statistics::*;
3use crate::Result;
4use rand::Rng;
5
6/// Implements the [Erlang](https://en.wikipedia.org/wiki/Erlang_distribution)
7/// distribution
8/// which is a special case of the
9/// [Gamma](https://en.wikipedia.org/wiki/Gamma_distribution)
10/// distribution
11///
12/// # Examples
13///
14/// ```
15/// use statrs::distribution::{Erlang, Continuous};
16/// use statrs::statistics::Distribution;
17/// use statrs::prec;
18///
19/// let n = Erlang::new(3, 1.0).unwrap();
20/// assert_eq!(n.mean().unwrap(), 3.0);
21/// assert!(prec::almost_eq(n.pdf(2.0), 0.270670566473225383788, 1e-15));
22/// ```
23#[derive(Debug, Copy, Clone, PartialEq)]
24pub struct Erlang {
25 g: Gamma,
26}
27
28impl Erlang {
29 /// Constructs a new erlang distribution with a shape (k)
30 /// of `shape` and a rate (λ) of `rate`
31 ///
32 /// # Errors
33 ///
34 /// Returns an error if `shape` or `rate` are `NaN`.
35 /// Also returns an error if `shape == 0` or `rate <= 0.0`
36 ///
37 /// # Examples
38 ///
39 /// ```
40 /// use statrs::distribution::Erlang;
41 ///
42 /// let mut result = Erlang::new(3, 1.0);
43 /// assert!(result.is_ok());
44 ///
45 /// result = Erlang::new(0, 0.0);
46 /// assert!(result.is_err());
47 /// ```
48 pub fn new(shape: u64, rate: f64) -> Result<Erlang> {
49 Gamma::new(shape as f64, rate).map(|g| Erlang { g })
50 }
51
52 /// Returns the shape (k) of the erlang distribution
53 ///
54 /// # Examples
55 ///
56 /// ```
57 /// use statrs::distribution::Erlang;
58 ///
59 /// let n = Erlang::new(3, 1.0).unwrap();
60 /// assert_eq!(n.shape(), 3);
61 /// ```
62 pub fn shape(&self) -> u64 {
63 self.g.shape() as u64
64 }
65
66 /// Returns the rate (λ) of the erlang distribution
67 ///
68 /// # Examples
69 ///
70 /// ```
71 /// use statrs::distribution::Erlang;
72 ///
73 /// let n = Erlang::new(3, 1.0).unwrap();
74 /// assert_eq!(n.rate(), 1.0);
75 /// ```
76 pub fn rate(&self) -> f64 {
77 self.g.rate()
78 }
79}
80
81impl ::rand::distributions::Distribution<f64> for Erlang {
82 fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> f64 {
83 ::rand::distributions::Distribution::sample(&self.g, rng)
84 }
85}
86
87impl ContinuousCDF<f64, f64> for Erlang {
88 /// Calculates the cumulative distribution function for the erlang
89 /// distribution
90 /// at `x`
91 ///
92 /// # Formula
93 ///
94 /// ```ignore
95 /// γ(k, λx) (k - 1)!
96 /// ```
97 ///
98 /// where `k` is the shape, `λ` is the rate, and `γ` is the lower
99 /// incomplete gamma function
100 fn cdf(&self, x: f64) -> f64 {
101 self.g.cdf(x)
102 }
103
104 /// Calculates the cumulative distribution function for the erlang
105 /// distribution
106 /// at `x`
107 ///
108 /// # Formula
109 ///
110 /// ```ignore
111 /// γ(k, λx) (k - 1)!
112 /// ```
113 ///
114 /// where `k` is the shape, `λ` is the rate, and `γ` is the upper
115 /// incomplete gamma function
116 fn sf(&self, x: f64) -> f64 {
117 self.g.sf(x)
118 }
119}
120
121impl Min<f64> for Erlang {
122 /// Returns the minimum value in the domain of the
123 /// erlang distribution representable by a double precision
124 /// float
125 ///
126 /// # Formula
127 ///
128 /// ```ignore
129 /// 0
130 /// ```
131 fn min(&self) -> f64 {
132 self.g.min()
133 }
134}
135
136impl Max<f64> for Erlang {
137 /// Returns the maximum value in the domain of the
138 /// erlang distribution representable by a double precision
139 /// float
140 ///
141 /// # Formula
142 ///
143 /// ```ignore
144 /// INF
145 /// ```
146 fn max(&self) -> f64 {
147 self.g.max()
148 }
149}
150
151impl Distribution<f64> for Erlang {
152 /// Returns the mean of the erlang distribution
153 ///
154 /// # Remarks
155 ///
156 /// Returns `shape` if `rate == f64::INFINITY`. This behavior
157 /// is borrowed from the Math.NET implementation
158 ///
159 /// # Formula
160 ///
161 /// ```ignore
162 /// k / λ
163 /// ```
164 ///
165 /// where `k` is the shape and `λ` is the rate
166 fn mean(&self) -> Option<f64> {
167 self.g.mean()
168 }
169 /// Returns the variance of the erlang distribution
170 ///
171 /// # Formula
172 ///
173 /// ```ignore
174 /// k / λ^2
175 /// ```
176 ///
177 /// where `α` is the shape and `λ` is the rate
178 fn variance(&self) -> Option<f64> {
179 self.g.variance()
180 }
181 /// Returns the entropy of the erlang distribution
182 ///
183 /// # Formula
184 ///
185 /// ```ignore
186 /// k - ln(λ) + ln(Γ(k)) + (1 - k) * ψ(k)
187 /// ```
188 ///
189 /// where `k` is the shape, `λ` is the rate, `Γ` is the gamma function,
190 /// and `ψ` is the digamma function
191 fn entropy(&self) -> Option<f64> {
192 self.g.entropy()
193 }
194 /// Returns the skewness of the erlang distribution
195 ///
196 /// # Formula
197 ///
198 /// ```ignore
199 /// 2 / sqrt(k)
200 /// ```
201 ///
202 /// where `k` is the shape
203 fn skewness(&self) -> Option<f64> {
204 self.g.skewness()
205 }
206}
207
208impl Mode<Option<f64>> for Erlang {
209 /// Returns the mode for the erlang distribution
210 ///
211 /// # Remarks
212 ///
213 /// Returns `shape` if `rate ==f64::INFINITY`. This behavior
214 /// is borrowed from the Math.NET implementation
215 ///
216 /// # Formula
217 ///
218 /// ```ignore
219 /// (k - 1) / λ
220 /// ```
221 ///
222 /// where `k` is the shape and `λ` is the rate
223 fn mode(&self) -> Option<f64> {
224 self.g.mode()
225 }
226}
227
228impl Continuous<f64, f64> for Erlang {
229 /// Calculates the probability density function for the erlang distribution
230 /// at `x`
231 ///
232 /// # Remarks
233 ///
234 /// Returns `NAN` if any of `shape` or `rate` are `INF`
235 /// or if `x` is `INF`
236 ///
237 /// # Formula
238 ///
239 /// ```ignore
240 /// (λ^k / Γ(k)) * x^(k - 1) * e^(-λ * x)
241 /// ```
242 ///
243 /// where `k` is the shape, `λ` is the rate, and `Γ` is the gamma function
244 fn pdf(&self, x: f64) -> f64 {
245 self.g.pdf(x)
246 }
247
248 /// Calculates the log probability density function for the erlang
249 /// distribution
250 /// at `x`
251 ///
252 /// # Remarks
253 ///
254 /// Returns `NAN` if any of `shape` or `rate` are `INF`
255 /// or if `x` is `INF`
256 ///
257 /// # Formula
258 ///
259 /// ```ignore
260 /// ln((λ^k / Γ(k)) * x^(k - 1) * e ^(-λ * x))
261 /// ```
262 ///
263 /// where `k` is the shape, `λ` is the rate, and `Γ` is the gamma function
264 fn ln_pdf(&self, x: f64) -> f64 {
265 self.g.ln_pdf(x)
266 }
267}
268
269#[rustfmt::skip]
270#[cfg(all(test, feature = "nightly"))]
271mod tests {
272 use crate::distribution::Erlang;
273 use crate::distribution::internal::*;
274 use crate::consts::ACC;
275
276 fn try_create(shape: u64, rate: f64) -> Erlang {
277 let n = Erlang::new(shape, rate);
278 assert!(n.is_ok());
279 n.unwrap()
280 }
281
282 fn create_case(shape: u64, rate: f64) {
283 let n = try_create(shape, rate);
284 assert_eq!(shape, n.shape());
285 assert_eq!(rate, n.rate());
286 }
287
288 fn bad_create_case(shape: u64, rate: f64) {
289 let n = Erlang::new(shape, rate);
290 assert!(n.is_err());
291 }
292
293 #[test]
294 fn test_create() {
295 create_case(1, 0.1);
296 create_case(1, 1.0);
297 create_case(10, 10.0);
298 create_case(10, 1.0);
299 create_case(10, f64::INFINITY);
300 }
301
302 #[test]
303 fn test_bad_create() {
304 bad_create_case(0, 1.0);
305 bad_create_case(1, 0.0);
306 bad_create_case(1, f64::NAN);
307 bad_create_case(1, -1.0);
308 }
309
310 #[test]
311 fn test_continuous() {
312 test::check_continuous_distribution(&try_create(1, 2.5), 0.0, 20.0);
313 test::check_continuous_distribution(&try_create(2, 1.5), 0.0, 20.0);
314 test::check_continuous_distribution(&try_create(3, 0.5), 0.0, 20.0);
315 }
316}