pearson 1.0.0

A Rust library implementing the Pearson correlation coefficient (PCC)
Documentation
#![no_std]

use core::marker::PhantomData;
use core::ops::{AddAssign, Mul};
use num_traits::{AsPrimitive, ConstZero, Float};

#[derive(Debug, Clone, Default, Copy)]
pub struct Pearson<T, O> {
    x: T,
    y: T,
    xy: T,
    x2: T,
    y2: T,
    len: usize,
    output: PhantomData<O>,
}

impl<T, O> Pearson<T, O>
where
    T: ConstZero,
{
    pub const fn new() -> Self {
        Self {
            x: T::ZERO,
            y: T::ZERO,
            xy: T::ZERO,
            x2: T::ZERO,
            y2: T::ZERO,
            len: 0,
            output: PhantomData,
        }
    }
}

impl<T, O> Pearson<T, O>
where
    T: AddAssign + Mul<Output = T> + Copy,
{
    pub fn record(&mut self, x: T, y: T) {
        self.x += x;
        self.y += y;
        self.xy += x * y;
        self.x2 += x * x;
        self.y2 += y * y;
        self.len += 1;
    }
}

impl<T, O> Pearson<T, O>
where
    O: Float + 'static,
    T: AsPrimitive<O> + ConstZero,
    usize: AsPrimitive<O>,
{
    pub fn correlate(&mut self) -> O {
        let x_prod_y = self.x.as_() * self.y.as_();
        let x_prod_x = self.x.as_() * self.x.as_();
        let y_prod_y = self.y.as_() * self.y.as_();
        let n = self.len.as_();

        let numerator = self.xy.as_() - (x_prod_y / n);
        let denominator =
            ((self.x2.as_() - (x_prod_x / n)) * (self.y2.as_() - (y_prod_y / n))).sqrt();

        *self = Self::new();

        numerator / denominator
    }
}

pub fn correlate<T, O>(x: &[T], y: &[T]) -> O
where
    O: Float + 'static,
    T: AsPrimitive<O> + ConstZero + AddAssign + Mul<Output = T> + Copy,
    usize: AsPrimitive<O>,
{
    assert_eq!(x.len(), y.len());
    let mut p = Pearson::new();
    x.iter().zip(y.iter()).for_each(|(&x, &y)| p.record(x, y));
    p.correlate()
}

#[test]
fn simple() {
    let x1 = [1, 2, 3, 4, 5];
    let x2 = [5, 4, 3, 2, 1];
    let y = [1, 2, 3, 4, 5];

    let c1: f32 = correlate(&x1, &y);
    let c2: f32 = correlate(&x2, &y);

    assert_eq!(c1, 1.0);
    assert_eq!(c2, -1.0);
}