Statistics
| Revision:

root / tmp / org.txm.statsengine.r.core.win32 / res / win32 / library / BH / include / boost / accumulators / statistics / weighted_density.hpp @ 2486

History | View | Annotate | Download (8.8 kB)

1 2486 sjacqu01
///////////////////////////////////////////////////////////////////////////////
2 2486 sjacqu01
// weighted_density.hpp
3 2486 sjacqu01
//
4 2486 sjacqu01
//  Copyright 2006 Daniel Egloff, Olivier Gygi. Distributed under the Boost
5 2486 sjacqu01
//  Software License, Version 1.0. (See accompanying file
6 2486 sjacqu01
//  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
7 2486 sjacqu01
8 2486 sjacqu01
#ifndef BOOST_ACCUMULATORS_STATISTICS_WEIGHTED_DENSITY_HPP_DE_01_01_2006
9 2486 sjacqu01
#define BOOST_ACCUMULATORS_STATISTICS_WEIGHTED_DENSITY_HPP_DE_01_01_2006
10 2486 sjacqu01
11 2486 sjacqu01
#include <vector>
12 2486 sjacqu01
#include <limits>
13 2486 sjacqu01
#include <functional>
14 2486 sjacqu01
#include <boost/range.hpp>
15 2486 sjacqu01
#include <boost/parameter/keyword.hpp>
16 2486 sjacqu01
#include <boost/mpl/placeholders.hpp>
17 2486 sjacqu01
#include <boost/accumulators/framework/accumulator_base.hpp>
18 2486 sjacqu01
#include <boost/accumulators/framework/extractor.hpp>
19 2486 sjacqu01
#include <boost/accumulators/numeric/functional.hpp>
20 2486 sjacqu01
#include <boost/accumulators/framework/parameters/sample.hpp>
21 2486 sjacqu01
#include <boost/accumulators/statistics_fwd.hpp>
22 2486 sjacqu01
#include <boost/accumulators/statistics/sum.hpp>
23 2486 sjacqu01
#include <boost/accumulators/statistics/max.hpp>
24 2486 sjacqu01
#include <boost/accumulators/statistics/min.hpp>
25 2486 sjacqu01
#include <boost/accumulators/statistics/density.hpp> // for named parameters density_cache_size and density_num_bins
26 2486 sjacqu01
27 2486 sjacqu01
namespace boost { namespace accumulators
28 2486 sjacqu01
{
29 2486 sjacqu01
30 2486 sjacqu01
namespace impl
31 2486 sjacqu01
{
32 2486 sjacqu01
    ///////////////////////////////////////////////////////////////////////////////
33 2486 sjacqu01
    // weighted_density_impl
34 2486 sjacqu01
    //  density histogram for weighted samples
35 2486 sjacqu01
    /**
36 2486 sjacqu01
        @brief Histogram density estimator for weighted samples
37 2486 sjacqu01

38 2486 sjacqu01
        The histogram density estimator returns a histogram of the sample distribution. The positions and sizes of the bins
39 2486 sjacqu01
        are determined using a specifiable number of cached samples (cache_size). The range between the minimum and the
40 2486 sjacqu01
        maximum of the cached samples is subdivided into a specifiable number of bins (num_bins) of same size. Additionally,
41 2486 sjacqu01
        an under- and an overflow bin is added to capture future under- and overflow samples. Once the bins are determined,
42 2486 sjacqu01
        the cached samples and all subsequent samples are added to the correct bins. At the end, a range of std::pair is
43 2486 sjacqu01
        returned, where each pair contains the position of the bin (lower bound) and the sum of the weights (normalized with the
44 2486 sjacqu01
        sum of all weights).
45 2486 sjacqu01

46 2486 sjacqu01
        @param density_cache_size Number of first samples used to determine min and max.
47 2486 sjacqu01
        @param density_num_bins Number of bins (two additional bins collect under- and overflow samples).
48 2486 sjacqu01
    */
49 2486 sjacqu01
    template<typename Sample, typename Weight>
50 2486 sjacqu01
    struct weighted_density_impl
51 2486 sjacqu01
      : accumulator_base
52 2486 sjacqu01
    {
53 2486 sjacqu01
        typedef typename numeric::functional::fdiv<Weight, std::size_t>::result_type float_type;
54 2486 sjacqu01
        typedef std::vector<std::pair<float_type, float_type> > histogram_type;
55 2486 sjacqu01
        typedef std::vector<float_type> array_type;
56 2486 sjacqu01
        // for boost::result_of
57 2486 sjacqu01
        typedef iterator_range<typename histogram_type::iterator> result_type;
58 2486 sjacqu01
59 2486 sjacqu01
        template<typename Args>
60 2486 sjacqu01
        weighted_density_impl(Args const &args)
61 2486 sjacqu01
            : cache_size(args[density_cache_size])
62 2486 sjacqu01
            , cache(cache_size)
63 2486 sjacqu01
            , num_bins(args[density_num_bins])
64 2486 sjacqu01
            , samples_in_bin(num_bins + 2, 0.)
65 2486 sjacqu01
            , bin_positions(num_bins + 2)
66 2486 sjacqu01
            , histogram(
67 2486 sjacqu01
                num_bins + 2
68 2486 sjacqu01
              , std::make_pair(
69 2486 sjacqu01
                    numeric::fdiv(args[sample | Sample()],(std::size_t)1)
70 2486 sjacqu01
                  , numeric::fdiv(args[sample | Sample()],(std::size_t)1)
71 2486 sjacqu01
                )
72 2486 sjacqu01
              )
73 2486 sjacqu01
            , is_dirty(true)
74 2486 sjacqu01
        {
75 2486 sjacqu01
        }
76 2486 sjacqu01
77 2486 sjacqu01
        template<typename Args>
78 2486 sjacqu01
        void operator ()(Args const &args)
79 2486 sjacqu01
        {
80 2486 sjacqu01
            this->is_dirty = true;
81 2486 sjacqu01
82 2486 sjacqu01
            std::size_t cnt = count(args);
83 2486 sjacqu01
84 2486 sjacqu01
            // Fill up cache with cache_size first samples
85 2486 sjacqu01
            if (cnt <= this->cache_size)
86 2486 sjacqu01
            {
87 2486 sjacqu01
                this->cache[cnt - 1] = std::make_pair(args[sample], args[weight]);
88 2486 sjacqu01
            }
89 2486 sjacqu01
90 2486 sjacqu01
            // Once cache_size samples have been accumulated, create num_bins bins of same size between
91 2486 sjacqu01
            // the minimum and maximum of the cached samples as well as an under- and an overflow bin.
92 2486 sjacqu01
            // Store their lower bounds (bin_positions) and fill the bins with the cached samples (samples_in_bin).
93 2486 sjacqu01
            if (cnt == this->cache_size)
94 2486 sjacqu01
            {
95 2486 sjacqu01
                float_type minimum = numeric::fdiv((min)(args),(std::size_t)1);
96 2486 sjacqu01
                float_type maximum = numeric::fdiv((max)(args),(std::size_t)1);
97 2486 sjacqu01
                float_type bin_size = numeric::fdiv(maximum - minimum, this->num_bins);
98 2486 sjacqu01
99 2486 sjacqu01
                // determine bin positions (their lower bounds)
100 2486 sjacqu01
                for (std::size_t i = 0; i < this->num_bins + 2; ++i)
101 2486 sjacqu01
                {
102 2486 sjacqu01
                    this->bin_positions[i] = minimum + (i - 1.) * bin_size;
103 2486 sjacqu01
                }
104 2486 sjacqu01
105 2486 sjacqu01
                for (typename histogram_type::const_iterator iter = this->cache.begin(); iter != this->cache.end(); ++iter)
106 2486 sjacqu01
                {
107 2486 sjacqu01
                    if (iter->first < this->bin_positions[1])
108 2486 sjacqu01
                    {
109 2486 sjacqu01
                        this->samples_in_bin[0] += iter->second;
110 2486 sjacqu01
                    }
111 2486 sjacqu01
                    else if (iter->first >= this->bin_positions[this->num_bins + 1])
112 2486 sjacqu01
                    {
113 2486 sjacqu01
                        this->samples_in_bin[this->num_bins + 1] += iter->second;
114 2486 sjacqu01
                    }
115 2486 sjacqu01
                    else
116 2486 sjacqu01
                    {
117 2486 sjacqu01
                        typename array_type::iterator it = std::upper_bound(
118 2486 sjacqu01
                            this->bin_positions.begin()
119 2486 sjacqu01
                          , this->bin_positions.end()
120 2486 sjacqu01
                          , iter->first
121 2486 sjacqu01
                        );
122 2486 sjacqu01
123 2486 sjacqu01
                        std::size_t d = std::distance(this->bin_positions.begin(), it);
124 2486 sjacqu01
                        this->samples_in_bin[d - 1] += iter->second;
125 2486 sjacqu01
                    }
126 2486 sjacqu01
                }
127 2486 sjacqu01
            }
128 2486 sjacqu01
            // Add each subsequent sample to the correct bin
129 2486 sjacqu01
            else if (cnt > this->cache_size)
130 2486 sjacqu01
            {
131 2486 sjacqu01
                if (args[sample] < this->bin_positions[1])
132 2486 sjacqu01
                {
133 2486 sjacqu01
                    this->samples_in_bin[0] += args[weight];
134 2486 sjacqu01
                }
135 2486 sjacqu01
                else if (args[sample] >= this->bin_positions[this->num_bins + 1])
136 2486 sjacqu01
                {
137 2486 sjacqu01
                    this->samples_in_bin[this->num_bins + 1] += args[weight];
138 2486 sjacqu01
                }
139 2486 sjacqu01
                else
140 2486 sjacqu01
                {
141 2486 sjacqu01
                    typename array_type::iterator it = std::upper_bound(
142 2486 sjacqu01
                        this->bin_positions.begin()
143 2486 sjacqu01
                      , this->bin_positions.end()
144 2486 sjacqu01
                      , args[sample]
145 2486 sjacqu01
                    );
146 2486 sjacqu01
147 2486 sjacqu01
                    std::size_t d = std::distance(this->bin_positions.begin(), it);
148 2486 sjacqu01
                    this->samples_in_bin[d - 1] += args[weight];
149 2486 sjacqu01
                }
150 2486 sjacqu01
            }
151 2486 sjacqu01
        }
152 2486 sjacqu01
153 2486 sjacqu01
        template<typename Args>
154 2486 sjacqu01
        result_type result(Args const &args) const
155 2486 sjacqu01
        {
156 2486 sjacqu01
            if (this->is_dirty)
157 2486 sjacqu01
            {
158 2486 sjacqu01
                this->is_dirty = false;
159 2486 sjacqu01
160 2486 sjacqu01
                // creates a vector of std::pair where each pair i holds
161 2486 sjacqu01
                // the values bin_positions[i] (x-axis of histogram) and
162 2486 sjacqu01
                // samples_in_bin[i] / cnt (y-axis of histogram).
163 2486 sjacqu01
164 2486 sjacqu01
                for (std::size_t i = 0; i < this->num_bins + 2; ++i)
165 2486 sjacqu01
                {
166 2486 sjacqu01
                    this->histogram[i] = std::make_pair(this->bin_positions[i], numeric::fdiv(this->samples_in_bin[i], sum_of_weights(args)));
167 2486 sjacqu01
                }
168 2486 sjacqu01
            }
169 2486 sjacqu01
170 2486 sjacqu01
            // returns a range of pairs
171 2486 sjacqu01
            return make_iterator_range(this->histogram);
172 2486 sjacqu01
        }
173 2486 sjacqu01
174 2486 sjacqu01
    private:
175 2486 sjacqu01
        std::size_t            cache_size;      // number of cached samples
176 2486 sjacqu01
        histogram_type         cache;           // cache to store the first cache_size samples with their weights as std::pair
177 2486 sjacqu01
        std::size_t            num_bins;        // number of bins
178 2486 sjacqu01
        array_type             samples_in_bin;  // number of samples in each bin
179 2486 sjacqu01
        array_type             bin_positions;   // lower bounds of bins
180 2486 sjacqu01
        mutable histogram_type histogram;       // histogram
181 2486 sjacqu01
        mutable bool is_dirty;
182 2486 sjacqu01
    };
183 2486 sjacqu01
184 2486 sjacqu01
} // namespace impl
185 2486 sjacqu01
186 2486 sjacqu01
///////////////////////////////////////////////////////////////////////////////
187 2486 sjacqu01
// tag::weighted_density
188 2486 sjacqu01
//
189 2486 sjacqu01
namespace tag
190 2486 sjacqu01
{
191 2486 sjacqu01
    struct weighted_density
192 2486 sjacqu01
      : depends_on<count, sum_of_weights, min, max>
193 2486 sjacqu01
      , density_cache_size
194 2486 sjacqu01
      , density_num_bins
195 2486 sjacqu01
    {
196 2486 sjacqu01
        /// INTERNAL ONLY
197 2486 sjacqu01
        ///
198 2486 sjacqu01
        typedef accumulators::impl::weighted_density_impl<mpl::_1, mpl::_2> impl;
199 2486 sjacqu01
200 2486 sjacqu01
        #ifdef BOOST_ACCUMULATORS_DOXYGEN_INVOKED
201 2486 sjacqu01
        static boost::parameter::keyword<density_cache_size> const cache_size;
202 2486 sjacqu01
        static boost::parameter::keyword<density_num_bins> const num_bins;
203 2486 sjacqu01
        #endif
204 2486 sjacqu01
    };
205 2486 sjacqu01
}
206 2486 sjacqu01
207 2486 sjacqu01
///////////////////////////////////////////////////////////////////////////////
208 2486 sjacqu01
// extract::weighted_density
209 2486 sjacqu01
//
210 2486 sjacqu01
namespace extract
211 2486 sjacqu01
{
212 2486 sjacqu01
    extractor<tag::density> const weighted_density = {};
213 2486 sjacqu01
214 2486 sjacqu01
    BOOST_ACCUMULATORS_IGNORE_GLOBAL(weighted_density)
215 2486 sjacqu01
}
216 2486 sjacqu01
217 2486 sjacqu01
using extract::weighted_density;
218 2486 sjacqu01
219 2486 sjacqu01
}} // namespace boost::accumulators
220 2486 sjacqu01
221 2486 sjacqu01
#endif