rand_distr/
gumbel.rs

1// Copyright 2021 Developers of the Rand project.
2//
3// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
4// https://www.apache.org/licenses/LICENSE-2.0> or the MIT license
5// <LICENSE-MIT or https://opensource.org/licenses/MIT>, at your
6// option. This file may not be copied, modified, or distributed
7// except according to those terms.
8
9//! The Gumbel distribution `Gumbel(μ, β)`.
10
11use crate::{Distribution, OpenClosed01};
12use core::fmt;
13use num_traits::Float;
14use rand::Rng;
15
16/// The [Gumbel distribution](https://en.wikipedia.org/wiki/Gumbel_distribution) `Gumbel(μ, β)`.
17///
18/// The Gumbel distribution is a continuous probability distribution
19/// with location parameter `μ` (`mu`) and scale parameter `β` (`beta`).
20/// It is used to model the distribution of the maximum (or minimum)
21/// of a number of samples of various distributions.
22///
23/// # Density function
24///
25/// `f(x) = exp(-(z + exp(-z))) / β`, where `z = (x - μ) / β`.
26///
27/// # Plot
28///
29/// The following plot illustrates the Gumbel distribution with various values of `μ` and `β`.
30/// Note how the location parameter `μ` shifts the distribution along the x-axis,
31/// and the scale parameter `β` changes the density around `μ`.
32/// Note also the asymptotic behavior of the distribution towards the right.
33///
34/// ![Gumbel distribution](https://raw.githubusercontent.com/rust-random/charts/main/charts/gumbel.svg)
35///
36/// # Example
37/// ```
38/// use rand::prelude::*;
39/// use rand_distr::Gumbel;
40///
41/// let val: f64 = rand::rng().sample(Gumbel::new(0.0, 1.0).unwrap());
42/// println!("{}", val);
43/// ```
44#[derive(Clone, Copy, Debug, PartialEq)]
45#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
46pub struct Gumbel<F>
47where
48    F: Float,
49    OpenClosed01: Distribution<F>,
50{
51    location: F,
52    scale: F,
53}
54
55/// Error type returned from [`Gumbel::new`].
56#[derive(Clone, Copy, Debug, PartialEq, Eq)]
57pub enum Error {
58    /// location is infinite or NaN
59    LocationNotFinite,
60    /// scale is not finite positive number
61    ScaleNotPositive,
62}
63
64impl fmt::Display for Error {
65    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
66        f.write_str(match self {
67            Error::ScaleNotPositive => "scale is not positive and finite in Gumbel distribution",
68            Error::LocationNotFinite => "location is not finite in Gumbel distribution",
69        })
70    }
71}
72
73#[cfg(feature = "std")]
74impl std::error::Error for Error {}
75
76impl<F> Gumbel<F>
77where
78    F: Float,
79    OpenClosed01: Distribution<F>,
80{
81    /// Construct a new `Gumbel` distribution with given `location` and `scale`.
82    pub fn new(location: F, scale: F) -> Result<Gumbel<F>, Error> {
83        if scale <= F::zero() || scale.is_infinite() || scale.is_nan() {
84            return Err(Error::ScaleNotPositive);
85        }
86        if location.is_infinite() || location.is_nan() {
87            return Err(Error::LocationNotFinite);
88        }
89        Ok(Gumbel { location, scale })
90    }
91}
92
93impl<F> Distribution<F> for Gumbel<F>
94where
95    F: Float,
96    OpenClosed01: Distribution<F>,
97{
98    fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> F {
99        let x: F = rng.sample(OpenClosed01);
100        self.location - self.scale * (-x.ln()).ln()
101    }
102}
103
104#[cfg(test)]
105mod tests {
106    use super::*;
107
108    #[test]
109    #[should_panic]
110    fn test_zero_scale() {
111        Gumbel::new(0.0, 0.0).unwrap();
112    }
113
114    #[test]
115    #[should_panic]
116    fn test_infinite_scale() {
117        Gumbel::new(0.0, f64::INFINITY).unwrap();
118    }
119
120    #[test]
121    #[should_panic]
122    fn test_nan_scale() {
123        Gumbel::new(0.0, f64::NAN).unwrap();
124    }
125
126    #[test]
127    #[should_panic]
128    fn test_infinite_location() {
129        Gumbel::new(f64::INFINITY, 1.0).unwrap();
130    }
131
132    #[test]
133    #[should_panic]
134    fn test_nan_location() {
135        Gumbel::new(f64::NAN, 1.0).unwrap();
136    }
137
138    #[test]
139    fn test_sample_against_cdf() {
140        fn neg_log_log(x: f64) -> f64 {
141            -(-x.ln()).ln()
142        }
143        let location = 0.0;
144        let scale = 1.0;
145        let iterations = 100_000;
146        let increment = 1.0 / iterations as f64;
147        let probabilities = [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9];
148        let mut quantiles = [0.0; 9];
149        for (i, p) in probabilities.iter().enumerate() {
150            quantiles[i] = neg_log_log(*p);
151        }
152        let mut proportions = [0.0; 9];
153        let d = Gumbel::new(location, scale).unwrap();
154        let mut rng = crate::test::rng(1);
155        for _ in 0..iterations {
156            let replicate = d.sample(&mut rng);
157            for (i, q) in quantiles.iter().enumerate() {
158                if replicate < *q {
159                    proportions[i] += increment;
160                }
161            }
162        }
163        assert!(proportions
164            .iter()
165            .zip(&probabilities)
166            .all(|(p_hat, p)| (p_hat - p).abs() < 0.003))
167    }
168
169    #[test]
170    fn gumbel_distributions_can_be_compared() {
171        assert_eq!(Gumbel::new(1.0, 2.0), Gumbel::new(1.0, 2.0));
172    }
173}
rand_distr/gumbel.rs

rand_distr/
gumbel.rs
