variance.hpp 7.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247
  1. ///////////////////////////////////////////////////////////////////////////////
  2. // variance.hpp
  3. //
  4. // Copyright 2005 Daniel Egloff, Eric Niebler. Distributed under the Boost
  5. // Software License, Version 1.0. (See accompanying file
  6. // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
  7. #ifndef BOOST_ACCUMULATORS_STATISTICS_VARIANCE_HPP_EAN_28_10_2005
  8. #define BOOST_ACCUMULATORS_STATISTICS_VARIANCE_HPP_EAN_28_10_2005
  9. #include <boost/mpl/placeholders.hpp>
  10. #include <boost/accumulators/framework/accumulator_base.hpp>
  11. #include <boost/accumulators/framework/extractor.hpp>
  12. #include <boost/accumulators/numeric/functional.hpp>
  13. #include <boost/accumulators/framework/parameters/sample.hpp>
  14. #include <boost/accumulators/framework/depends_on.hpp>
  15. #include <boost/accumulators/statistics_fwd.hpp>
  16. #include <boost/accumulators/statistics/count.hpp>
  17. #include <boost/accumulators/statistics/sum.hpp>
  18. #include <boost/accumulators/statistics/mean.hpp>
  19. #include <boost/accumulators/statistics/moment.hpp>
  20. namespace boost { namespace accumulators
  21. {
  22. namespace impl
  23. {
  24. //! Lazy calculation of variance.
  25. /*!
  26. Default sample variance implementation based on the second moment \f$ M_n^{(2)} \f$ moment<2>, mean and count.
  27. \f[
  28. \sigma_n^2 = M_n^{(2)} - \mu_n^2.
  29. \f]
  30. where
  31. \f[
  32. \mu_n = \frac{1}{n} \sum_{i = 1}^n x_i.
  33. \f]
  34. is the estimate of the sample mean and \f$n\f$ is the number of samples.
  35. */
  36. template<typename Sample, typename MeanFeature>
  37. struct lazy_variance_impl
  38. : accumulator_base
  39. {
  40. // for boost::result_of
  41. typedef typename numeric::functional::fdiv<Sample, std::size_t>::result_type result_type;
  42. lazy_variance_impl(dont_care) {}
  43. template<typename Args>
  44. result_type result(Args const &args) const
  45. {
  46. extractor<MeanFeature> mean;
  47. result_type tmp = mean(args);
  48. return accumulators::moment<2>(args) - tmp * tmp;
  49. }
  50. // serialization is done by accumulators it depends on
  51. template<class Archive>
  52. void serialize(Archive & ar, const unsigned int file_version) {}
  53. };
  54. //! Iterative calculation of variance.
  55. /*!
  56. Iterative calculation of sample variance \f$\sigma_n^2\f$ according to the formula
  57. \f[
  58. \sigma_n^2 = \frac{1}{n} \sum_{i = 1}^n (x_i - \mu_n)^2 = \frac{n-1}{n} \sigma_{n-1}^2 + \frac{1}{n-1}(x_n - \mu_n)^2.
  59. \f]
  60. where
  61. \f[
  62. \mu_n = \frac{1}{n} \sum_{i = 1}^n x_i.
  63. \f]
  64. is the estimate of the sample mean and \f$n\f$ is the number of samples.
  65. Note that the sample variance is not defined for \f$n <= 1\f$.
  66. A simplification can be obtained by the approximate recursion
  67. \f[
  68. \sigma_n^2 \approx \frac{n-1}{n} \sigma_{n-1}^2 + \frac{1}{n}(x_n - \mu_n)^2.
  69. \f]
  70. because the difference
  71. \f[
  72. \left(\frac{1}{n-1} - \frac{1}{n}\right)(x_n - \mu_n)^2 = \frac{1}{n(n-1)}(x_n - \mu_n)^2.
  73. \f]
  74. converges to zero as \f$n \rightarrow \infty\f$. However, for small \f$ n \f$ the difference
  75. can be non-negligible.
  76. */
  77. template<typename Sample, typename MeanFeature, typename Tag>
  78. struct variance_impl
  79. : accumulator_base
  80. {
  81. // for boost::result_of
  82. typedef typename numeric::functional::fdiv<Sample, std::size_t>::result_type result_type;
  83. template<typename Args>
  84. variance_impl(Args const &args)
  85. : variance(numeric::fdiv(args[sample | Sample()], numeric::one<std::size_t>::value))
  86. {
  87. }
  88. template<typename Args>
  89. void operator ()(Args const &args)
  90. {
  91. std::size_t cnt = count(args);
  92. if(cnt > 1)
  93. {
  94. extractor<MeanFeature> mean;
  95. result_type tmp = args[parameter::keyword<Tag>::get()] - mean(args);
  96. this->variance =
  97. numeric::fdiv(this->variance * (cnt - 1), cnt)
  98. + numeric::fdiv(tmp * tmp, cnt - 1);
  99. }
  100. }
  101. result_type result(dont_care) const
  102. {
  103. return this->variance;
  104. }
  105. // make this accumulator serializeable
  106. template<class Archive>
  107. void serialize(Archive & ar, const unsigned int file_version)
  108. {
  109. ar & variance;
  110. }
  111. private:
  112. result_type variance;
  113. };
  114. } // namespace impl
  115. ///////////////////////////////////////////////////////////////////////////////
  116. // tag::variance
  117. // tag::immediate_variance
  118. //
  119. namespace tag
  120. {
  121. struct lazy_variance
  122. : depends_on<moment<2>, mean>
  123. {
  124. /// INTERNAL ONLY
  125. ///
  126. typedef accumulators::impl::lazy_variance_impl<mpl::_1, mean> impl;
  127. };
  128. struct variance
  129. : depends_on<count, immediate_mean>
  130. {
  131. /// INTERNAL ONLY
  132. ///
  133. typedef accumulators::impl::variance_impl<mpl::_1, mean, sample> impl;
  134. };
  135. }
  136. ///////////////////////////////////////////////////////////////////////////////
  137. // extract::lazy_variance
  138. // extract::variance
  139. //
  140. namespace extract
  141. {
  142. extractor<tag::lazy_variance> const lazy_variance = {};
  143. extractor<tag::variance> const variance = {};
  144. BOOST_ACCUMULATORS_IGNORE_GLOBAL(lazy_variance)
  145. BOOST_ACCUMULATORS_IGNORE_GLOBAL(variance)
  146. }
  147. using extract::lazy_variance;
  148. using extract::variance;
  149. // variance(lazy) -> lazy_variance
  150. template<>
  151. struct as_feature<tag::variance(lazy)>
  152. {
  153. typedef tag::lazy_variance type;
  154. };
  155. // variance(immediate) -> variance
  156. template<>
  157. struct as_feature<tag::variance(immediate)>
  158. {
  159. typedef tag::variance type;
  160. };
  161. // for the purposes of feature-based dependency resolution,
  162. // immediate_variance provides the same feature as variance
  163. template<>
  164. struct feature_of<tag::lazy_variance>
  165. : feature_of<tag::variance>
  166. {
  167. };
  168. // So that variance can be automatically substituted with
  169. // weighted_variance when the weight parameter is non-void.
  170. template<>
  171. struct as_weighted_feature<tag::variance>
  172. {
  173. typedef tag::weighted_variance type;
  174. };
  175. // for the purposes of feature-based dependency resolution,
  176. // weighted_variance provides the same feature as variance
  177. template<>
  178. struct feature_of<tag::weighted_variance>
  179. : feature_of<tag::variance>
  180. {
  181. };
  182. // So that immediate_variance can be automatically substituted with
  183. // immediate_weighted_variance when the weight parameter is non-void.
  184. template<>
  185. struct as_weighted_feature<tag::lazy_variance>
  186. {
  187. typedef tag::lazy_weighted_variance type;
  188. };
  189. // for the purposes of feature-based dependency resolution,
  190. // immediate_weighted_variance provides the same feature as immediate_variance
  191. template<>
  192. struct feature_of<tag::lazy_weighted_variance>
  193. : feature_of<tag::lazy_variance>
  194. {
  195. };
  196. ////////////////////////////////////////////////////////////////////////////
  197. //// droppable_accumulator<variance_impl>
  198. //// need to specialize droppable lazy variance to cache the result at the
  199. //// point the accumulator is dropped.
  200. ///// INTERNAL ONLY
  201. /////
  202. //template<typename Sample, typename MeanFeature>
  203. //struct droppable_accumulator<impl::variance_impl<Sample, MeanFeature> >
  204. // : droppable_accumulator_base<
  205. // with_cached_result<impl::variance_impl<Sample, MeanFeature> >
  206. // >
  207. //{
  208. // template<typename Args>
  209. // droppable_accumulator(Args const &args)
  210. // : droppable_accumulator::base(args)
  211. // {
  212. // }
  213. //};
  214. }} // namespace boost::accumulators
  215. #endif