variance.hpp 7.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236
  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. };
  51. //! Iterative calculation of variance.
  52. /*!
  53. Iterative calculation of sample variance \f$\sigma_n^2\f$ according to the formula
  54. \f[
  55. \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.
  56. \f]
  57. where
  58. \f[
  59. \mu_n = \frac{1}{n} \sum_{i = 1}^n x_i.
  60. \f]
  61. is the estimate of the sample mean and \f$n\f$ is the number of samples.
  62. Note that the sample variance is not defined for \f$n <= 1\f$.
  63. A simplification can be obtained by the approximate recursion
  64. \f[
  65. \sigma_n^2 \approx \frac{n-1}{n} \sigma_{n-1}^2 + \frac{1}{n}(x_n - \mu_n)^2.
  66. \f]
  67. because the difference
  68. \f[
  69. \left(\frac{1}{n-1} - \frac{1}{n}\right)(x_n - \mu_n)^2 = \frac{1}{n(n-1)}(x_n - \mu_n)^2.
  70. \f]
  71. converges to zero as \f$n \rightarrow \infty\f$. However, for small \f$ n \f$ the difference
  72. can be non-negligible.
  73. */
  74. template<typename Sample, typename MeanFeature, typename Tag>
  75. struct variance_impl
  76. : accumulator_base
  77. {
  78. // for boost::result_of
  79. typedef typename numeric::functional::fdiv<Sample, std::size_t>::result_type result_type;
  80. template<typename Args>
  81. variance_impl(Args const &args)
  82. : variance(numeric::fdiv(args[sample | Sample()], numeric::one<std::size_t>::value))
  83. {
  84. }
  85. template<typename Args>
  86. void operator ()(Args const &args)
  87. {
  88. std::size_t cnt = count(args);
  89. if(cnt > 1)
  90. {
  91. extractor<MeanFeature> mean;
  92. result_type tmp = args[parameter::keyword<Tag>::get()] - mean(args);
  93. this->variance =
  94. numeric::fdiv(this->variance * (cnt - 1), cnt)
  95. + numeric::fdiv(tmp * tmp, cnt - 1);
  96. }
  97. }
  98. result_type result(dont_care) const
  99. {
  100. return this->variance;
  101. }
  102. private:
  103. result_type variance;
  104. };
  105. } // namespace impl
  106. ///////////////////////////////////////////////////////////////////////////////
  107. // tag::variance
  108. // tag::immediate_variance
  109. //
  110. namespace tag
  111. {
  112. struct lazy_variance
  113. : depends_on<moment<2>, mean>
  114. {
  115. /// INTERNAL ONLY
  116. ///
  117. typedef accumulators::impl::lazy_variance_impl<mpl::_1, mean> impl;
  118. };
  119. struct variance
  120. : depends_on<count, immediate_mean>
  121. {
  122. /// INTERNAL ONLY
  123. ///
  124. typedef accumulators::impl::variance_impl<mpl::_1, mean, sample> impl;
  125. };
  126. }
  127. ///////////////////////////////////////////////////////////////////////////////
  128. // extract::lazy_variance
  129. // extract::variance
  130. //
  131. namespace extract
  132. {
  133. extractor<tag::lazy_variance> const lazy_variance = {};
  134. extractor<tag::variance> const variance = {};
  135. BOOST_ACCUMULATORS_IGNORE_GLOBAL(lazy_variance)
  136. BOOST_ACCUMULATORS_IGNORE_GLOBAL(variance)
  137. }
  138. using extract::lazy_variance;
  139. using extract::variance;
  140. // variance(lazy) -> lazy_variance
  141. template<>
  142. struct as_feature<tag::variance(lazy)>
  143. {
  144. typedef tag::lazy_variance type;
  145. };
  146. // variance(immediate) -> variance
  147. template<>
  148. struct as_feature<tag::variance(immediate)>
  149. {
  150. typedef tag::variance type;
  151. };
  152. // for the purposes of feature-based dependency resolution,
  153. // immediate_variance provides the same feature as variance
  154. template<>
  155. struct feature_of<tag::lazy_variance>
  156. : feature_of<tag::variance>
  157. {
  158. };
  159. // So that variance can be automatically substituted with
  160. // weighted_variance when the weight parameter is non-void.
  161. template<>
  162. struct as_weighted_feature<tag::variance>
  163. {
  164. typedef tag::weighted_variance type;
  165. };
  166. // for the purposes of feature-based dependency resolution,
  167. // weighted_variance provides the same feature as variance
  168. template<>
  169. struct feature_of<tag::weighted_variance>
  170. : feature_of<tag::variance>
  171. {
  172. };
  173. // So that immediate_variance can be automatically substituted with
  174. // immediate_weighted_variance when the weight parameter is non-void.
  175. template<>
  176. struct as_weighted_feature<tag::lazy_variance>
  177. {
  178. typedef tag::lazy_weighted_variance type;
  179. };
  180. // for the purposes of feature-based dependency resolution,
  181. // immediate_weighted_variance provides the same feature as immediate_variance
  182. template<>
  183. struct feature_of<tag::lazy_weighted_variance>
  184. : feature_of<tag::lazy_variance>
  185. {
  186. };
  187. ////////////////////////////////////////////////////////////////////////////
  188. //// droppable_accumulator<variance_impl>
  189. //// need to specialize droppable lazy variance to cache the result at the
  190. //// point the accumulator is dropped.
  191. ///// INTERNAL ONLY
  192. /////
  193. //template<typename Sample, typename MeanFeature>
  194. //struct droppable_accumulator<impl::variance_impl<Sample, MeanFeature> >
  195. // : droppable_accumulator_base<
  196. // with_cached_result<impl::variance_impl<Sample, MeanFeature> >
  197. // >
  198. //{
  199. // template<typename Args>
  200. // droppable_accumulator(Args const &args)
  201. // : droppable_accumulator::base(args)
  202. // {
  203. // }
  204. //};
  205. }} // namespace boost::accumulators
  206. #endif