plotters-statistical 0.1.0

//! Histogram binning with several standard bin-count rules.

use super::{quartiles, sorted_finite, StatsError};

/// How to choose the number/width of histogram bins.
#[derive(Debug, Clone, Copy, PartialEq)]
pub enum BinRule {
    /// A fixed number of equal-width bins.
    Count(usize),
    /// A fixed bin width (bin count derived from the data range).
    Width(f64),
    /// Sturges' rule: `ceil(log2(n)) + 1`.
    Sturges,
    /// Freedman–Diaconis rule: width `= 2 * IQR * n^(-1/3)`.
    FreedmanDiaconis,
    /// Scott's rule: width `= 3.49 * std * n^(-1/3)`.
    Scott,
}

/// A computed histogram: bin edges (length `bins + 1`), per-bin counts, and
/// per-bin density (`count / (n * width)`, so the areas sum to 1).
#[derive(Debug, Clone, PartialEq)]
pub struct Histogram {
    /// Bin edges, ascending; `edges.len() == counts.len() + 1`.
    pub edges: Vec<f64>,
    /// Count of samples in each bin.
    pub counts: Vec<usize>,
    /// Density of each bin (`count / (n * bin_width)`).
    pub density: Vec<f64>,
}

impl Histogram {
    /// The center of each bin.
    pub fn centers(&self) -> Vec<f64> {
        self.edges.windows(2).map(|w| (w[0] + w[1]) / 2.0).collect()
    }

    /// The uniform bin width.
    pub fn bin_width(&self) -> f64 {
        if self.edges.len() < 2 {
            0.0
        } else {
            self.edges[1] - self.edges[0]
        }
    }
}

/// Build a [`Histogram`] from `data` using the chosen [`BinRule`]. Non-finite
/// values are ignored.
///
/// # Errors
/// * [`StatsError::EmptyInput`] if no finite values remain.
/// * [`StatsError::ZeroVariance`] if all values are identical (no range to bin).
/// * [`StatsError::InvalidBandwidth`] if an explicit [`BinRule::Width`] or
///   [`BinRule::Count`] is non-positive.
pub fn histogram(data: &[f64], rule: BinRule) -> Result<Histogram, StatsError> {
    let sorted = sorted_finite(data);
    if sorted.is_empty() {
        return Err(StatsError::EmptyInput);
    }
    let n = sorted.len();
    let min = sorted[0];
    let max = sorted[n - 1];
    if max <= min {
        return Err(StatsError::ZeroVariance);
    }
    let range = max - min;
    let nf = n as f64;

    let bins = match rule {
        BinRule::Count(c) => {
            if c == 0 {
                return Err(StatsError::InvalidBandwidth);
            }
            c
        }
        BinRule::Width(w) => {
            if !(w.is_finite() && w > 0.0) {
                return Err(StatsError::InvalidBandwidth);
            }
            (range / w).ceil().max(1.0) as usize
        }
        BinRule::Sturges => (nf.log2().ceil() as usize) + 1,
        BinRule::FreedmanDiaconis => {
            let q = quartiles(&sorted)?;
            let w = 2.0 * q.iqr * nf.powf(-1.0 / 3.0);
            if w > 0.0 {
                (range / w).ceil().max(1.0) as usize
            } else {
                (nf.log2().ceil() as usize) + 1 // fall back to Sturges on zero IQR
            }
        }
        BinRule::Scott => {
            let mean = sorted.iter().sum::<f64>() / nf;
            let std = (sorted.iter().map(|x| (x - mean).powi(2)).sum::<f64>()
                / (nf - 1.0).max(1.0))
            .sqrt();
            let w = 3.49 * std * nf.powf(-1.0 / 3.0);
            if w > 0.0 {
                (range / w).ceil().max(1.0) as usize
            } else {
                (nf.log2().ceil() as usize) + 1
            }
        }
    }
    .max(1);

    let width = range / bins as f64;
    let edges: Vec<f64> = (0..=bins).map(|i| min + width * i as f64).collect();
    let mut counts = vec![0usize; bins];
    for &v in &sorted {
        let mut idx = ((v - min) / width).floor() as usize;
        if idx >= bins {
            idx = bins - 1; // the maximum lands in the last bin
        }
        counts[idx] += 1;
    }
    let density = counts.iter().map(|&c| c as f64 / (nf * width)).collect();

    Ok(Histogram {
        edges,
        counts,
        density,
    })
}