1//! Kernels
2
3use crate::stats::float::Float;
4
5/// Kernel function
6pub trait Kernel<A>: Copy + Sync
7where
8 A: Float,
9{
10 /// Apply the kernel function to the given x-value.
11 fn evaluate(&self, x: A) -> A;
12}
13
14/// Gaussian kernel
15#[derive(Clone, Copy)]
16pub struct Gaussian;
17
18impl<A> Kernel<A> for Gaussian
19where
20 A: Float,
21{
22 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)]
30macro_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)]
81mod test {
82 test!(f32);
83 test!(f64);
84}
85