From 22a2884f2d333097e363758def387f0719d2dabf Mon Sep 17 00:00:00 2001 From: julianknodt Date: Mon, 18 May 2020 02:07:56 -0700 Subject: [PATCH 1/2] Created Dirac distribution The dirac distribution is a trivial distribution but is useful for modelling certainty in situations. This implementation is fairly simple in that most outputs are either 0, 1, Infinity, or the point of the distribution. --- src/distribution/dirac.rs | 324 ++++++++++++++++++++++++++++++++++++++ src/distribution/mod.rs | 2 + 2 files changed, 326 insertions(+) create mode 100644 src/distribution/dirac.rs diff --git a/src/distribution/dirac.rs b/src/distribution/dirac.rs new file mode 100644 index 00000000..524c83ed --- /dev/null +++ b/src/distribution/dirac.rs @@ -0,0 +1,324 @@ +use crate::distribution::{Continuous, Univariate}; +use crate::statistics::*; +use crate::{Result, StatsError}; +use rand::distributions::Distribution; +use rand::Rng; +use std::f64; + +/// Implements the [Dirac Delta](https://en.wikipedia.org/wiki/Dirac_delta_function#As_a_distribution) +/// distribution +/// +/// # Examples +/// +/// ``` +/// use statrs::distribution::{Dirac}; +/// +/// let n = Dirac::new(3.0).unwrap(); +/// assert_eq!(n.mean(), 3.0); +/// assert_eq!(n.pdf(1.0), 0.0); +/// ``` +#[derive(Debug, Copy, Clone, PartialEq)] +pub struct Dirac(f64); + +impl Dirac { + /// Constructs a new dirac distribution function at value `v`. + /// + /// # Errors + /// + /// Returns an error if `v` is not-a-number. + /// + /// # Examples + /// + /// ``` + /// use statrs::distribution::Dirac; + /// + /// let mut result = Dirac::new(0.0); + /// assert!(result.is_ok()); + /// + /// result = Dirac::new(f64::NAN); + /// assert!(result.is_err()); + /// ``` + pub fn new(v: f64) -> Result { + if v.is_nan() { + Err(StatsError::BadParams) + } else { + Ok(Dirac(v)) + } + } +} + +impl Distribution for Dirac { + fn sample(&self, _: &mut R) -> f64 { + self.0 + } +} + +impl Univariate for Dirac { + /// Calculates the cumulative distribution function for the + /// dirac distribution at `x` + /// + /// Where the value is 1 if x > `v`, 0 otherwise. + /// + fn cdf(&self, x: f64) -> f64 { + if x < self.0 { 0.0 } else { 1.0 } + } +} + +impl Min for Dirac { + /// Returns the minimum value in the domain of the + /// dirac distribution representable by a double precision float + /// + /// # Formula + /// + /// ```ignore + /// v + /// ``` + fn min(&self) -> f64 { self.0 } +} + +impl Max for Dirac { + /// Returns the maximum value in the domain of the + /// dirac distribution representable by a double precision float + /// + /// # Formula + /// + /// ```ignore + /// v + /// ``` + fn max(&self) -> f64 { self.0 } +} + +impl Mean for Dirac { + /// Returns the mean of the dirac distribution + /// + /// # Remarks + /// + /// Since the only value that can be produced by this distribution is `v` with probability + /// 1, it is just `v`. + fn mean(&self) -> f64 { self.0 } +} + +impl Variance for Dirac { + /// Returns the variance of the normal distribution + /// + /// # Formula + /// + /// ```ignore + /// 0 + /// ``` + /// + /// Since only one value can be produced there is no variance. + fn variance(&self) -> f64 { 0.0 } + + /// Returns the standard deviation of the normal distribution + /// + /// # Remarks + /// + /// Since there is no variance in draws from this distribution the standard deviation is + /// also 0. + fn std_dev(&self) -> f64 { 0.0 } +} + +impl Entropy for Dirac { + /// Returns the entropy of the dirac distribution + /// + /// # Formula + /// + /// ```ignore + /// 0 + /// ``` + /// + /// Since this distribution has full certainty, it encodes no information + fn entropy(&self) -> f64 { 0.0 } +} + +impl Skewness for Dirac { + /// Returns the skewness of the dirac distribution + /// + /// # Formula + /// + /// ```ignore + /// 0 + /// ``` + fn skewness(&self) -> f64 { 0.0 } +} + +impl Median for Dirac { + /// Returns the median of the dirac distribution + /// + /// # Formula + /// + /// ```ignore + /// v + /// ``` + /// + /// where `v` is the point of the dirac distribution + fn median(&self) -> f64 { + self.0 + } +} + +impl Mode for Dirac { + /// Returns the mode of the normal distribution + /// + /// # Formula + /// + /// ```ignore + /// v + /// ``` + /// + /// where `v` is the point of the dirac distribution + fn mode(&self) -> f64 { + self.0 + } +} + +impl Continuous for Dirac { + /// Calculates the probability density function for the dirac distribution + /// at `x` + /// + /// # Formula + /// + /// ```ignore + /// 1 if x = v, 0 otherwise + /// ``` + /// + /// where `v` is point of this dirac distribution + fn pdf(&self, x: f64) -> f64 { + if x == self.0 { 1.0 } else { 0.0 } + } + + /// Calculates the log probability density function for the normal + /// distribution + /// at `x` + /// + /// # Formula + /// + /// ```ignore + /// ln(1 if x = v, 0 otherwise) + /// ``` + /// + /// where `v` is the point of this dirac distribution + /// + /// # Remarks + /// This distribution is usually negative infinity everywhere except at `v`. + fn ln_pdf(&self, x: f64) -> f64 { + if self.0 == x { 0.0 } else { f64::NEG_INFINITY } + } +} + +#[cfg_attr(rustfmt, rustfmt_skip)] +#[cfg(test)] +mod test { + use std::f64; + use crate::statistics::*; + use crate::distribution::{Univariate, Continuous, Dirac}; + use crate::distribution::internal::*; + + fn try_create(v: f64) -> Dirac { + let d = Dirac::new(v); + assert!(d.is_ok()); + d.unwrap() + } + + fn create_case(v: f64) { + let d = try_create(v); + assert_eq!(v, d.mean()); + } + + fn bad_create_case(v: f64) { + let d = Dirac::new(v); + assert!(d.is_err()); + } + + fn test_case(v: f64, expected: f64, eval: F) + where F: Fn(Dirac) -> f64 + { + let x = eval(try_create(v)); + assert_eq!(expected, x); + } + + #[test] + fn test_create() { + create_case(10.0); + create_case(-5.0); + create_case(10.0); + create_case(100.0); + create_case(f64::INFINITY); + } + + #[test] + fn test_bad_create() { + bad_create_case(f64::NAN); + } + + #[test] + fn test_variance() { + test_case(0.0, 0.0, |x| x.variance()); + test_case(-5.0, 0.0, |x| x.variance()); + test_case(f64::INFINITY, 0.0, |x| x.variance()); + } + + #[test] + fn test_entropy() { + test_case(0.0, 0.0, |x| x.entropy()); + test_case(f64::INFINITY, 0.0, |x| x.entropy()); + } + + #[test] + fn test_skewness() { + test_case(0.0, 0.0, |x| x.skewness()); + test_case(4.0, 0.0, |x| x.skewness()); + test_case(0.3, 0.0, |x| x.skewness()); + test_case(f64::INFINITY, 0.0, |x| x.skewness()); + } + + #[test] + fn test_mode() { + test_case(0.0, 0.0, |x| x.mode()); + test_case(3.0, 3.0, |x| x.mode()); + test_case(f64::INFINITY, f64::INFINITY, |x| x.mode()); + } + + #[test] + fn test_median() { + test_case(0.0, 0.0, |x| x.median()); + test_case(3.0, 3.0, |x| x.median()); + test_case(f64::INFINITY, f64::INFINITY, |x| x.median()); + } + + #[test] + fn test_min_max() { + test_case(0.0, 0.0, |x| x.min()); + test_case(3.0, 3.0, |x| x.min()); + test_case(f64::INFINITY, f64::INFINITY, |x| x.min()); + + test_case(0.0, 0.0, |x| x.max()); + test_case(3.0, 3.0, |x| x.max()); + test_case(f64::NEG_INFINITY, f64::NEG_INFINITY, |x| x.max()); + } + + #[test] + fn test_pdf() { + test_case(0.0, 0.0, |x| x.pdf(1.0)); + test_case(3.0, 1.0, |x| x.pdf(3.0)); + test_case(f64::NEG_INFINITY, 0.0, |x| x.pdf(1.0)); + test_case(f64::NEG_INFINITY, 1.0, |x| x.pdf(f64::NEG_INFINITY)); + } + + #[test] + fn test_ln_pdf() { + test_case(0.0, 0.0, |x| x.ln_pdf(0.0)); + test_case(3.0, 0.0, |x| x.ln_pdf(3.0)); + test_case(f64::INFINITY, f64::NEG_INFINITY, |x| x.ln_pdf(1.0)); + test_case(f64::INFINITY, 0.0, |x| x.ln_pdf(f64::INFINITY)); + } + + #[test] + fn test_cdf() { + test_case(0.0, 1.0, |x| x.cdf(0.0)); + test_case(3.0, 1.0, |x| x.cdf(3.0)); + test_case(f64::INFINITY, 0.0, |x| x.cdf(1.0)); + test_case(f64::INFINITY, 1.0, |x| x.cdf(f64::INFINITY)); + } +} diff --git a/src/distribution/mod.rs b/src/distribution/mod.rs index 0c30279e..83a02ae9 100644 --- a/src/distribution/mod.rs +++ b/src/distribution/mod.rs @@ -9,6 +9,7 @@ pub use self::categorical::Categorical; pub use self::cauchy::Cauchy; pub use self::chi::Chi; pub use self::chi_squared::ChiSquared; +pub use self::dirac::Dirac; pub use self::dirichlet::Dirichlet; pub use self::discrete_uniform::DiscreteUniform; pub use self::erlang::Erlang; @@ -37,6 +38,7 @@ mod categorical; mod cauchy; mod chi; mod chi_squared; +mod dirac; mod dirichlet; mod discrete_uniform; mod erlang; From a9a41aabeb9da0004be8829b801ad8113fb974f2 Mon Sep 17 00:00:00 2001 From: julianknodt Date: Thu, 21 May 2020 17:29:40 -0700 Subject: [PATCH 2/2] Changed some leftover normals from documentation --- src/distribution/dirac.rs | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/distribution/dirac.rs b/src/distribution/dirac.rs index 524c83ed..62aba625 100644 --- a/src/distribution/dirac.rs +++ b/src/distribution/dirac.rs @@ -99,7 +99,7 @@ impl Mean for Dirac { } impl Variance for Dirac { - /// Returns the variance of the normal distribution + /// Returns the variance of the dirac distribution /// /// # Formula /// @@ -110,7 +110,7 @@ impl Variance for Dirac { /// Since only one value can be produced there is no variance. fn variance(&self) -> f64 { 0.0 } - /// Returns the standard deviation of the normal distribution + /// Returns the standard deviation of the dirac distribution /// /// # Remarks /// @@ -159,7 +159,7 @@ impl Median for Dirac { } impl Mode for Dirac { - /// Returns the mode of the normal distribution + /// Returns the mode of the dirac distribution /// /// # Formula /// @@ -188,7 +188,7 @@ impl Continuous for Dirac { if x == self.0 { 1.0 } else { 0.0 } } - /// Calculates the log probability density function for the normal + /// Calculates the log probability density function for the dirac /// distribution /// at `x` ///