A toolkit for working with phylogenetic data.
v0.20.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
statistics.hpp
Go to the documentation of this file.
1 #ifndef GENESIS_UTILS_MATH_STATISTICS_H_
2 #define GENESIS_UTILS_MATH_STATISTICS_H_
3 
4 /*
5  Genesis - A toolkit for working with phylogenetic data.
6  Copyright (C) 2014-2018 Lucas Czech and HITS gGmbH
7 
8  This program is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  This program is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with this program. If not, see <http://www.gnu.org/licenses/>.
20 
21  Contact:
22  Lucas Czech <lucas.czech@h-its.org>
23  Exelixis Lab, Heidelberg Institute for Theoretical Studies
24  Schloss-Wolfsbrunnenweg 35, D-69118 Heidelberg, Germany
25 */
26 
35 
36 #include <algorithm>
37 #include <cassert>
38 #include <cmath>
39 #include <cstddef>
40 #include <limits>
41 #include <stdexcept>
42 #include <utility>
43 #include <vector>
44 
45 namespace genesis {
46 namespace utils {
47 
48 // =================================================================================================
49 // Forward Declarations
50 // =================================================================================================
51 
52 // Needed for spearmans_rank_correlation_coefficient()
53 std::vector<double> ranking_fractional( std::vector<double> const& vec );
54 
55 // =================================================================================================
56 // Structures and Classes
57 // =================================================================================================
58 
64 template< typename T >
65 struct MinMaxPair
66 {
67  T min;
68  T max;
69 };
70 
78 {
79  double mean;
80  double stddev;
81 };
82 
87 struct Quartiles
88 {
89  double q0 = 0.0;
90  double q1 = 0.0;
91  double q2 = 0.0;
92  double q3 = 0.0;
93  double q4 = 0.0;
94 };
95 
96 // =================================================================================================
97 // Mean Stddev
98 // =================================================================================================
99 
113 template <class ForwardIterator>
114 MeanStddevPair mean_stddev( ForwardIterator first, ForwardIterator last, double epsilon = -1.0 )
115 {
116  // Prepare result.
117  MeanStddevPair result;
118  result.mean = 0.0;
119  result.stddev = 0.0;
120  size_t count = 0;
121 
122  // Sum up elements.
123  auto it = first;
124  while( it != last ) {
125  if( std::isfinite( *it ) ) {
126  result.mean += *it;
127  ++count;
128  }
129  ++it;
130  }
131 
132  // If there are no valid elements, return an all-zero result.
133  if( count == 0 ) {
134  return result;
135  }
136 
137  // Calculate mean.
138  result.mean /= static_cast<double>( count );
139 
140  // Calculate column std dev.
141  it = first;
142  while( it != last ) {
143  if( std::isfinite( *it ) ) {
144  result.stddev += (( *it - result.mean ) * ( *it - result.mean ));
145  }
146  ++it;
147  }
148  assert( count > 0 );
149  result.stddev /= static_cast<double>( count );
150  result.stddev = std::sqrt( result.stddev );
151 
152  // The following in an inelegant (but usual) way to handle near-zero values,
153  // which later would cause a division by zero.
154  assert( result.stddev >= 0.0 );
155  if( result.stddev <= epsilon ){
156  result.stddev = 1.0;
157  }
158 
159  return result;
160 }
161 
167 inline MeanStddevPair mean_stddev( std::vector<double> const& vec, double epsilon = -1.0 )
168 {
169  return mean_stddev( vec.begin(), vec.end(), epsilon );
170 }
171 
172 // =================================================================================================
173 // Median
174 // =================================================================================================
175 
176 // TODO this weird range plus inidces implementation comes from before the whole functionw as a template.
177 // now, using just the range should be enough, so no need for the helper function!
178 
182 template <class RandomAccessIterator>
183 double median( RandomAccessIterator first, RandomAccessIterator last, size_t l, size_t r )
184 {
185  auto const size = static_cast<size_t>( std::distance( first, last ));
186  assert( l < size && r < size && l <= r );
187  (void) size;
188 
189  // Size of the interval.
190  size_t const sz = r - l + 1;
191 
192  // Even or odd size? Median is calculated differently.
193  if( sz % 2 == 0 ) {
194 
195  // Get the two middle positions.
196  size_t pl = l + sz / 2 - 1;
197  size_t pu = l + sz / 2;
198  assert( pl < size && pu < size );
199 
200  return ( *(first + pl) + *(first + pu) ) / 2.0;
201 
202  } else {
203 
204  // Int division, rounds down. This is what we want.
205  size_t p = l + sz / 2;
206  assert( p < size );
207 
208  return *(first + p);
209  }
210 }
211 
217 template <class RandomAccessIterator>
218 double median( RandomAccessIterator first, RandomAccessIterator last )
219 {
220  // Checks.
221  if( ! std::is_sorted( first, last )) {
222  throw std::runtime_error( "Range has to be sorted for median calculation." );
223  }
224  auto const size = static_cast<size_t>( std::distance( first, last ));
225  if( size == 0 ) {
226  return 0.0;
227  }
228 
229  // Use helper function, which takes the range inclusively.
230  return median( first, last, 0, size - 1 );
231 }
232 
238 inline double median( std::vector<double> const& vec )
239 {
240  return median( vec.begin(), vec.end() );
241 }
242 
243 // =================================================================================================
244 // Quartiles
245 // =================================================================================================
246 
247 template <class RandomAccessIterator>
248 Quartiles quartiles( RandomAccessIterator first, RandomAccessIterator last )
249 {
250  // Prepare result.
251  Quartiles result;
252 
253  // Checks.
254  if( ! std::is_sorted( first, last )) {
255  throw std::runtime_error( "Range has to be sorted for quartiles calculation." );
256  }
257  auto const size = static_cast<size_t>( std::distance( first, last ));
258  if( size == 0 ) {
259  return result;
260  }
261 
262  // Set min, 50% and max.
263  result.q0 = *first;
264  result.q2 = median( first, last, 0, size - 1 );
265  result.q4 = *(first + size - 1);
266 
267  // Even or odd size? Quartiles are calculated differently.
268  // This could be done shorter, but this way feels more expressive.
269  if( size % 2 == 0 ) {
270 
271  // Even: Split exaclty in halves.
272  result.q1 = median( first, last, 0, size / 2 - 1 );
273  result.q3 = median( first, last, size / 2, size - 1 );
274 
275  } else {
276 
277  // Odd: Do not include the median value itself.
278  result.q1 = median( first, last, 0, size / 2 - 1 );
279  result.q3 = median( first, last, size / 2 + 1, size - 1 );
280  }
281 
282  return result;
283 }
284 
290 inline Quartiles quartiles( std::vector<double> const& vec )
291 {
292  return quartiles( vec.begin(), vec.end() );
293 }
294 
295 // =================================================================================================
296 // Dispersion
297 // =================================================================================================
298 
307 inline double coefficient_of_variation( MeanStddevPair const& ms )
308 {
309  return ms.stddev / ms.mean;
310 }
311 
315 inline std::vector<double> coefficient_of_variation( std::vector<MeanStddevPair> const& ms )
316 {
317  auto res = std::vector<double>( ms.size() );
318  for( size_t i = 0; i < ms.size(); ++i ) {
319  res[ i ] = coefficient_of_variation( ms[i] );
320  }
321  return res;
322 }
323 
333 inline double index_of_dispersion( MeanStddevPair const& ms )
334 {
335  return ms.stddev * ms.stddev / ms.mean;
336 }
337 
341 inline std::vector<double> index_of_dispersion( std::vector<MeanStddevPair> const& ms )
342 {
343  auto res = std::vector<double>( ms.size() );
344  for( size_t i = 0; i < ms.size(); ++i ) {
345  res[ i ] = index_of_dispersion( ms[i] );
346  }
347  return res;
348 }
349 
358 {
359  return ( q.q3 - q.q1 ) / ( q.q3 + q.q1 );
360 }
361 
365 inline std::vector<double> quartile_coefficient_of_dispersion( std::vector<Quartiles> const& q )
366 {
367  auto res = std::vector<double>( q.size() );
368  for( size_t i = 0; i < q.size(); ++i ) {
369  res[ i ] = quartile_coefficient_of_dispersion( q[i] );
370  }
371  return res;
372 }
373 
374 // =================================================================================================
375 // Correlation Coefficients
376 // =================================================================================================
377 
385 template <class ForwardIteratorA, class ForwardIteratorB>
386 std::pair<std::vector<double>, std::vector<double>> finite_pairs(
387  ForwardIteratorA first_a, ForwardIteratorA last_a,
388  ForwardIteratorB first_b, ForwardIteratorB last_b
389 ) {
390  // Prepare result.
391  std::vector<double> vec_a;
392  std::vector<double> vec_b;
393 
394  // Iterate in parallel.
395  while( first_a != last_a && first_b != last_b ) {
396  if( std::isfinite( *first_a ) && std::isfinite( *first_b ) ) {
397  vec_a.push_back( *first_a );
398  vec_b.push_back( *first_b );
399  }
400  ++first_a;
401  ++first_b;
402  }
403  if( first_a != last_a || first_b != last_b ) {
404  throw std::runtime_error(
405  "Ranges need to have same length."
406  );
407  }
408 
409  assert( vec_a.size() == vec_b.size() );
410  return { vec_a, vec_b };
411 }
412 
425 template <class ForwardIteratorA, class ForwardIteratorB>
427  ForwardIteratorA first_a, ForwardIteratorA last_a,
428  ForwardIteratorB first_b, ForwardIteratorB last_b
429 ) {
430  // Calculate means.
431  double mean_a = 0.0;
432  double mean_b = 0.0;
433  size_t count = 0;
434  auto it_a = first_a;
435  auto it_b = first_b;
436  while( it_a != last_a && it_b != last_b ) {
437  if( std::isfinite( *it_a ) && std::isfinite( *it_b ) ) {
438  mean_a += *it_a;
439  mean_b += *it_b;
440  ++count;
441  }
442  ++it_a;
443  ++it_b;
444  }
445  if( it_a != last_a || it_b != last_b ) {
446  throw std::runtime_error(
447  "Ranges need to have same length to calculate their Pearson Correlation Coefficient."
448  );
449  }
450  if( count == 0 ) {
451  return std::numeric_limits<double>::quiet_NaN();
452  }
453  assert( count > 0 );
454  mean_a /= static_cast<double>( count );
455  mean_b /= static_cast<double>( count );
456 
457  // Calculate PCC parts.
458  double numerator = 0.0;
459  double stddev_a = 0.0;
460  double stddev_b = 0.0;
461  it_a = first_a;
462  it_b = first_b;
463  while( it_a != last_a && it_b != last_b ) {
464  if( std::isfinite( *it_a ) && std::isfinite( *it_b ) ) {
465  double const d1 = *it_a - mean_a;
466  double const d2 = *it_b - mean_b;
467  numerator += d1 * d2;
468  stddev_a += d1 * d1;
469  stddev_b += d2 * d2;
470  }
471  ++it_a;
472  ++it_b;
473  }
474  assert( it_a == last_a && it_b == last_b );
475 
476  // Calcualte PCC, and assert that it is in the correct range
477  // (or not a number, which can happen if the std dev is 0.0, e.g. in all-zero vectors).
478  auto const pcc = numerator / ( std::sqrt( stddev_a ) * std::sqrt( stddev_b ) );
479  assert(( -1.0 <= pcc && pcc <= 1.0 ) || ( ! std::isfinite( pcc ) ));
480  return pcc;
481 }
482 
489  std::vector<double> const& vec_a,
490  std::vector<double> const& vec_b
491 ) {
493  vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end()
494  );
495 }
496 
505 template <class RandomAccessIteratorA, class RandomAccessIteratorB>
507  RandomAccessIteratorA first_a, RandomAccessIteratorA last_a,
508  RandomAccessIteratorB first_b, RandomAccessIteratorB last_b
509 ) {
510  // Get cleaned results.
511  auto const cleaned = finite_pairs( first_a, last_a, first_b, last_b );
512 
513  // Get the ranking of both vectors.
514  auto ranks_a = ranking_fractional( cleaned.first );
515  auto ranks_b = ranking_fractional( cleaned.second );
516  assert( ranks_a.size() == ranks_b.size() );
517 
518  // Nice try. But removing them later does not work, as the ranges would be calculated
519  // differently... So we have to live with making copies of the data :-(
520  //
521  // if( ranks_a.size() != ranks_b.size() ) {
522  // throw std::runtime_error(
523  // "Ranges need to have same length to calculate their "
524  // "Spearman's Rank Correlation Coefficient."
525  // );
526  // }
527  //
528  // // Remove non finite values.
529  // size_t ins = 0;
530  // size_t pos = 0;
531  // while( first_a != last_a && first_b != last_b ) {
532  // if( std::isfinite( *first_a ) && std::isfinite( *first_b ) ) {
533  // ranks_a[ ins ] = ranks_a[ pos ];
534  // ranks_b[ ins ] = ranks_b[ pos ];
535  // ++ins;
536  // }
537  // ++first_a;
538  // ++first_b;
539  // ++pos;
540  // }
541  //
542  // // We already checked that the ranks have the same length. Assert this.
543  // assert( first_a == last_a && first_b == last_b );
544  // assert( pos == ranks_a.size() && pos == ranks_b.size() );
545  //
546  // // Erase the removed elements.
547  // assert( ins <= ranks_a.size() && ins <= ranks_b.size() );
548  // ranks_a.resize( ins );
549  // ranks_b.resize( ins );
550 
551  return pearson_correlation_coefficient( ranks_a, ranks_b );
552 }
553 
560  std::vector<double> const& vec_a,
561  std::vector<double> const& vec_b
562 ) {
564  vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end()
565  );
566 }
567 
579 inline double fisher_transformation( double correlation_coefficient )
580 {
581  auto const r = correlation_coefficient;
582  if( r < -1.0 || r > 1.0 ) {
583  throw std::invalid_argument(
584  "Cannot apply fisher transformation to value " + std::to_string( r ) +
585  " outside of [ -1.0, 1.0 ]."
586  );
587  }
588 
589  // LOG_DBG << "formula " << 0.5 * log( ( 1.0 + r ) / ( 1.0 - r ) );
590  // LOG_DBG << "simple " << std::atanh( r );
591  return std::atanh( r );
592 }
593 
599 inline std::vector<double> fisher_transformation( std::vector<double> const& correlation_coefficients )
600 {
601  auto res = correlation_coefficients;
602  for( auto& elem : res ) {
603  elem = fisher_transformation( elem );
604  }
605  return res;
606 }
607 
608 // =================================================================================================
609 // Ranking Standard
610 // =================================================================================================
611 
621 template <class RandomAccessIterator>
622 std::vector<size_t> ranking_standard( RandomAccessIterator first, RandomAccessIterator last )
623 {
624  // Prepare result, and get the sorting order of the vector.
625  auto const size = static_cast<size_t>( std::distance( first, last ));
626  auto result = std::vector<size_t>( size, 1 );
627  auto const order = stable_sort_indices( first, last );
628 
629  // Shortcuts for better readability.
630  auto ordered_value = [&]( size_t i ){
631  return *( first + order[i] );
632  };
633  auto ordered_result = [&]( size_t i ) -> size_t& {
634  return result[ order[i] ];
635  };
636 
637  // Calculate ranks.
638  for( size_t i = 1; i < size; ++i ) {
639 
640  // Same values get the same rank. The next bigger one continues at the current i.
641  if( ordered_value( i ) == ordered_value( i - 1 ) ) {
642  ordered_result( i ) = ordered_result( i - 1 );
643  } else {
644  ordered_result( i ) = i + 1;
645  }
646  }
647 
648  return result;
649 }
650 
654 inline std::vector<size_t> ranking_standard( std::vector<double> const& vec )
655 {
656  return ranking_standard( vec.begin(), vec.end() );
657 }
658 
659 // =================================================================================================
660 // Ranking Modified
661 // =================================================================================================
662 
672 template <class RandomAccessIterator>
673 std::vector<size_t> ranking_modified( RandomAccessIterator first, RandomAccessIterator last )
674 {
675  // Prepare result, and get the sorting order of the vector.
676  auto const size = static_cast<size_t>( std::distance( first, last ));
677  auto result = std::vector<size_t>( size, 1 );
678  auto const order = stable_sort_indices( first, last );
679 
680  // Shortcuts for better readability.
681  auto ordered_value = [&]( size_t i ){
682  return *( first + order[i] );
683  };
684  auto ordered_result = [&]( size_t i ) -> size_t& {
685  return result[ order[i] ];
686  };
687 
688  // Calculate ranks. The loop variable is incremented at the end.
689  for( size_t i = 0; i < size; ) {
690 
691  // Look ahead: How often does the value occur?
692  size_t j = 1;
693  while( i+j < size && ordered_value(i+j) == ordered_value(i) ) {
694  ++j;
695  }
696 
697  // Set the j-next entries.
698  for( size_t k = 0; k < j; ++k ) {
699  ordered_result( i + k ) = i + j;
700  }
701 
702  // We can skip the j-next loop iterations, as we just set their values
703  i += j;
704  }
705 
706  return result;
707 }
708 
712 inline std::vector<size_t> ranking_modified( std::vector<double> const& vec )
713 {
714  return ranking_modified( vec.begin(), vec.end() );
715 }
716 
717 // =================================================================================================
718 // Ranking Dense
719 // =================================================================================================
720 
729 template <class RandomAccessIterator>
730 std::vector<size_t> ranking_dense( RandomAccessIterator first, RandomAccessIterator last )
731 {
732  // Prepare result, and get the sorting order of the vector.
733  auto const size = static_cast<size_t>( std::distance( first, last ));
734  auto result = std::vector<size_t>( size, 1 );
735  auto const order = stable_sort_indices( first, last );
736 
737  // Shortcuts for better readability.
738  auto ordered_value = [&]( size_t i ){
739  return *( first + order[i] );
740  };
741  auto ordered_result = [&]( size_t i ) -> size_t& {
742  return result[ order[i] ];
743  };
744 
745  // Calculate ranks.
746  for( size_t i = 1; i < size; ++i ) {
747 
748  // Same values get the same rank. The next bigger one continues by incrementing.
749  if( ordered_value( i ) == ordered_value( i - 1 ) ) {
750  ordered_result( i ) = ordered_result( i - 1 );
751  } else {
752  ordered_result( i ) = ordered_result( i - 1 ) + 1;
753  }
754  }
755 
756  return result;
757 }
758 
762 inline std::vector<size_t> ranking_dense( std::vector<double> const& vec )
763 {
764  return ranking_dense( vec.begin(), vec.end() );
765 }
766 
767 // =================================================================================================
768 // Ranking Ordinal
769 // =================================================================================================
770 
779 template <class RandomAccessIterator>
780 std::vector<size_t> ranking_ordinal( RandomAccessIterator first, RandomAccessIterator last )
781 {
782  // Prepare result, and get the sorting order of the vector.
783  auto const size = static_cast<size_t>( std::distance( first, last ));
784  auto result = std::vector<size_t>( size, 1 );
785  auto const order = stable_sort_indices( first, last );
786 
787  // Shortcuts for better readability.
788  auto ordered_result = [&]( size_t i ) -> size_t& {
789  return result[ order[i] ];
790  };
791 
792  // Calculate ranks. This is simply the order plus 1 (as ranks are 1-based).
793  for( size_t i = 0; i < size; ++i ) {
794  ordered_result( i ) = i + 1;
795  }
796 
797  return result;
798 }
799 
803 inline std::vector<size_t> ranking_ordinal( std::vector<double> const& vec )
804 {
805  return ranking_ordinal( vec.begin(), vec.end() );
806 }
807 
808 // =================================================================================================
809 // Ranking Fractional
810 // =================================================================================================
811 
822 template <class RandomAccessIterator>
823 std::vector<double> ranking_fractional( RandomAccessIterator first, RandomAccessIterator last )
824 {
825  // Prepare result, and get the sorting order of the vector.
826  auto const size = static_cast<size_t>( std::distance( first, last ));
827  auto result = std::vector<double>( size, 1 );
828  auto const order = stable_sort_indices( first, last );
829 
830  // Shortcuts for better readability.
831  auto ordered_value = [&]( size_t i ){
832  return *( first + order[i] );
833  };
834  auto ordered_result = [&]( size_t i ) -> double& {
835  return result[ order[i] ];
836  };
837 
838  // Calculate the average of the sum of numbers in the given inclusive range.
839  auto sum_avg = []( size_t l, size_t r )
840  {
841  assert( l <= r );
842 
843  // Example: l == 7, r == 9
844  // We want: (7 + 8 + 9) / 3 = 8.0
845  // Upper: 1+2+3+4+5+6+7+8+9 = 45
846  // Lower: 1+2+3+4+5+6 = 21
847  // Diff: 45 - 21 = 24
848  // Count: 9 - 7 + 1 = 3
849  // Result: 24 / 3 = 8
850  auto const upper = r * ( r + 1 ) / 2;
851  auto const lower = ( l - 1 ) * l / 2;
852  return static_cast<double>( upper - lower ) / static_cast<double>( r - l + 1 );
853  };
854 
855  // Calculate ranks. The loop variable is incremented at the end.
856  for( size_t i = 0; i < size; ) {
857 
858  // Look ahead: How often does the value occur?
859  size_t j = 1;
860  while( i+j < size && ordered_value(i+j) == ordered_value(i) ) {
861  ++j;
862  }
863 
864  // Set the j-next entries.
865  auto entry = sum_avg( i + 1, i + j );
866  for( size_t k = 0; k < j; ++k ) {
867  ordered_result( i + k ) = entry;
868  }
869 
870  // We can skip the j-next loop iterations, as we just set their values
871  i += j;
872  }
873 
874  return result;
875 }
876 
880 inline std::vector<double> ranking_fractional( std::vector<double> const& vec )
881 {
882  return ranking_fractional( vec.begin(), vec.end() );
883 }
884 
885 } // namespace utils
886 } // namespace genesis
887 
888 #endif // include guard
std::vector< size_t > ranking_modified(RandomAccessIterator first, RandomAccessIterator last)
Return the ranking of the values in the given range, using Modified competition ranking ("1334" ranki...
Definition: statistics.hpp:673
Provides some valuable algorithms that are not part of the C++ 11 STL.
Store a pair of min and max values.
Definition: statistics.hpp:65
double pearson_correlation_coefficient(ForwardIteratorA first_a, ForwardIteratorA last_a, ForwardIteratorB first_b, ForwardIteratorB last_b)
Calculate the Pearson Correlation Coefficient between two ranges of double.
Definition: statistics.hpp:426
std::vector< size_t > ranking_ordinal(RandomAccessIterator first, RandomAccessIterator last)
Return the ranking of the values in the given range, using Ordinal ranking ("1234" ranking)...
Definition: statistics.hpp:780
double fisher_transformation(double correlation_coefficient)
Apply Fisher z-transformation to a correlation coefficient.
Definition: statistics.hpp:579
std::pair< std::vector< double >, std::vector< double > > finite_pairs(ForwardIteratorA first_a, ForwardIteratorA last_a, ForwardIteratorB first_b, ForwardIteratorB last_b)
Helper function that cleans two ranges of double of the same length from non-finite values...
Definition: statistics.hpp:386
double median(const Histogram &h)
Quartiles quartiles(RandomAccessIterator first, RandomAccessIterator last)
Definition: statistics.hpp:248
double quartile_coefficient_of_dispersion(Quartiles const &q)
Calculate the quartile_coefficient_of_dispersion.
Definition: statistics.hpp:357
std::string to_string(T const &v)
Return a string representation of a given value.
Definition: string.hpp:381
std::vector< size_t > ranking_dense(RandomAccessIterator first, RandomAccessIterator last)
Return the ranking of the values in the given range, using Dense ranking ("1223" ranking).
Definition: statistics.hpp:730
double coefficient_of_variation(MeanStddevPair const &ms)
Calculate the index of dispersion.
Definition: statistics.hpp:307
Store a mean and a standard deviation value.
Definition: statistics.hpp:77
std::vector< size_t > ranking_standard(RandomAccessIterator first, RandomAccessIterator last)
Return the ranking of the values in the given range, using Standard competition ranking ("1224" ranki...
Definition: statistics.hpp:622
double index_of_dispersion(MeanStddevPair const &ms)
Calculate the index of dispersion.
Definition: statistics.hpp:333
std::vector< size_t > stable_sort_indices(RandomAccessIterator first, RandomAccessIterator last, Comparator comparator)
Get the indices to the stable sorted order of the given range.
Definition: algorithm.hpp:181
std::vector< double > ranking_fractional(std::vector< double > const &vec)
Return the ranking of the values in the given range, using Fractional ranking ("1 2...
Definition: statistics.hpp:880
Store the values of quartiles: q0 == min, q1 == 25%, q2 == 50%, q3 == 75%, q4 == max.
Definition: statistics.hpp:87
double spearmans_rank_correlation_coefficient(RandomAccessIteratorA first_a, RandomAccessIteratorA last_a, RandomAccessIteratorB first_b, RandomAccessIteratorB last_b)
Calculate Spearman's Rank Correlation Coefficient between two ranges of double.
Definition: statistics.hpp:506
MeanStddevPair mean_stddev(ForwardIterator first, ForwardIterator last, double epsilon=-1.0)
Calculate the mean and standard deviation of a range of double elements.
Definition: statistics.hpp:114