mean.hpp 4.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147
  1. // Copyright 2015-2018 Hans Dembinski
  2. //
  3. // Distributed under the Boost Software License, version 1.0.
  4. // (See accompanying file LICENSE_1_0.txt
  5. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  6. #ifndef BOOST_HISTOGRAM_ACCUMULATORS_MEAN_HPP
  7. #define BOOST_HISTOGRAM_ACCUMULATORS_MEAN_HPP
  8. #include <boost/assert.hpp>
  9. #include <boost/core/nvp.hpp>
  10. #include <boost/histogram/fwd.hpp> // for mean<>
  11. #include <boost/throw_exception.hpp>
  12. #include <stdexcept>
  13. #include <type_traits>
  14. namespace boost {
  15. namespace histogram {
  16. namespace accumulators {
  17. /** Calculates mean and variance of sample.
  18. Uses Welfords's incremental algorithm to improve the numerical
  19. stability of mean and variance computation.
  20. */
  21. template <class ValueType>
  22. class mean {
  23. public:
  24. using value_type = ValueType;
  25. using const_reference = const value_type&;
  26. mean() = default;
  27. /// Allow implicit conversion from mean<T>
  28. template <class T>
  29. mean(const mean<T>& o) noexcept
  30. : sum_{o.sum_}, mean_{o.mean_}, sum_of_deltas_squared_{o.sum_of_deltas_squared_} {}
  31. /// Initialize to external count, mean, and variance
  32. mean(const_reference n, const_reference mean, const_reference variance) noexcept
  33. : sum_(n), mean_(mean), sum_of_deltas_squared_(variance * (n - 1)) {}
  34. /// Insert sample x
  35. void operator()(const_reference x) noexcept {
  36. sum_ += static_cast<value_type>(1);
  37. const auto delta = x - mean_;
  38. mean_ += delta / sum_;
  39. sum_of_deltas_squared_ += delta * (x - mean_);
  40. }
  41. /// Insert sample x with weight w
  42. void operator()(const weight_type<value_type>& w, const_reference x) noexcept {
  43. sum_ += w.value;
  44. const auto delta = x - mean_;
  45. mean_ += w.value * delta / sum_;
  46. sum_of_deltas_squared_ += w.value * delta * (x - mean_);
  47. }
  48. /// Add another mean accumulator
  49. mean& operator+=(const mean& rhs) noexcept {
  50. if (sum_ != 0 || rhs.sum_ != 0) {
  51. const auto tmp = mean_ * sum_ + rhs.mean_ * rhs.sum_;
  52. sum_ += rhs.sum_;
  53. mean_ = tmp / sum_;
  54. }
  55. sum_of_deltas_squared_ += rhs.sum_of_deltas_squared_;
  56. return *this;
  57. }
  58. /** Scale by value
  59. This acts as if all samples were scaled by the value.
  60. */
  61. mean& operator*=(const_reference s) noexcept {
  62. mean_ *= s;
  63. sum_of_deltas_squared_ *= s * s;
  64. return *this;
  65. }
  66. bool operator==(const mean& rhs) const noexcept {
  67. return sum_ == rhs.sum_ && mean_ == rhs.mean_ &&
  68. sum_of_deltas_squared_ == rhs.sum_of_deltas_squared_;
  69. }
  70. bool operator!=(const mean& rhs) const noexcept { return !operator==(rhs); }
  71. /// Return how many samples were accumulated
  72. const_reference count() const noexcept { return sum_; }
  73. /// Return mean value of accumulated samples
  74. const_reference value() const noexcept { return mean_; }
  75. /// Return variance of accumulated samples
  76. value_type variance() const noexcept { return sum_of_deltas_squared_ / (sum_ - 1); }
  77. template <class Archive>
  78. void serialize(Archive& ar, unsigned version) {
  79. if (version == 0) {
  80. // read only
  81. std::size_t sum;
  82. ar& make_nvp("sum", sum);
  83. sum_ = static_cast<value_type>(sum);
  84. } else {
  85. ar& make_nvp("sum", sum_);
  86. }
  87. ar& make_nvp("mean", mean_);
  88. ar& make_nvp("sum_of_deltas_squared", sum_of_deltas_squared_);
  89. }
  90. private:
  91. value_type sum_{};
  92. value_type mean_{};
  93. value_type sum_of_deltas_squared_{};
  94. };
  95. } // namespace accumulators
  96. } // namespace histogram
  97. } // namespace boost
  98. #ifndef BOOST_HISTOGRAM_DOXYGEN_INVOKED
  99. namespace boost {
  100. namespace serialization {
  101. template <class T>
  102. struct version;
  103. // version 1 for boost::histogram::accumulators::mean<T>
  104. template <class T>
  105. struct version<boost::histogram::accumulators::mean<T>> : std::integral_constant<int, 1> {
  106. };
  107. } // namespace serialization
  108. } // namespace boost
  109. namespace std {
  110. template <class T, class U>
  111. /// Specialization for boost::histogram::accumulators::mean.
  112. struct common_type<boost::histogram::accumulators::mean<T>,
  113. boost::histogram::accumulators::mean<U>> {
  114. using type = boost::histogram::accumulators::mean<common_type_t<T, U>>;
  115. };
  116. } // namespace std
  117. #endif
  118. #endif