boost/math/distributions/empirical_cumulative_distribution_function.hpp

   1 //  Copyright Nick Thompson 2019.
   2 //  Use, modification and distribution are subject to the
   3 //  Boost Software License, Version 1.0. (See accompanying file
   4 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   5
   6 #ifndef BOOST_MATH_DISTRIBUTIONS_EMPIRICAL_CUMULATIVE_DISTRIBUTION_FUNCTION_HPP
   7 #define BOOST_MATH_DISTRIBUTIONS_EMPIRICAL_CUMULATIVE_DISTRIBUTION_FUNCTION_HPP
   8 #include <algorithm>
   9 #include <iterator>
  10 #include <stdexcept>
  11
  12 namespace boost { namespace math{
  13
  14 template<class RandomAccessContainer>
  15 class empirical_cumulative_distribution_function {
  16     using Real = typename RandomAccessContainer::value_type;
  17 public:
  18     empirical_cumulative_distribution_function(RandomAccessContainer && v, bool sorted = false)
  19     {
  20         if (v.size() == 0) {
  21             throw std::domain_error("At least one sample is required to compute an empirical CDF.");
  22         }
  23         m_v = std::move(v);
  24         if (!sorted) {
  25             std::sort(m_v.begin(), m_v.end());
  26         }
  27     }
  28
  29     auto operator()(Real x) const {
  30        if constexpr (std::is_integral_v<Real>)
  31        {
  32          if (x < m_v[0]) {
  33            return double(0);
  34          }
  35          if (x >= m_v[m_v.size()-1]) {
  36            return double(1);
  37          }
  38          auto it = std::upper_bound(m_v.begin(), m_v.end(), x);
  39          return static_cast<double>(std::distance(m_v.begin(), it))/static_cast<double>(m_v.size());
  40        }
  41        else
  42        {
  43          if (x < m_v[0]) {
  44            return Real(0);
  45          }
  46          if (x >= m_v[m_v.size()-1]) {
  47            return Real(1);
  48          }
  49          auto it = std::upper_bound(m_v.begin(), m_v.end(), x);
  50          return static_cast<Real>(std::distance(m_v.begin(), it))/static_cast<Real>(m_v.size());
  51       }
  52     }
  53
  54     RandomAccessContainer&& return_data() {
  55         return std::move(m_v);
  56     }
  57
  58 private:
  59     RandomAccessContainer m_v;
  60 };
  61
  62 }}
  63 #endif