probability_rs/dist/
lognormal.rs

1use crate::dist::normal::Normal;
2use crate::dist::{Continuous, DistError, Distribution, Moments};
3use crate::rng::RngCore;
4
5/// Lognormal with parameters (mu, sigma) where ln(X) ~ Normal(mu, sigma).
6#[derive(Debug, Clone, Copy)]
7pub struct LogNormal {
8    mu: f64,
9    sigma: f64,
10    normal: Normal,
11}
12
13impl LogNormal {
14    pub fn new(mu: f64, sigma: f64) -> Result<Self, DistError> {
15        if !(sigma > 0.0 && sigma.is_finite() && mu.is_finite()) {
16            return Err(DistError::InvalidParameter);
17        }
18        let normal = Normal::new(mu, sigma)?;
19        Ok(Self { mu, sigma, normal })
20    }
21    #[inline]
22    pub fn mu(&self) -> f64 {
23        self.mu
24    }
25    #[inline]
26    pub fn sigma(&self) -> f64 {
27        self.sigma
28    }
29}
30
31impl Distribution for LogNormal {
32    type Value = f64;
33    fn cdf(&self, x: f64) -> f64 {
34        if x <= 0.0 || !x.is_finite() {
35            return 0.0;
36        }
37        self.normal.cdf(x.ln())
38    }
39    fn in_support(&self, x: f64) -> bool {
40        x > 0.0 && x.is_finite()
41    }
42    fn sample<R: RngCore>(&self, rng: &mut R) -> f64 {
43        (self.normal.sample(rng)).exp()
44    }
45}
46
47impl Continuous for LogNormal {
48    fn pdf(&self, x: f64) -> f64 {
49        if !self.in_support(x) {
50            return 0.0;
51        }
52        let z = x.ln();
53        self.normal.pdf(z) / x
54    }
55    fn inv_cdf(&self, p: f64) -> f64 {
56        debug_assert!((0.0..=1.0).contains(&p));
57        self.normal.inv_cdf(p).exp()
58    }
59}
60
61impl Moments for LogNormal {
62    fn mean(&self) -> f64 {
63        (self.mu + 0.5 * self.sigma * self.sigma).exp()
64    }
65    fn variance(&self) -> f64 {
66        let s2 = self.sigma * self.sigma;
67        ((2.0 * s2).exp() - (s2).exp()) * (self.mu + s2).exp()
68    }
69        fn skewness(&self) -> f64 {
70            let s2 = self.sigma * self.sigma;
71            let es2 = s2.exp();
72            (es2 + 2.0) * (es2 - 1.0).sqrt()
73        }
74        fn kurtosis(&self) -> f64 {
75            let s2 = self.sigma * self.sigma;
76            let es2 = s2.exp();
77            es2.powi(4) + 2.0 * es2.powi(3) + 3.0 * es2.powi(2) - 6.0
78        }
79        fn entropy(&self) -> f64 {
80            // H = mu + 0.5 * ln(2π e σ^2)
81            self.mu + 0.5 * (2.0 * std::f64::consts::PI * std::f64::consts::E * self.sigma * self.sigma).ln()
82        }
83}
84
85#[cfg(test)]
86mod tests {
87    use super::*;
88    #[test]
89    fn pdf_cdf_basic() {
90        let ln = LogNormal::new(0.0, 1.0).unwrap();
91        assert!(ln.pdf(1.0) > 0.0);
92        assert!((ln.cdf(1.0) - 0.5).abs() < 2e-6);
93    }
94    #[test]
95    fn moments_higher() {
96        let ln = LogNormal::new(0.0_f64, 1.0_f64).unwrap();
97        let s2: f64 = 1.0;
98        let es2: f64 = s2.exp();
99        let skew = (es2 + 2.0) * (es2 - 1.0).sqrt();
100        let kurt = es2.powi(4) + 2.0*es2.powi(3) + 3.0*es2.powi(2) - 6.0;
101        assert!((ln.skewness() - skew).abs() < 1e-12);
102        assert!((ln.kurtosis() - kurt).abs() < 1e-12);
103    }
104}
probability_rs/dist/lognormal.rs

probability_rs/dist/
lognormal.rs
