Support KDE in no_std

.github/workflows/test.yml

@@ -87,6 +87,9 @@ jobs:
       - name: Check no_std with nalgebra and rand
         run: cargo check --target riscv32imac-unknown-none-elf --no-default-features --features nalgebra,rand --lib
 
+      - name: Check no_std with KDE
+        run: cargo check --target riscv32imac-unknown-none-elf --no-default-features --features kde --lib
+
   features:
     needs: [clippy, fmt]
     runs-on: ubuntu-latest

Cargo.toml

@@ -31,11 +31,11 @@ required-features = ["rand", "std", "kde"]
 
 [features]
 default = ["std", "nalgebra", "rand"]
-std = ["approx/std", "num-traits/std", "nalgebra?/std", "rand?/std"]
+std = ["approx/std", "num-traits/std", "nalgebra?/std", "rand?/std", "kdtree?/std"]
 nalgebra = ["dep:nalgebra", "nalgebra/alloc", "nalgebra/libm"]
 rand = ["dep:rand", "nalgebra?/rand-no-std", "rand?/std_rng"]
 # kd-tree backed density estimation (src/density), implemented in terms of nalgebra vectors
-kde = ["dep:kdtree", "nalgebra", "std"]
+kde = ["dep:kdtree", "kdtree/libm", "nalgebra"]
 
 [dependencies]
 approx = { version = "0.5.0", default-features = false }
@@ -53,8 +53,11 @@ optional = true
 default-features = false
 
 [dependencies.kdtree]
-version = "0.7.0"
+version = "0.8.1"
+git = "https://github.com/mrhooray/kdtree-rs"
+rev = "db9af2f8e987c4f1c21251998f032b89a8bd2d52"
 optional = true
+default-features = false
 
 [dev-dependencies]
 criterion = "0.8"

README.md

@@ -55,6 +55,13 @@ matrix-backed distributions without `std`, enable `nalgebra` directly:
 statrs = { version = "*", default-features = false, features = ["nalgebra"] }
 ```
 
+To enable kernel density estimation without `std`, enable `kde`:
+
+```toml
+[dependencies]
+statrs = { version = "*", default-features = false, features = ["kde"] }
+```
+
 Heap-backed APIs use Rust's `alloc` crate. A `no_std` application that calls
 these APIs must provide and initialize a global allocator suitable for its
 target:

src/density/kde.rs

@@ -1,7 +1,9 @@
 use kdtree::distance::squared_euclidean;
+#[cfg(not(feature = "std"))]
+use num_traits::Float as _;
 
 use crate::{
-    density::{Container, DensityError, nearest_neighbors},
+    density::{Container, DensityError, nearest_neighbors, neighborhood_radius},
     function::kernel::{Gaussian, Kernel},
 };
 
@@ -30,7 +32,7 @@ where
     if neighbors.is_empty() {
         Err(DensityError::EmptyNeighborhood)
     } else {
-        let radius = neighbors.last().unwrap().sqrt(); // safe to unwrap here since `neighbors` is not empty
+        let radius = neighborhood_radius(&neighbors).unwrap();
         let d = x.length() as i32;
         Ok((1. / (n_samples * radius.powi(d)))
             * samples

src/density/knn.rs

@@ -1,9 +1,11 @@
 use super::Container;
 use crate::{
-    density::{DensityError, nearest_neighbors},
+    density::{DensityError, nearest_neighbors, neighborhood_radius},
     function::gamma::gamma,
 };
 use core::f64::consts::PI;
+#[cfg(not(feature = "std"))]
+use num_traits::Float as _;
 
 /// Computes the `k`-nearest neighbor density estimate for a given point `x`
 /// using the samples provided.
@@ -30,7 +32,7 @@ where
     if neighbors.is_empty() {
         Err(DensityError::EmptyNeighborhood)
     } else {
-        let radius = neighbors.last().unwrap().sqrt();
+        let radius = neighborhood_radius(&neighbors).unwrap();
         let d = x.length() as f64;
         Ok((k / n_samples) * (gamma(d / 2. + 1.) / (PI.powf(d / 2.) * radius.powf(d))))
     }

src/density/mod.rs

@@ -16,6 +16,8 @@ pub mod kde;
 pub mod knn;
 use alloc::vec::Vec;
 use kdtree::{ErrorKind, KdTree, distance::squared_euclidean};
+#[cfg(not(feature = "std"))]
+use num_traits::Float as _;
 use thiserror::Error;
 
 /// Errors that can occur when estimating a density from a sample.
@@ -79,6 +81,14 @@ impl_container!(
 );
 pub type NearestNeighbors = (Vec<f64>, f64);
 
+pub(crate) fn neighborhood_radius(neighbors: &[f64]) -> Option<f64> {
+    neighbors
+        .iter()
+        .copied()
+        .max_by(f64::total_cmp)
+        .map(f64::sqrt)
+}
+
 pub(crate) fn nearest_neighbors<S, X>(
     x: &X,
     samples: &S,
@@ -118,6 +128,11 @@ mod tests {
 
     use super::*;
 
+    #[test]
+    fn neighborhood_radius_uses_farthest_distance() {
+        assert_eq!(neighborhood_radius(&[4.0, 1.0, 9.0, 2.0]), Some(3.0));
+    }
+
     #[test]
     fn test_vec_container() {
         let v1 = vec![1.0, 2.0, 3.0];

tests/no_std/Cargo.toml

@@ -12,7 +12,7 @@ crate-type = ["cdylib"]
 
 [dependencies]
 dlmalloc = { version = "0.2.14", features = ["global"] }
-statrs = { path = "../..", default-features = false, features = ["nalgebra"] }
+statrs = { path = "../..", default-features = false, features = ["kde"] }
 
 [profile.dev]
 panic = "abort"

tests/no_std/src/lib.rs

@@ -3,6 +3,7 @@
 extern crate alloc;
 
 use alloc::vec;
+use statrs::density::{kde::kde_pdf, knn::knn_pdf};
 use statrs::distribution::{Categorical, Continuous, Empirical, Multinomial, MultivariateNormal};
 use statrs::generate::log_spaced;
 use statrs::statistics::{Data, Distribution, MeanN, OrderStatistics, RankTieBreaker};
@@ -22,6 +23,10 @@ fn assert_close(actual: f64, expected: f64, tolerance: f64) {
 
 #[unsafe(no_mangle)]
 pub extern "C" fn verify() {
+    let kde_samples = vec![[-1.0], [0.0], [1.0]];
+    assert!(kde_pdf(&[0.0], &kde_samples, Some(1.0)).unwrap() > 0.0);
+    assert!(knn_pdf(&[0.0], &kde_samples, Some(1.0)).unwrap() > 0.0);
+
     let categorical = Categorical::new(&[1.0, 2.0, 3.0]).unwrap();
     assert_close(categorical.mean().unwrap(), 4.0 / 3.0, 1e-12);