Skip to main content

numeris/stats/
exponential.rs

1use crate::FloatScalar;
2use super::{ContinuousDistribution, StatsError};
3
4/// Exponential distribution with rate λ.
5///
6/// f(x) = λ exp(−λx) for x ≥ 0.
7///
8/// # Example
9///
10/// ```
11/// use numeris::stats::{Exponential, ContinuousDistribution};
12///
13/// let e = Exponential::new(2.0_f64).unwrap();
14/// assert!((e.mean() - 0.5).abs() < 1e-14);
15/// assert!((e.cdf(0.0)).abs() < 1e-14);
16/// ```
17#[derive(Debug, Clone, Copy)]
18pub struct Exponential<T> {
19    lambda: T,
20}
21
22impl<T: FloatScalar> Exponential<T> {
23    /// Create an exponential distribution with rate `lambda`. Requires `lambda > 0`.
24    pub fn new(lambda: T) -> Result<Self, StatsError> {
25        if lambda <= T::zero() {
26            return Err(StatsError::InvalidParameter);
27        }
28        Ok(Self { lambda })
29    }
30}
31
32impl<T: FloatScalar> Exponential<T> {
33    /// Draw a random sample from this distribution (inverse CDF method).
34    pub fn sample(&self, rng: &mut super::Rng) -> T {
35        let u: T = rng.next_float();
36        -(T::one() - u).ln() / self.lambda
37    }
38
39    /// Fill a fixed-size array with independent samples.
40    pub fn sample_array<const K: usize>(&self, rng: &mut super::Rng) -> [T; K] {
41        let mut out = [T::zero(); K];
42        for v in out.iter_mut() {
43            *v = self.sample(rng);
44        }
45        out
46    }
47}
48
49impl<T: FloatScalar> ContinuousDistribution<T> for Exponential<T> {
50    fn pdf(&self, x: T) -> T {
51        if x < T::zero() {
52            T::zero()
53        } else {
54            self.lambda * (-self.lambda * x).exp()
55        }
56    }
57
58    fn ln_pdf(&self, x: T) -> T {
59        if x < T::zero() {
60            T::neg_infinity()
61        } else {
62            self.lambda.ln() - self.lambda * x
63        }
64    }
65
66    fn cdf(&self, x: T) -> T {
67        if x <= T::zero() {
68            T::zero()
69        } else {
70            T::one() - (-self.lambda * x).exp()
71        }
72    }
73
74    fn quantile(&self, p: T) -> T {
75        -(T::one() - p).ln() / self.lambda
76    }
77
78    fn mean(&self) -> T {
79        T::one() / self.lambda
80    }
81
82    fn variance(&self) -> T {
83        T::one() / (self.lambda * self.lambda)
84    }
85}