histogram/utility/wald_interval.hpp

0001 // Copyright 2022 Jay Gohil, Hans Dembinski
0002 //
0003 // Distributed under the Boost Software License, version 1.0.
0004 // (See accompanying file LICENSE_1_0.txt
0005 // or copy at http://www.boost.org/LICENSE_1_0.txt)
0006
0007 #ifndef BOOST_HISTOGRAM_UTILITY_WALD_INTERVAL_HPP
0008 #define BOOST_HISTOGRAM_UTILITY_WALD_INTERVAL_HPP
0009
0010 #include <boost/histogram/fwd.hpp>
0011 #include <boost/histogram/utility/binomial_proportion_interval.hpp>
0012 #include <cmath>
0013 #include <utility>
0014
0015 namespace boost {
0016 namespace histogram {
0017 namespace utility {
0018
0019 /**
0020   Wald interval or normal approximation interval.
0021
0022   The Wald interval is a symmetric interval. It is simple to compute, but has poor
0023   statistical properties and is universally rejected by statisticians. It should always be
0024   replaced by another iternal, for example, the Wilson interval.
0025
0026   The Wald interval can be derived easily using the plug-in estimate of the variance for
0027   the binomial distribution, which is likely a reason for its omnipresence. Without
0028   further insight into statistical theory, it is not obvious that this derivation is
0029   flawed and that better alternatives exist.
0030
0031   The Wald interval undercovers on average. It is unsuitable when the sample size is small
0032   or when the fraction is close to 0 or 1. e. Its limits are not naturally bounded by 0
0033   or 1. It produces empty intervals if the number of successes or failures is zero.
0034
0035   For a critique of the Wald interval, see (a selection):
0036
0037   L.D. Brown, T.T. Cai, A. DasGupta, Statistical Science 16 (2001) 101-133.
0038   R. D. Cousins, K. E. Hymes, J. Tucker, Nucl. Instrum. Meth. A 612 (2010) 388-398.
0039 */
0040 template <class ValueType>
0041 class wald_interval : public binomial_proportion_interval<ValueType> {
0042 public:
0043   using value_type = typename wald_interval::value_type;
0044   using interval_type = typename wald_interval::interval_type;
0045
0046   /** Construct Wald interval computer.
0047
0048     @param d Number of standard deviations for the interval. The default value 1
0049     corresponds to a confidence level of 68 %. Both `deviation` and `confidence_level`
0050     objects can be used to initialize the interval.
0051   */
0052   explicit wald_interval(deviation d = deviation{1.0}) noexcept
0053       : z_{static_cast<value_type>(d)} {}
0054
0055   using binomial_proportion_interval<ValueType>::operator();
0056
0057   /** Compute interval for given number of successes and failures.
0058
0059     @param successes Number of successful trials.
0060     @param failures Number of failed trials.
0061   */
0062   interval_type operator()(value_type successes, value_type failures) const noexcept {
0063     // See https://en.wikipedia.org/wiki/
0064     //   Binomial_proportion_confidence_interval
0065     //   #Normal_approximation_interval_or_Wald_interval
0066     const value_type total_inv = 1 / (successes + failures);
0067     const value_type a = successes * total_inv;
0068     const value_type b = (z_ * total_inv) * std::sqrt(successes * failures * total_inv);
0069     return {a - b, a + b};
0070   }
0071
0072 private:
0073   value_type z_;
0074 };
0075
0076 } // namespace utility
0077 } // namespace histogram
0078 } // namespace boost
0079
0080 #endif