1 //! Kernels 2 3 use crate::stats::float::Float; 4 5 /// Kernel function 6 pub trait Kernel<A>: Copy + Sync 7 where 8 A: Float, 9 { 10 /// Apply the kernel function to the given x-value. evaluate(&self, x: A) -> A11 fn evaluate(&self, x: A) -> A; 12 } 13 14 /// Gaussian kernel 15 #[derive(Clone, Copy)] 16 pub struct Gaussian; 17 18 impl<A> Kernel<A> for Gaussian 19 where 20 A: Float, 21 { evaluate(&self, x: A) -> A22 fn evaluate(&self, x: A) -> A { 23 use std::f32::consts::PI; 24 25 (x.powi(2).exp() * A::cast(2. * PI)).sqrt().recip() 26 } 27 } 28 29 #[cfg(test)] 30 macro_rules! test { 31 ($ty:ident) => { 32 mod $ty { 33 mod gaussian { 34 use approx::relative_eq; 35 use quickcheck::quickcheck; 36 use quickcheck::TestResult; 37 38 use crate::stats::univariate::kde::kernel::{Gaussian, Kernel}; 39 40 quickcheck! { 41 fn symmetric(x: $ty) -> bool { 42 x.is_nan() || relative_eq!(Gaussian.evaluate(-x), Gaussian.evaluate(x)) 43 } 44 } 45 46 // Any [a b] integral should be in the range [0 1] 47 quickcheck! { 48 fn integral(a: $ty, b: $ty) -> TestResult { 49 let a = a.sin().abs(); // map the value to [0 1] 50 let b = b.sin().abs(); // map the value to [0 1] 51 const DX: $ty = 1e-3; 52 53 if a > b { 54 TestResult::discard() 55 } else { 56 let mut acc = 0.; 57 let mut x = a; 58 let mut y = Gaussian.evaluate(a); 59 60 while x < b { 61 acc += DX * y / 2.; 62 63 x += DX; 64 y = Gaussian.evaluate(x); 65 66 acc += DX * y / 2.; 67 } 68 69 TestResult::from_bool( 70 (acc > 0. || relative_eq!(acc, 0.)) && 71 (acc < 1. || relative_eq!(acc, 1.))) 72 } 73 } 74 } 75 } 76 } 77 }; 78 } 79 80 #[cfg(test)] 81 mod test { 82 test!(f32); 83 test!(f64); 84 } 85