Skip to content

Commit

Permalink
feat: enhance Levy distribution with comprehensive documentation and …
Browse files Browse the repository at this point in the history
…methods
  • Loading branch information
Qazalbash committed Jan 16, 2025
1 parent f82ee4c commit 347b4b6
Showing 1 changed file with 171 additions and 0 deletions.
171 changes: 171 additions & 0 deletions src/distribution/levy.rs
Original file line number Diff line number Diff line change
Expand Up @@ -4,12 +4,24 @@ use crate::function::erf::{erf, erfc, erfc_inv};
use crate::statistics::*;
use core::f64;

/// Implements the [Levy](https://en.wikipedia.org/wiki/L%C3%A9vy_distribution) distribution.
///
/// # Example
///
/// ```
/// use statrs::distribution::{Levy, Continuous};
/// use statrs::statistics::Distribution;
///
/// let n = Levy::new(1.0, 1.0).unwrap();
/// assert_eq!(n.pdf(0.0), 0.0);
/// ```
#[derive(Copy, Clone, PartialEq, Debug)]
pub struct Levy {
mu: f64,
c: f64,
}

/// Represents the errors that can occur when creating a [`Levy`].
#[derive(Copy, Clone, PartialEq, Eq, Debug, Hash)]
#[non_exhaustive]
pub enum LevyError {
Expand All @@ -32,6 +44,23 @@ impl std::fmt::Display for LevyError {
impl std::error::Error for LevyError {}

impl Levy {
/// Constructs a new Levy distribution with a location (μ) and dispersion (c)
///
/// # Errors
///
/// Returns and error if `mu` is NaN or infinite or if `c` is NaN, infinite or nonpositive
///
/// # Example
///
/// ```
/// use statrs::distribution::Levy;
///
/// let mut result = Levy::new(0.0, 1.0);
/// assert!(result.is_ok());
///
/// result = Levy::new(0.0, 0.0);
/// assert!(result.is_err());
/// ```
pub fn new(mu: f64, c: f64) -> Result<Levy, LevyError> {
if mu.is_nan() || mu.is_infinite() {
return Err(LevyError::LocationInvalid);
Expand All @@ -42,10 +71,30 @@ impl Levy {
Ok(Levy { mu, c })
}

/// Returns the location (μ) of the Levy distribution
///
/// # Example
///
/// ```
/// use statrs::distribution::Levy;
///
/// let n = Levy::new(1.0, 1.0).unwrap();
/// assert_eq!(n.mu(), 1.0);
/// ```
pub fn mu(&self) -> f64 {
self.mu
}

/// Returns the dispersion (c) of the Levy distribution
///
/// # Example
///
/// ```
/// use statrs::distribution::Levy;
///
/// let n = Levy::new(1.0, 1.0).unwrap();
/// assert_eq!(n.c(), 1.0);
/// ```
pub fn c(&self) -> f64 {
self.c
}
Expand All @@ -70,6 +119,17 @@ impl ::rand::distributions::Distribution<f64> for Levy {
}

impl ContinuousCDF<f64, f64> for Levy {
/// Calculates the cumulative distribution function for the Levy distribution at `x`
///
/// # Formula
///
/// ```text
/// 0 if x <= μ
/// erfc(sqrt(c / (2 * (x - μ)))) if x > μ
/// ```
///
/// where `μ` is the location, `c` is the dispersion, and `erfc` is the
/// complementary error function.
fn cdf(&self, x: f64) -> f64 {
if x <= self.mu {
0.0
Expand All @@ -80,6 +140,17 @@ impl ContinuousCDF<f64, f64> for Levy {
}
}

/// Calculates the survival function for the Levy distribution at `x`
///
/// # Formula
///
/// ```text
/// 1 if x <= μ
/// erf(sqrt(c / (2 * (x - μ)))) if x > μ
/// ```
///
/// where `μ` is the location, `c` is the dispersion, and `erf` is the error
/// function.
fn sf(&self, x: f64) -> f64 {
if x <= self.mu {
1.0
Expand All @@ -90,6 +161,21 @@ impl ContinuousCDF<f64, f64> for Levy {
}
}

/// Calculates the inverse cumulative distribution function for the
/// normal distribution at `x`.
///
/// # Panics
///
/// If `x < 0.0` or `x > 1.0`
///
/// # Formula
///
/// ```text
/// μ + c * (erfc_inv(x)^2)/2
/// ```
///
/// where `μ` is the mean, `σ` is the standard deviation and `erfc_inv` is
/// the inverse of the complementary error function
fn inverse_cdf(&self, x: f64) -> f64 {
if !(0.0..=1.0).contains(&x) {
panic!("x must be in [0, 1]");
Expand All @@ -100,26 +186,74 @@ impl ContinuousCDF<f64, f64> for Levy {
}

impl Min<f64> for Levy {
/// Returns the minimum value in the domain of the
/// Levy distribution representable by a double precision float
///
/// # Formula
///
/// ```text
/// μ
/// ```
fn min(&self) -> f64 {
self.mu
}
}

impl Max<f64> for Levy {
/// Returns the maximum value in the domain of the
/// Levy distribution representable by a double precision float
///
/// # Formula
///
/// ```text
/// f64::INFINITY
/// ```
fn max(&self) -> f64 {
f64::INFINITY
}
}

impl Distribution<f64> for Levy {
/// Returns the mean of the Levy distribution
///
/// # Formula
///
/// ```text
/// f64::INFINITY
/// ```
fn mean(&self) -> Option<f64> {
Some(f64::INFINITY)
}

/// Returns the variance of the Levy distribution
///
/// # Formula
///
/// ```text
/// f64::INFINITY
/// ```
fn variance(&self) -> Option<f64> {
Some(f64::INFINITY)
}

/// Returns the standard deviation of the Levy distribution
///
/// # Formula
///
/// ```text
/// f64::INFINITY
/// ```
fn std_dev(&self) -> Option<f64> {
Some(f64::INFINITY)
}

/// Returns the entropy of the Levy distribution
///
/// # Formula
///
/// ```text
/// (1 + 3γ + ln(16πc^2))/2
/// ```
fn entropy(&self) -> Option<f64> {
/// CONSTANT_PART = 1.5 * EULER_MASCHERONI + 0.5 * (1.0 + LN_PI + 16.0_f64.ln())
const CONSTANT_PART: f64 = 3.32448280139688989720525569282472133636474609375;
Expand All @@ -128,18 +262,46 @@ impl Distribution<f64> for Levy {
}

impl Median<f64> for Levy {
/// Returns the median of the Levy distribution
///
/// # Formula
///
/// ```text
/// μ + c/(2 * erfc_inv(0.5)^2)
/// ```
///
/// where `μ` is the mean, `c` is the dispersion and `erfc_inv` is
/// the inverse of the complementary error function.
fn median(&self) -> f64 {
self.mu + self.c * 0.5 * erfc_inv(0.5).powf(-2.0)
}
}

impl Mode<Option<f64>> for Levy {
/// Returns the mode of the Levy distribution
///
/// # Formula
///
/// ```text
/// μ + c/3
/// ```
///
/// where `μ` is the mean and `c` is the dispersion.
fn mode(&self) -> Option<f64> {
Some(self.mu + self.c / 3.0)
}
}

impl Continuous<f64, f64> for Levy {
/// Calculates the probability density function for the Levy distribution at `x`
///
/// # Formula
///
/// ```text
/// (sqrt(c / 2 * π)) * e^(-1/2 * (c / (x - μ))) * (1 / (x - μ)^(3/2))
/// ```
///
/// where `μ` is the mean and `c` is the dispersion.
fn pdf(&self, x: f64) -> f64 {
if x <= self.mu {
0.0
Expand All @@ -149,6 +311,15 @@ impl Continuous<f64, f64> for Levy {
}
}

/// Calculates the log probability density function for the Levy distribution at `x`
///
/// # Formula
///
/// ```text
/// 1/2 * (ln(c) - ln(2 * π) - (c / (x - μ))) - 3/2 * ln(x - μ)
/// ```
///
/// where `μ` is the mean and `σ` is the standard deviation
fn ln_pdf(&self, x: f64) -> f64 {
use crate::consts::LN_SQRT_2PI;

Expand Down

0 comments on commit 347b4b6

Please sign in to comment.