src/libtest/stats.rs

   1 // Copyright 2012 The Rust Project Developers. See the COPYRIGHT
   2 // file at the top-level directory of this distribution and at
   3 // http://rust-lang.org/COPYRIGHT.
   4 //
   5 // Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
   6 // http://www.apache.org/licenses/LICENSE-2.0> or the MIT license
   7 // <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your
   8 // option. This file may not be copied, modified, or distributed
   9 // except according to those terms.
  10
  11 #![allow(missing_docs)]
  12
  13 use std::cmp::Ordering::{mod, Less, Greater, Equal};
  14 use std::collections::hash_map::Entry::{Occupied, Vacant};
  15 use std::collections::hash_map;
  16 use std::fmt::Show;
  17 use std::hash::Hash;
  18 use std::io;
  19 use std::mem;
  20 use std::num::{Float, FloatMath, FromPrimitive};
  21
  22 fn local_cmp<T:Float>(x: T, y: T) -> Ordering {
  23     // arbitrarily decide that NaNs are larger than everything.
  24     if y.is_nan() {
  25         Less
  26     } else if x.is_nan() {
  27         Greater
  28     } else if x < y {
  29         Less
  30     } else if x == y {
  31         Equal
  32     } else {
  33         Greater
  34     }
  35 }
  36
  37 fn local_sort<T: Float>(v: &mut [T]) {
  38     v.sort_by(|x: &T, y: &T| local_cmp(*x, *y));
  39 }
  40
  41 /// Trait that provides simple descriptive statistics on a univariate set of numeric samples.
  42 pub trait Stats <T: FloatMath + FromPrimitive> for Sized? {
  43
  44     /// Sum of the samples.
  45     ///
  46     /// Note: this method sacrifices performance at the altar of accuracy
  47     /// Depends on IEEE-754 arithmetic guarantees. See proof of correctness at:
  48     /// ["Adaptive Precision Floating-Point Arithmetic and Fast Robust Geometric Predicates"]
  49     /// (http://www.cs.cmu.edu/~quake-papers/robust-arithmetic.ps)
  50     /// *Discrete & Computational Geometry 18*, 3 (Oct 1997), 305-363, Shewchuk J.R.
  51     fn sum(&self) -> T;
  52
  53     /// Minimum value of the samples.
  54     fn min(&self) -> T;
  55
  56     /// Maximum value of the samples.
  57     fn max(&self) -> T;
  58
  59     /// Arithmetic mean (average) of the samples: sum divided by sample-count.
  60     ///
  61     /// See: https://en.wikipedia.org/wiki/Arithmetic_mean
  62     fn mean(&self) -> T;
  63
  64     /// Median of the samples: value separating the lower half of the samples from the higher half.
  65     /// Equal to `self.percentile(50.0)`.
  66     ///
  67     /// See: https://en.wikipedia.org/wiki/Median
  68     fn median(&self) -> T;
  69
  70     /// Variance of the samples: bias-corrected mean of the squares of the differences of each
  71     /// sample from the sample mean. Note that this calculates the _sample variance_ rather than the
  72     /// population variance, which is assumed to be unknown. It therefore corrects the `(n-1)/n`
  73     /// bias that would appear if we calculated a population variance, by dividing by `(n-1)` rather
  74     /// than `n`.
  75     ///
  76     /// See: https://en.wikipedia.org/wiki/Variance
  77     fn var(&self) -> T;
  78
  79     /// Standard deviation: the square root of the sample variance.
  80     ///
  81     /// Note: this is not a robust statistic for non-normal distributions. Prefer the
  82     /// `median_abs_dev` for unknown distributions.
  83     ///
  84     /// See: https://en.wikipedia.org/wiki/Standard_deviation
  85     fn std_dev(&self) -> T;
  86
  87     /// Standard deviation as a percent of the mean value. See `std_dev` and `mean`.
  88     ///
  89     /// Note: this is not a robust statistic for non-normal distributions. Prefer the
  90     /// `median_abs_dev_pct` for unknown distributions.
  91     fn std_dev_pct(&self) -> T;
  92
  93     /// Scaled median of the absolute deviations of each sample from the sample median. This is a
  94     /// robust (distribution-agnostic) estimator of sample variability. Use this in preference to
  95     /// `std_dev` if you cannot assume your sample is normally distributed. Note that this is scaled
  96     /// by the constant `1.4826` to allow its use as a consistent estimator for the standard
  97     /// deviation.
  98     ///
  99     /// See: http://en.wikipedia.org/wiki/Median_absolute_deviation
 100     fn median_abs_dev(&self) -> T;
 101
 102     /// Median absolute deviation as a percent of the median. See `median_abs_dev` and `median`.
 103     fn median_abs_dev_pct(&self) -> T;
 104
 105     /// Percentile: the value below which `pct` percent of the values in `self` fall. For example,
 106     /// percentile(95.0) will return the value `v` such that 95% of the samples `s` in `self`
 107     /// satisfy `s <= v`.
 108     ///
 109     /// Calculated by linear interpolation between closest ranks.
 110     ///
 111     /// See: http://en.wikipedia.org/wiki/Percentile
 112     fn percentile(&self, pct: T) -> T;
 113
 114     /// Quartiles of the sample: three values that divide the sample into four equal groups, each
 115     /// with 1/4 of the data. The middle value is the median. See `median` and `percentile`. This
 116     /// function may calculate the 3 quartiles more efficiently than 3 calls to `percentile`, but
 117     /// is otherwise equivalent.
 118     ///
 119     /// See also: https://en.wikipedia.org/wiki/Quartile
 120     fn quartiles(&self) -> (T,T,T);
 121
 122     /// Inter-quartile range: the difference between the 25th percentile (1st quartile) and the 75th
 123     /// percentile (3rd quartile). See `quartiles`.
 124     ///
 125     /// See also: https://en.wikipedia.org/wiki/Interquartile_range
 126     fn iqr(&self) -> T;
 127 }
 128
 129 /// Extracted collection of all the summary statistics of a sample set.
 130 #[deriving(Clone, PartialEq)]
 131 #[allow(missing_docs)]
 132 pub struct Summary<T> {
 133     pub sum: T,
 134     pub min: T,
 135     pub max: T,
 136     pub mean: T,
 137     pub median: T,
 138     pub var: T,
 139     pub std_dev: T,
 140     pub std_dev_pct: T,
 141     pub median_abs_dev: T,
 142     pub median_abs_dev_pct: T,
 143     pub quartiles: (T,T,T),
 144     pub iqr: T,
 145 }
 146
 147 impl<T: FloatMath + FromPrimitive> Summary<T> {
 148     /// Construct a new summary of a sample set.
 149     pub fn new(samples: &[T]) -> Summary<T> {
 150         Summary {
 151             sum: samples.sum(),
 152             min: samples.min(),
 153             max: samples.max(),
 154             mean: samples.mean(),
 155             median: samples.median(),
 156             var: samples.var(),
 157             std_dev: samples.std_dev(),
 158             std_dev_pct: samples.std_dev_pct(),
 159             median_abs_dev: samples.median_abs_dev(),
 160             median_abs_dev_pct: samples.median_abs_dev_pct(),
 161             quartiles: samples.quartiles(),
 162             iqr: samples.iqr()
 163         }
 164     }
 165 }
 166
 167 impl<T: FloatMath + FromPrimitive> Stats<T> for [T] {
 168     // FIXME #11059 handle NaN, inf and overflow
 169     fn sum(&self) -> T {
 170         let mut partials = vec![];
 171
 172         for &mut x in self.iter() {
 173             let mut j = 0;
 174             // This inner loop applies `hi`/`lo` summation to each
 175             // partial so that the list of partial sums remains exact.
 176             for i in range(0, partials.len()) {
 177                 let mut y: T = partials[i];
 178                 if x.abs() < y.abs() {
 179                     mem::swap(&mut x, &mut y);
 180                 }
 181                 // Rounded `x+y` is stored in `hi` with round-off stored in
 182                 // `lo`. Together `hi+lo` are exactly equal to `x+y`.
 183                 let hi = x + y;
 184                 let lo = y - (hi - x);
 185                 if lo != Float::zero() {
 186                     partials[j] = lo;
 187                     j += 1;
 188                 }
 189                 x = hi;
 190             }
 191             if j >= partials.len() {
 192                 partials.push(x);
 193             } else {
 194                 partials[j] = x;
 195                 partials.truncate(j+1);
 196             }
 197         }
 198         let zero: T = Float::zero();
 199         partials.iter().fold(zero, |p, q| p + *q)
 200     }
 201
 202     fn min(&self) -> T {
 203         assert!(self.len() != 0);
 204         self.iter().fold(self[0], |p, q| p.min(*q))
 205     }
 206
 207     fn max(&self) -> T {
 208         assert!(self.len() != 0);
 209         self.iter().fold(self[0], |p, q| p.max(*q))
 210     }
 211
 212     fn mean(&self) -> T {
 213         assert!(self.len() != 0);
 214         self.sum() / FromPrimitive::from_uint(self.len()).unwrap()
 215     }
 216
 217     fn median(&self) -> T {
 218         self.percentile(FromPrimitive::from_uint(50).unwrap())
 219     }
 220
 221     fn var(&self) -> T {
 222         if self.len() < 2 {
 223             Float::zero()
 224         } else {
 225             let mean = self.mean();
 226             let mut v: T = Float::zero();
 227             for s in self.iter() {
 228                 let x = *s - mean;
 229                 v = v + x*x;
 230             }
 231             // NB: this is _supposed to be_ len-1, not len. If you
 232             // change it back to len, you will be calculating a
 233             // population variance, not a sample variance.
 234             let denom = FromPrimitive::from_uint(self.len()-1).unwrap();
 235             v/denom
 236         }
 237     }
 238
 239     fn std_dev(&self) -> T {
 240         self.var().sqrt()
 241     }
 242
 243     fn std_dev_pct(&self) -> T {
 244         let hundred = FromPrimitive::from_uint(100).unwrap();
 245         (self.std_dev() / self.mean()) * hundred
 246     }
 247
 248     fn median_abs_dev(&self) -> T {
 249         let med = self.median();
 250         let abs_devs: Vec<T> = self.iter().map(|&v| (med - v).abs()).collect();
 251         // This constant is derived by smarter statistics brains than me, but it is
 252         // consistent with how R and other packages treat the MAD.
 253         let number = FromPrimitive::from_f64(1.4826).unwrap();
 254         abs_devs.median() * number
 255     }
 256
 257     fn median_abs_dev_pct(&self) -> T {
 258         let hundred = FromPrimitive::from_uint(100).unwrap();
 259         (self.median_abs_dev() / self.median()) * hundred
 260     }
 261
 262     fn percentile(&self, pct: T) -> T {
 263         let mut tmp = self.to_vec();
 264         local_sort(tmp.as_mut_slice());
 265         percentile_of_sorted(tmp.as_slice(), pct)
 266     }
 267
 268     fn quartiles(&self) -> (T,T,T) {
 269         let mut tmp = self.to_vec();
 270         local_sort(tmp.as_mut_slice());
 271         let first = FromPrimitive::from_uint(25).unwrap();
 272         let a = percentile_of_sorted(tmp.as_slice(), first);
 273         let secound = FromPrimitive::from_uint(50).unwrap();
 274         let b = percentile_of_sorted(tmp.as_slice(), secound);
 275         let third = FromPrimitive::from_uint(75).unwrap();
 276         let c = percentile_of_sorted(tmp.as_slice(), third);
 277         (a,b,c)
 278     }
 279
 280     fn iqr(&self) -> T {
 281         let (a,_,c) = self.quartiles();
 282         c - a
 283     }
 284 }
 285
 286
 287 // Helper function: extract a value representing the `pct` percentile of a sorted sample-set, using
 288 // linear interpolation. If samples are not sorted, return nonsensical value.
 289 fn percentile_of_sorted<T: Float + FromPrimitive>(sorted_samples: &[T],
 290                                                              pct: T) -> T {
 291     assert!(sorted_samples.len() != 0);
 292     if sorted_samples.len() == 1 {
 293         return sorted_samples[0];
 294     }
 295     let zero: T = Float::zero();
 296     assert!(zero <= pct);
 297     let hundred = FromPrimitive::from_uint(100).unwrap();
 298     assert!(pct <= hundred);
 299     if pct == hundred {
 300         return sorted_samples[sorted_samples.len() - 1];
 301     }
 302     let length = FromPrimitive::from_uint(sorted_samples.len() - 1).unwrap();
 303     let rank = (pct / hundred) * length;
 304     let lrank = rank.floor();
 305     let d = rank - lrank;
 306     let n = lrank.to_uint().unwrap();
 307     let lo = sorted_samples[n];
 308     let hi = sorted_samples[n+1];
 309     lo + (hi - lo) * d
 310 }
 311
 312
 313 /// Winsorize a set of samples, replacing values above the `100-pct` percentile and below the `pct`
 314 /// percentile with those percentiles themselves. This is a way of minimizing the effect of
 315 /// outliers, at the cost of biasing the sample. It differs from trimming in that it does not
 316 /// change the number of samples, just changes the values of those that are outliers.
 317 ///
 318 /// See: http://en.wikipedia.org/wiki/Winsorising
 319 pub fn winsorize<T: Float + FromPrimitive>(samples: &mut [T], pct: T) {
 320     let mut tmp = samples.to_vec();
 321     local_sort(tmp.as_mut_slice());
 322     let lo = percentile_of_sorted(tmp.as_slice(), pct);
 323     let hundred: T = FromPrimitive::from_uint(100).unwrap();
 324     let hi = percentile_of_sorted(tmp.as_slice(), hundred-pct);
 325     for samp in samples.iter_mut() {
 326         if *samp > hi {
 327             *samp = hi
 328         } else if *samp < lo {
 329             *samp = lo
 330         }
 331     }
 332 }
 333
 334 /// Render writes the min, max and quartiles of the provided `Summary` to the provided `Writer`.
 335 pub fn write_5_number_summary<W: Writer, T: Float + Show>(w: &mut W,
 336                                                           s: &Summary<T>) -> io::IoResult<()> {
 337     let (q1,q2,q3) = s.quartiles;
 338     write!(w, "(min={}, q1={}, med={}, q3={}, max={})",
 339                      s.min,
 340                      q1,
 341                      q2,
 342                      q3,
 343                      s.max)
 344 }
 345
 346 /// Render a boxplot to the provided writer. The boxplot shows the min, max and quartiles of the
 347 /// provided `Summary` (thus includes the mean) and is scaled to display within the range of the
 348 /// nearest multiple-of-a-power-of-ten above and below the min and max of possible values, and
 349 /// target `width_hint` characters of display (though it will be wider if necessary).
 350 ///
 351 /// As an example, the summary with 5-number-summary `(min=15, q1=17, med=20, q3=24, max=31)` might
 352 /// display as:
 353 ///
 354 /// ```{.ignore}
 355 ///   10 |        [--****#******----------]          | 40
 356 /// ```
 357 pub fn write_boxplot<W: Writer, T: Float + Show + FromPrimitive>(
 358                      w: &mut W,
 359                      s: &Summary<T>,
 360                      width_hint: uint)
 361                       -> io::IoResult<()> {
 362
 363     let (q1,q2,q3) = s.quartiles;
 364
 365     // the .abs() handles the case where numbers are negative
 366     let ten: T = FromPrimitive::from_uint(10).unwrap();
 367     let lomag = ten.powf(s.min.abs().log10().floor());
 368     let himag = ten.powf(s.max.abs().log10().floor());
 369
 370     // need to consider when the limit is zero
 371     let zero: T = Float::zero();
 372     let lo = if lomag == Float::zero() {
 373         zero
 374     } else {
 375         (s.min / lomag).floor() * lomag
 376     };
 377
 378     let hi = if himag == Float::zero() {
 379         zero
 380     } else {
 381         (s.max / himag).ceil() * himag
 382     };
 383
 384     let range = hi - lo;
 385
 386     let lostr = lo.to_string();
 387     let histr = hi.to_string();
 388
 389     let overhead_width = lostr.len() + histr.len() + 4;
 390     let range_width = width_hint - overhead_width;
 391     let range_float = FromPrimitive::from_uint(range_width).unwrap();
 392     let char_step = range / range_float;
 393
 394     try!(write!(w, "{} |", lostr));
 395
 396     let mut c = 0;
 397     let mut v = lo;
 398
 399     while c < range_width && v < s.min {
 400         try!(write!(w, " "));
 401         v = v + char_step;
 402         c += 1;
 403     }
 404     try!(write!(w, "["));
 405     c += 1;
 406     while c < range_width && v < q1 {
 407         try!(write!(w, "-"));
 408         v = v + char_step;
 409         c += 1;
 410     }
 411     while c < range_width && v < q2 {
 412         try!(write!(w, "*"));
 413         v = v + char_step;
 414         c += 1;
 415     }
 416     try!(write!(w, "#"));
 417     c += 1;
 418     while c < range_width && v < q3 {
 419         try!(write!(w, "*"));
 420         v = v + char_step;
 421         c += 1;
 422     }
 423     while c < range_width && v < s.max {
 424         try!(write!(w, "-"));
 425         v = v + char_step;
 426         c += 1;
 427     }
 428     try!(write!(w, "]"));
 429     while c < range_width {
 430         try!(write!(w, " "));
 431         v = v + char_step;
 432         c += 1;
 433     }
 434
 435     try!(write!(w, "| {}", histr));
 436     Ok(())
 437 }
 438
 439 /// Returns a HashMap with the number of occurrences of every element in the
 440 /// sequence that the iterator exposes.
 441 pub fn freq_count<T: Iterator<Item=U>, U: Eq+Hash>(mut iter: T) -> hash_map::HashMap<U, uint> {
 442     let mut map: hash_map::HashMap<U,uint> = hash_map::HashMap::new();
 443     for elem in iter {
 444         match map.entry(elem) {
 445             Occupied(mut entry) => { *entry.get_mut() += 1; },
 446             Vacant(entry) => { entry.set(1); },
 447         }
 448     }
 449     map
 450 }
 451
 452 // Test vectors generated from R, using the script src/etc/stat-test-vectors.r.
 453
 454 #[cfg(test)]
 455 mod tests {
 456     use stats::Stats;
 457     use stats::Summary;
 458     use stats::write_5_number_summary;
 459     use stats::write_boxplot;
 460     use std::io;
 461     use std::f64;
 462
 463     macro_rules! assert_approx_eq {
 464         ($a:expr, $b:expr) => ({
 465             use std::num::Float;
 466             let (a, b) = (&$a, &$b);
 467             assert!((*a - *b).abs() < 1.0e-6,
 468                     "{} is not approximately equal to {}", *a, *b);
 469         })
 470     }
 471
 472     fn check(samples: &[f64], summ: &Summary<f64>) {
 473
 474         let summ2 = Summary::new(samples);
 475
 476         let mut w = io::stdout();
 477         let w = &mut w;
 478         (write!(w, "\n")).unwrap();
 479         write_5_number_summary(w, &summ2).unwrap();
 480         (write!(w, "\n")).unwrap();
 481         write_boxplot(w, &summ2, 50).unwrap();
 482         (write!(w, "\n")).unwrap();
 483
 484         assert_eq!(summ.sum, summ2.sum);
 485         assert_eq!(summ.min, summ2.min);
 486         assert_eq!(summ.max, summ2.max);
 487         assert_eq!(summ.mean, summ2.mean);
 488         assert_eq!(summ.median, summ2.median);
 489
 490         // We needed a few more digits to get exact equality on these
 491         // but they're within float epsilon, which is 1.0e-6.
 492         assert_approx_eq!(summ.var, summ2.var);
 493         assert_approx_eq!(summ.std_dev, summ2.std_dev);
 494         assert_approx_eq!(summ.std_dev_pct, summ2.std_dev_pct);
 495         assert_approx_eq!(summ.median_abs_dev, summ2.median_abs_dev);
 496         assert_approx_eq!(summ.median_abs_dev_pct, summ2.median_abs_dev_pct);
 497
 498         assert_eq!(summ.quartiles, summ2.quartiles);
 499         assert_eq!(summ.iqr, summ2.iqr);
 500     }
 501
 502     #[test]
 503     fn test_min_max_nan() {
 504         let xs = &[1.0, 2.0, f64::NAN, 3.0, 4.0];
 505         let summary = Summary::new(xs);
 506         assert_eq!(summary.min, 1.0);
 507         assert_eq!(summary.max, 4.0);
 508     }
 509
 510     #[test]
 511     fn test_norm2() {
 512         let val = &[
 513             958.0000000000,
 514             924.0000000000,
 515         ];
 516         let summ = &Summary {
 517             sum: 1882.0000000000,
 518             min: 924.0000000000,
 519             max: 958.0000000000,
 520             mean: 941.0000000000,
 521             median: 941.0000000000,
 522             var: 578.0000000000,
 523             std_dev: 24.0416305603,
 524             std_dev_pct: 2.5549022912,
 525             median_abs_dev: 25.2042000000,
 526             median_abs_dev_pct: 2.6784484591,
 527             quartiles: (932.5000000000,941.0000000000,949.5000000000),
 528             iqr: 17.0000000000,
 529         };
 530         check(val, summ);
 531     }
 532     #[test]
 533     fn test_norm10narrow() {
 534         let val = &[
 535             966.0000000000,
 536             985.0000000000,
 537             1110.0000000000,
 538             848.0000000000,
 539             821.0000000000,
 540             975.0000000000,
 541             962.0000000000,
 542             1157.0000000000,
 543             1217.0000000000,
 544             955.0000000000,
 545         ];
 546         let summ = &Summary {
 547             sum: 9996.0000000000,
 548             min: 821.0000000000,
 549             max: 1217.0000000000,
 550             mean: 999.6000000000,
 551             median: 970.5000000000,
 552             var: 16050.7111111111,
 553             std_dev: 126.6914010938,
 554             std_dev_pct: 12.6742097933,
 555             median_abs_dev: 102.2994000000,
 556             median_abs_dev_pct: 10.5408964451,
 557             quartiles: (956.7500000000,970.5000000000,1078.7500000000),
 558             iqr: 122.0000000000,
 559         };
 560         check(val, summ);
 561     }
 562     #[test]
 563     fn test_norm10medium() {
 564         let val = &[
 565             954.0000000000,
 566             1064.0000000000,
 567             855.0000000000,
 568             1000.0000000000,
 569             743.0000000000,
 570             1084.0000000000,
 571             704.0000000000,
 572             1023.0000000000,
 573             357.0000000000,
 574             869.0000000000,
 575         ];
 576         let summ = &Summary {
 577             sum: 8653.0000000000,
 578             min: 357.0000000000,
 579             max: 1084.0000000000,
 580             mean: 865.3000000000,
 581             median: 911.5000000000,
 582             var: 48628.4555555556,
 583             std_dev: 220.5186059170,
 584             std_dev_pct: 25.4846418487,
 585             median_abs_dev: 195.7032000000,
 586             median_abs_dev_pct: 21.4704552935,
 587             quartiles: (771.0000000000,911.5000000000,1017.2500000000),
 588             iqr: 246.2500000000,
 589         };
 590         check(val, summ);
 591     }
 592     #[test]
 593     fn test_norm10wide() {
 594         let val = &[
 595             505.0000000000,
 596             497.0000000000,
 597             1591.0000000000,
 598             887.0000000000,
 599             1026.0000000000,
 600             136.0000000000,
 601             1580.0000000000,
 602             940.0000000000,
 603             754.0000000000,
 604             1433.0000000000,
 605         ];
 606         let summ = &Summary {
 607             sum: 9349.0000000000,
 608             min: 136.0000000000,
 609             max: 1591.0000000000,
 610             mean: 934.9000000000,
 611             median: 913.5000000000,
 612             var: 239208.9888888889,
 613             std_dev: 489.0899599142,
 614             std_dev_pct: 52.3146817750,
 615             median_abs_dev: 611.5725000000,
 616             median_abs_dev_pct: 66.9482758621,
 617             quartiles: (567.2500000000,913.5000000000,1331.2500000000),
 618             iqr: 764.0000000000,
 619         };
 620         check(val, summ);
 621     }
 622     #[test]
 623     fn test_norm25verynarrow() {
 624         let val = &[
 625             991.0000000000,
 626             1018.0000000000,
 627             998.0000000000,
 628             1013.0000000000,
 629             974.0000000000,
 630             1007.0000000000,
 631             1014.0000000000,
 632             999.0000000000,
 633             1011.0000000000,
 634             978.0000000000,
 635             985.0000000000,
 636             999.0000000000,
 637             983.0000000000,
 638             982.0000000000,
 639             1015.0000000000,
 640             1002.0000000000,
 641             977.0000000000,
 642             948.0000000000,
 643             1040.0000000000,
 644             974.0000000000,
 645             996.0000000000,
 646             989.0000000000,
 647             1015.0000000000,
 648             994.0000000000,
 649             1024.0000000000,
 650         ];
 651         let summ = &Summary {
 652             sum: 24926.0000000000,
 653             min: 948.0000000000,
 654             max: 1040.0000000000,
 655             mean: 997.0400000000,
 656             median: 998.0000000000,
 657             var: 393.2066666667,
 658             std_dev: 19.8294393937,
 659             std_dev_pct: 1.9888308788,
 660             median_abs_dev: 22.2390000000,
 661             median_abs_dev_pct: 2.2283567134,
 662             quartiles: (983.0000000000,998.0000000000,1013.0000000000),
 663             iqr: 30.0000000000,
 664         };
 665         check(val, summ);
 666     }
 667     #[test]
 668     fn test_exp10a() {
 669         let val = &[
 670             23.0000000000,
 671             11.0000000000,
 672             2.0000000000,
 673             57.0000000000,
 674             4.0000000000,
 675             12.0000000000,
 676             5.0000000000,
 677             29.0000000000,
 678             3.0000000000,
 679             21.0000000000,
 680         ];
 681         let summ = &Summary {
 682             sum: 167.0000000000,
 683             min: 2.0000000000,
 684             max: 57.0000000000,
 685             mean: 16.7000000000,
 686             median: 11.5000000000,
 687             var: 287.7888888889,
 688             std_dev: 16.9643416875,
 689             std_dev_pct: 101.5828843560,
 690             median_abs_dev: 13.3434000000,
 691             median_abs_dev_pct: 116.0295652174,
 692             quartiles: (4.2500000000,11.5000000000,22.5000000000),
 693             iqr: 18.2500000000,
 694         };
 695         check(val, summ);
 696     }
 697     #[test]
 698     fn test_exp10b() {
 699         let val = &[
 700             24.0000000000,
 701             17.0000000000,
 702             6.0000000000,
 703             38.0000000000,
 704             25.0000000000,
 705             7.0000000000,
 706             51.0000000000,
 707             2.0000000000,
 708             61.0000000000,
 709             32.0000000000,
 710         ];
 711         let summ = &Summary {
 712             sum: 263.0000000000,
 713             min: 2.0000000000,
 714             max: 61.0000000000,
 715             mean: 26.3000000000,
 716             median: 24.5000000000,
 717             var: 383.5666666667,
 718             std_dev: 19.5848580967,
 719             std_dev_pct: 74.4671410520,
 720             median_abs_dev: 22.9803000000,
 721             median_abs_dev_pct: 93.7971428571,
 722             quartiles: (9.5000000000,24.5000000000,36.5000000000),
 723             iqr: 27.0000000000,
 724         };
 725         check(val, summ);
 726     }
 727     #[test]
 728     fn test_exp10c() {
 729         let val = &[
 730             71.0000000000,
 731             2.0000000000,
 732             32.0000000000,
 733             1.0000000000,
 734             6.0000000000,
 735             28.0000000000,
 736             13.0000000000,
 737             37.0000000000,
 738             16.0000000000,
 739             36.0000000000,
 740         ];
 741         let summ = &Summary {
 742             sum: 242.0000000000,
 743             min: 1.0000000000,
 744             max: 71.0000000000,
 745             mean: 24.2000000000,
 746             median: 22.0000000000,
 747             var: 458.1777777778,
 748             std_dev: 21.4050876611,
 749             std_dev_pct: 88.4507754589,
 750             median_abs_dev: 21.4977000000,
 751             median_abs_dev_pct: 97.7168181818,
 752             quartiles: (7.7500000000,22.0000000000,35.0000000000),
 753             iqr: 27.2500000000,
 754         };
 755         check(val, summ);
 756     }
 757     #[test]
 758     fn test_exp25() {
 759         let val = &[
 760             3.0000000000,
 761             24.0000000000,
 762             1.0000000000,
 763             19.0000000000,
 764             7.0000000000,
 765             5.0000000000,
 766             30.0000000000,
 767             39.0000000000,
 768             31.0000000000,
 769             13.0000000000,
 770             25.0000000000,
 771             48.0000000000,
 772             1.0000000000,
 773             6.0000000000,
 774             42.0000000000,
 775             63.0000000000,
 776             2.0000000000,
 777             12.0000000000,
 778             108.0000000000,
 779             26.0000000000,
 780             1.0000000000,
 781             7.0000000000,
 782             44.0000000000,
 783             25.0000000000,
 784             11.0000000000,
 785         ];
 786         let summ = &Summary {
 787             sum: 593.0000000000,
 788             min: 1.0000000000,
 789             max: 108.0000000000,
 790             mean: 23.7200000000,
 791             median: 19.0000000000,
 792             var: 601.0433333333,
 793             std_dev: 24.5161851301,
 794             std_dev_pct: 103.3565983562,
 795             median_abs_dev: 19.2738000000,
 796             median_abs_dev_pct: 101.4410526316,
 797             quartiles: (6.0000000000,19.0000000000,31.0000000000),
 798             iqr: 25.0000000000,
 799         };
 800         check(val, summ);
 801     }
 802     #[test]
 803     fn test_binom25() {
 804         let val = &[
 805             18.0000000000,
 806             17.0000000000,
 807             27.0000000000,
 808             15.0000000000,
 809             21.0000000000,
 810             25.0000000000,
 811             17.0000000000,
 812             24.0000000000,
 813             25.0000000000,
 814             24.0000000000,
 815             26.0000000000,
 816             26.0000000000,
 817             23.0000000000,
 818             15.0000000000,
 819             23.0000000000,
 820             17.0000000000,
 821             18.0000000000,
 822             18.0000000000,
 823             21.0000000000,
 824             16.0000000000,
 825             15.0000000000,
 826             31.0000000000,
 827             20.0000000000,
 828             17.0000000000,
 829             15.0000000000,
 830         ];
 831         let summ = &Summary {
 832             sum: 514.0000000000,
 833             min: 15.0000000000,
 834             max: 31.0000000000,
 835             mean: 20.5600000000,
 836             median: 20.0000000000,
 837             var: 20.8400000000,
 838             std_dev: 4.5650848842,
 839             std_dev_pct: 22.2037202539,
 840             median_abs_dev: 5.9304000000,
 841             median_abs_dev_pct: 29.6520000000,
 842             quartiles: (17.0000000000,20.0000000000,24.0000000000),
 843             iqr: 7.0000000000,
 844         };
 845         check(val, summ);
 846     }
 847     #[test]
 848     fn test_pois25lambda30() {
 849         let val = &[
 850             27.0000000000,
 851             33.0000000000,
 852             34.0000000000,
 853             34.0000000000,
 854             24.0000000000,
 855             39.0000000000,
 856             28.0000000000,
 857             27.0000000000,
 858             31.0000000000,
 859             28.0000000000,
 860             38.0000000000,
 861             21.0000000000,
 862             33.0000000000,
 863             36.0000000000,
 864             29.0000000000,
 865             37.0000000000,
 866             32.0000000000,
 867             34.0000000000,
 868             31.0000000000,
 869             39.0000000000,
 870             25.0000000000,
 871             31.0000000000,
 872             32.0000000000,
 873             40.0000000000,
 874             24.0000000000,
 875         ];
 876         let summ = &Summary {
 877             sum: 787.0000000000,
 878             min: 21.0000000000,
 879             max: 40.0000000000,
 880             mean: 31.4800000000,
 881             median: 32.0000000000,
 882             var: 26.5933333333,
 883             std_dev: 5.1568724372,
 884             std_dev_pct: 16.3814245145,
 885             median_abs_dev: 5.9304000000,
 886             median_abs_dev_pct: 18.5325000000,
 887             quartiles: (28.0000000000,32.0000000000,34.0000000000),
 888             iqr: 6.0000000000,
 889         };
 890         check(val, summ);
 891     }
 892     #[test]
 893     fn test_pois25lambda40() {
 894         let val = &[
 895             42.0000000000,
 896             50.0000000000,
 897             42.0000000000,
 898             46.0000000000,
 899             34.0000000000,
 900             45.0000000000,
 901             34.0000000000,
 902             49.0000000000,
 903             39.0000000000,
 904             28.0000000000,
 905             40.0000000000,
 906             35.0000000000,
 907             37.0000000000,
 908             39.0000000000,
 909             46.0000000000,
 910             44.0000000000,
 911             32.0000000000,
 912             45.0000000000,
 913             42.0000000000,
 914             37.0000000000,
 915             48.0000000000,
 916             42.0000000000,
 917             33.0000000000,
 918             42.0000000000,
 919             48.0000000000,
 920         ];
 921         let summ = &Summary {
 922             sum: 1019.0000000000,
 923             min: 28.0000000000,
 924             max: 50.0000000000,
 925             mean: 40.7600000000,
 926             median: 42.0000000000,
 927             var: 34.4400000000,
 928             std_dev: 5.8685603004,
 929             std_dev_pct: 14.3978417577,
 930             median_abs_dev: 5.9304000000,
 931             median_abs_dev_pct: 14.1200000000,
 932             quartiles: (37.0000000000,42.0000000000,45.0000000000),
 933             iqr: 8.0000000000,
 934         };
 935         check(val, summ);
 936     }
 937     #[test]
 938     fn test_pois25lambda50() {
 939         let val = &[
 940             45.0000000000,
 941             43.0000000000,
 942             44.0000000000,
 943             61.0000000000,
 944             51.0000000000,
 945             53.0000000000,
 946             59.0000000000,
 947             52.0000000000,
 948             49.0000000000,
 949             51.0000000000,
 950             51.0000000000,
 951             50.0000000000,
 952             49.0000000000,
 953             56.0000000000,
 954             42.0000000000,
 955             52.0000000000,
 956             51.0000000000,
 957             43.0000000000,
 958             48.0000000000,
 959             48.0000000000,
 960             50.0000000000,
 961             42.0000000000,
 962             43.0000000000,
 963             42.0000000000,
 964             60.0000000000,
 965         ];
 966         let summ = &Summary {
 967             sum: 1235.0000000000,
 968             min: 42.0000000000,
 969             max: 61.0000000000,
 970             mean: 49.4000000000,
 971             median: 50.0000000000,
 972             var: 31.6666666667,
 973             std_dev: 5.6273143387,
 974             std_dev_pct: 11.3913245723,
 975             median_abs_dev: 4.4478000000,
 976             median_abs_dev_pct: 8.8956000000,
 977             quartiles: (44.0000000000,50.0000000000,52.0000000000),
 978             iqr: 8.0000000000,
 979         };
 980         check(val, summ);
 981     }
 982     #[test]
 983     fn test_unif25() {
 984         let val = &[
 985             99.0000000000,
 986             55.0000000000,
 987             92.0000000000,
 988             79.0000000000,
 989             14.0000000000,
 990             2.0000000000,
 991             33.0000000000,
 992             49.0000000000,
 993             3.0000000000,
 994             32.0000000000,
 995             84.0000000000,
 996             59.0000000000,
 997             22.0000000000,
 998             86.0000000000,
 999             76.0000000000,
1000             31.0000000000,
1001             29.0000000000,
1002             11.0000000000,
1003             41.0000000000,
1004             53.0000000000,
1005             45.0000000000,
1006             44.0000000000,
1007             98.0000000000,
1008             98.0000000000,
1009             7.0000000000,
1010         ];
1011         let summ = &Summary {
1012             sum: 1242.0000000000,
1013             min: 2.0000000000,
1014             max: 99.0000000000,
1015             mean: 49.6800000000,
1016             median: 45.0000000000,
1017             var: 1015.6433333333,
1018             std_dev: 31.8691595957,
1019             std_dev_pct: 64.1488719719,
1020             median_abs_dev: 45.9606000000,
1021             median_abs_dev_pct: 102.1346666667,
1022             quartiles: (29.0000000000,45.0000000000,79.0000000000),
1023             iqr: 50.0000000000,
1024         };
1025         check(val, summ);
1026     }
1027
1028     #[test]
1029     fn test_boxplot_nonpositive() {
1030         fn t(s: &Summary<f64>, expected: String) {
1031             let mut m = Vec::new();
1032             write_boxplot(&mut m, s, 30).unwrap();
1033             let out = String::from_utf8(m).unwrap();
1034             assert_eq!(out, expected);
1035         }
1036
1037         t(&Summary::new(&[-2.0f64, -1.0f64]),
1038                         "-2 |[------******#*****---]| -1".to_string());
1039         t(&Summary::new(&[0.0f64, 2.0f64]),
1040                         "0 |[-------*****#*******---]| 2".to_string());
1041         t(&Summary::new(&[-2.0f64, 0.0f64]),
1042                         "-2 |[------******#******---]| 0".to_string());
1043
1044     }
1045     #[test]
1046     fn test_sum_f64s() {
1047         assert_eq!([0.5f64, 3.2321f64, 1.5678f64].sum(), 5.2999);
1048     }
1049     #[test]
1050     fn test_sum_f64_between_ints_that_sum_to_0() {
1051         assert_eq!([1e30f64, 1.2f64, -1e30f64].sum(), 1.2);
1052     }
1053 }
1054
1055 #[cfg(test)]
1056 mod bench {
1057     use Bencher;
1058     use stats::Stats;
1059
1060     #[bench]
1061     pub fn sum_three_items(b: &mut Bencher) {
1062         b.iter(|| {
1063             [1e20f64, 1.5f64, -1e20f64].sum();
1064         })
1065     }
1066     #[bench]
1067     pub fn sum_many_f64(b: &mut Bencher) {
1068         let nums = [-1e30f64, 1e60, 1e30, 1.0, -1e60];
1069         let v = Vec::from_fn(500, |i| nums[i%5]);
1070
1071         b.iter(|| {
1072             v.as_slice().sum();
1073         })
1074     }
1075 }