// Copyright 2012 The Rust Project Developers. See the COPYRIGHT // file at the top-level directory of this distribution and at // http://rust-lang.org/COPYRIGHT. // // Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or // http://www.apache.org/licenses/LICENSE-2.0> or the MIT license // <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your // option. This file may not be copied, modified, or distributed // except according to those terms.
use core::cmp::Ordering::{self, Equal, Greater, Less}; use core::mem;
fn local_cmp(x: f64, y: f64) -> Ordering { // arbitrarily decide that NaNs are larger than everything. if y.is_nan() {
Less
} elseif x.is_nan() {
Greater
} elseif x < y {
Less
} elseif x == y {
Equal
} else {
Greater
}
}
/// Trait that provides simple descriptive statistics on a univariate set of numeric samples. pubtrait Stats { /// Sum of the samples. /// /// Note: this method sacrifices performance at the altar of accuracy /// Depends on IEEE-754 arithmetic guarantees. See proof of correctness at: /// ["Adaptive Precision Floating-Point Arithmetic and Fast Robust Geometric Predicates"] /// (http://www.cs.cmu.edu/~quake-papers/robust-arithmetic.ps) fn sum(&self) -> f64;
/// Minimum value of the samples. fn min(&self) -> f64;
/// Maximum value of the samples. fn max(&self) -> f64;
/// Median of the samples: value separating the lower half of the samples from the higher half. /// Equal to `self.percentile(50.0)`. /// /// See: https://en.wikipedia.org/wiki/Median fn median(&self) -> f64;
/// Variance of the samples: bias-corrected mean of the squares of the differences of each /// sample from the sample mean. Note that this calculates the _sample variance_ rather than the /// population variance, which is assumed to be unknown. It therefore corrects the `(n-1)/n` /// bias that would appear if we calculated a population variance, by dividing by `(n-1)` rather /// than `n`. /// /// See: https://en.wikipedia.org/wiki/Variance fn var(&self) -> f64;
/// Standard deviation: the square root of the sample variance. /// /// Note: this is not a robust statistic for non-normal distributions. Prefer the /// `median_abs_dev` for unknown distributions. /// /// See: https://en.wikipedia.org/wiki/Standard_deviation fn std_dev(&self) -> f64;
/// Standard deviation as a percent of the mean value. See `std_dev` and `mean`. /// /// Note: this is not a robust statistic for non-normal distributions. Prefer the /// `median_abs_dev_pct` for unknown distributions. fn std_dev_pct(&self) -> f64;
/// Scaled median of the absolute deviations of each sample from the sample median. This is a /// robust (distribution-agnostic) estimator of sample variability. Use this in preference to /// `std_dev` if you cannot assume your sample is normally distributed. Note that this is scaled /// by the constant `1.4826` to allow its use as a consistent estimator for the standard /// deviation. /// /// See: http://en.wikipedia.org/wiki/Median_absolute_deviation fn median_abs_dev(&self) -> f64;
/// Median absolute deviation as a percent of the median. See `median_abs_dev` and `median`. fn median_abs_dev_pct(&self) -> f64;
/// Percentile: the value below which `pct` percent of the values in `self` fall. For example, /// percentile(95.0) will return the value `v` such that 95% of the samples `s` in `self` /// satisfy `s <= v`. /// /// Calculated by linear interpolation between closest ranks. /// /// See: http://en.wikipedia.org/wiki/Percentile fn percentile(&self, pct: f64) -> f64;
/// Quartiles of the sample: three values that divide the sample into four equal groups, each /// with 1/4 of the data. The middle value is the median. See `median` and `percentile`. This /// function may calculate the 3 quartiles more efficiently than 3 calls to `percentile`, but /// is otherwise equivalent. /// /// See also: https://en.wikipedia.org/wiki/Quartile fn quartiles(&self) -> (f64, f64, f64);
/// Inter-quartile range: the difference between the 25th percentile (1st quartile) and the 75th /// percentile (3rd quartile). See `quartiles`. /// /// See also: https://en.wikipedia.org/wiki/Interquartile_range fn iqr(&self) -> f64;
}
for &x inself { letmut x = x; letmut j = 0; // This inner loop applies `hi`/`lo` summation to each // partial so that the list of partial sums remains exact. for i in0..partials.len() { letmut y: f64 = partials[i]; if x.abs() < y.abs() {
mem::swap(&mut x, &mut y);
} // Rounded `x+y` is stored in `hi` with round-off stored in // `lo`. Together `hi+lo` are exactly equal to `x+y`. let hi = x + y; let lo = y - (hi - x); if lo != 0.0 {
partials[j] = lo;
j += 1;
}
x = hi;
} if j >= partials.len() {
partials.push(x);
} else {
partials[j] = x;
partials.truncate(j + 1);
}
} let zero: f64 = 0.0;
partials.iter().fold(zero, |p, q| p + *q)
}
fn var(&self) -> f64 { ifself.len() < 2 { 0.0
} else { let mean = self.mean(); letmut v: f64 = 0.0; for s inself { let x = *s - mean;
v += x * x;
} // NB: this is _supposed to be_ len-1, not len. If you // change it back to len, you will be calculating a // population variance, not a sample variance. let denom = (self.len() - 1) as f64;
v / denom
}
}
fn median_abs_dev(&self) -> f64 { let med = self.median(); let abs_devs: Vec<f64> = self.iter().map(|&v| (med - v).abs()).collect(); // This constant is derived by smarter statistics brains than me, but it is // consistent with how R and other packages treat the MAD. let number = 1.4826;
abs_devs.median() * number
}
fn quartiles(&self) -> (f64, f64, f64) { letmut tmp = self.to_vec();
local_sort(&mut tmp); let first = 25f64; let a = percentile_of_sorted(&tmp, first); let secound = 50f64; let b = percentile_of_sorted(&tmp, secound); let third = 75f64; let c = percentile_of_sorted(&tmp, third);
(a, b, c)
}
fn iqr(&self) -> f64 { let (a, _, c) = self.quartiles();
c - a
}
}
// Helper function: extract a value representing the `pct` percentile of a sorted sample-set, using // linear interpolation. If samples are not sorted, return nonsensical value. fn percentile_of_sorted(sorted_samples: &[f64], pct: f64) -> f64 {
assert!(!sorted_samples.is_empty()); if sorted_samples.len() == 1 { return sorted_samples[0];
} let zero: f64 = 0.0;
assert!(zero <= pct); let hundred = 100f64;
assert!(pct <= hundred); if pct == hundred { return sorted_samples[sorted_samples.len() - 1];
} let length = (sorted_samples.len() - 1) as f64; let rank = (pct / hundred) * length; let lrank = rank.floor(); let d = rank - lrank; let n = lrank as usize; let lo = sorted_samples[n]; let hi = sorted_samples[n + 1];
lo + (hi - lo) * d
}
/// Winsorize a set of samples, replacing values above the `100-pct` percentile /// and below the `pct` percentile with those percentiles themselves. This is a /// way of minimizing the effect of outliers, at the cost of biasing the sample. /// It differs from trimming in that it does not change the number of samples, /// just changes the values of those that are outliers. /// /// See: http://en.wikipedia.org/wiki/Winsorising pubfn winsorize(samples: &mut [f64], pct: f64) { letmut tmp = samples.to_vec();
local_sort(&mut tmp); let lo = percentile_of_sorted(&tmp, pct); let hundred = 100as f64; let hi = percentile_of_sorted(&tmp, hundred - pct); for samp in samples { if *samp > hi {
*samp = hi
} elseif *samp < lo {
*samp = lo
}
}
}
// Test vectors generated from R, using the script src/etc/stat-test-vectors.r.
#[cfg(test)] mod tests { use stats::Stats; use stats::Summary; use std::f64; use std::io; use std::io::prelude::*;
macro_rules! assert_approx_eq {
($a:expr, $b:expr) => {{ let (a, b) = (&$a, &$b);
assert!((*a - *b).abs() < 1.0e-6, "{} is not approximately equal to {}", *a, *b);
}};
}
fn check(samples: &[f64], summ: &Summary) { let summ2 = Summary::new(samples);
letmut w = io::sink(); let w = &mut w;
(write!(w, "\n")).unwrap();
// We needed a few more digits to get exact equality on these // but they're within float epsilon, which is 1.0e-6.
assert_approx_eq!(summ.var, summ2.var);
assert_approx_eq!(summ.std_dev, summ2.std_dev);
assert_approx_eq!(summ.std_dev_pct, summ2.std_dev_pct);
assert_approx_eq!(summ.median_abs_dev, summ2.median_abs_dev);
assert_approx_eq!(summ.median_abs_dev_pct, summ2.median_abs_dev_pct);
Die Informationen auf dieser Webseite wurden
nach bestem Wissen sorgfältig zusammengestellt. Es wird jedoch weder Vollständigkeit, noch Richtigkeit,
noch Qualität der bereit gestellten Informationen zugesichert.
Bemerkung:
Die farbliche Syntaxdarstellung und die Messung sind noch experimentell.