libs/math/doc/html/math_toolkit/stat_tut/weg/neg_binom_eg/neg_binom_conf.html

   1 <html>
   2 <head>
   3 <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
   4 <title>Calculating Confidence Limits on the Frequency of Occurrence for the Negative Binomial Distribution</title>
   5 <link rel="stylesheet" href="../../../../math.css" type="text/css">
   6 <meta name="generator" content="DocBook XSL Stylesheets V1.78.1">
   7 <link rel="home" href="../../../../index.html" title="Math Toolkit 2.1.0">
   8 <link rel="up" href="../neg_binom_eg.html" title="Negative Binomial Distribution Examples">
   9 <link rel="prev" href="../neg_binom_eg.html" title="Negative Binomial Distribution Examples">
  10 <link rel="next" href="neg_binom_size_eg.html" title="Estimating Sample Sizes for the Negative Binomial.">
  11 </head>
  12 <body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF">
  13 <table cellpadding="2" width="100%"><tr>
  14 <td valign="top"><img alt="Boost C++ Libraries" width="277" height="86" src="../../../../../../../../boost.png"></td>
  15 <td align="center"><a href="../../../../../../../../index.html">Home</a></td>
  16 <td align="center"><a href="../../../../../../../../libs/libraries.htm">Libraries</a></td>
  17 <td align="center"><a href="http://www.boost.org/users/people.html">People</a></td>
  18 <td align="center"><a href="http://www.boost.org/users/faq.html">FAQ</a></td>
  19 <td align="center"><a href="../../../../../../../../more/index.htm">More</a></td>
  20 </tr></table>
  21 <hr>
  22 <div class="spirit-nav">
  23 <a accesskey="p" href="../neg_binom_eg.html"><img src="../../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../neg_binom_eg.html"><img src="../../../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../../../index.html"><img src="../../../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="neg_binom_size_eg.html"><img src="../../../../../../../../doc/src/images/next.png" alt="Next"></a>
  24 </div>
  25 <div class="section">
  26 <div class="titlepage"><div><div><h5 class="title">
  27 <a name="math_toolkit.stat_tut.weg.neg_binom_eg.neg_binom_conf"></a><a class="link" href="neg_binom_conf.html" title="Calculating Confidence Limits on the Frequency of Occurrence for the Negative Binomial Distribution">Calculating
  28           Confidence Limits on the Frequency of Occurrence for the Negative Binomial
  29           Distribution</a>
  30 </h5></div></div></div>
  31 <p>
  32             Imagine you have a process that follows a negative binomial distribution:
  33             for each trial conducted, an event either occurs or does it does not,
  34             referred to as "successes" and "failures". The frequency
  35             with which successes occur is variously referred to as the success fraction,
  36             success ratio, success percentage, occurrence frequency, or probability
  37             of occurrence.
  38           </p>
  39 <p>
  40             If, by experiment, you want to measure the the best estimate of success
  41             fraction is given simply by <span class="emphasis"><em>k</em></span> / <span class="emphasis"><em>N</em></span>,
  42             for <span class="emphasis"><em>k</em></span> successes out of <span class="emphasis"><em>N</em></span> trials.
  43           </p>
  44 <p>
  45             However our confidence in that estimate will be shaped by how many trials
  46             were conducted, and how many successes were observed. The static member
  47             functions <code class="computeroutput"><span class="identifier">negative_binomial_distribution</span><span class="special">&lt;&gt;::</span><span class="identifier">find_lower_bound_on_p</span></code>
  48             and <code class="computeroutput"><span class="identifier">negative_binomial_distribution</span><span class="special">&lt;&gt;::</span><span class="identifier">find_upper_bound_on_p</span></code>
  49             allow you to calculate the confidence intervals for your estimate of
  50             the success fraction.
  51           </p>
  52 <p>
  53             The sample program <a href="../../../../../../example/neg_binom_confidence_limits.cpp" target="_top">neg_binom_confidence_limits.cpp</a>
  54             illustrates their use.
  55           </p>
  56 <p>
  57             First we need some includes to access the negative binomial distribution
  58             (and some basic std output of course).
  59           </p>
  60 <pre class="programlisting"><span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">distributions</span><span class="special">/</span><span class="identifier">negative_binomial</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
  61 <span class="keyword">using</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">negative_binomial</span><span class="special">;</span>
  62
  63 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iostream</span><span class="special">&gt;</span>
  64 <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
  65 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iomanip</span><span class="special">&gt;</span>
  66 <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">setprecision</span><span class="special">;</span>
  67 <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">setw</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">left</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">fixed</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">right</span><span class="special">;</span>
  68 </pre>
  69 <p>
  70             First define a table of significance levels: these are the probabilities
  71             that the true occurrence frequency lies outside the calculated interval:
  72           </p>
  73 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">alpha</span><span class="special">[]</span> <span class="special">=</span> <span class="special">{</span> <span class="number">0.5</span><span class="special">,</span> <span class="number">0.25</span><span class="special">,</span> <span class="number">0.1</span><span class="special">,</span> <span class="number">0.05</span><span class="special">,</span> <span class="number">0.01</span><span class="special">,</span> <span class="number">0.001</span><span class="special">,</span> <span class="number">0.0001</span><span class="special">,</span> <span class="number">0.00001</span> <span class="special">};</span>
  74 </pre>
  75 <p>
  76             Confidence value as % is (1 - alpha) * 100, so alpha 0.05 == 95% confidence
  77             that the true occurence frequency lies <span class="bold"><strong>inside</strong></span>
  78             the calculated interval.
  79           </p>
  80 <p>
  81             We need a function to calculate and print confidence limits for an observed
  82             frequency of occurrence that follows a negative binomial distribution.
  83           </p>
  84 <pre class="programlisting"><span class="keyword">void</span> <span class="identifier">confidence_limits_on_frequency</span><span class="special">(</span><span class="keyword">unsigned</span> <span class="identifier">trials</span><span class="special">,</span> <span class="keyword">unsigned</span> <span class="identifier">successes</span><span class="special">)</span>
  85 <span class="special">{</span>
  86    <span class="comment">// trials = Total number of trials.</span>
  87    <span class="comment">// successes = Total number of observed successes.</span>
  88    <span class="comment">// failures = trials - successes.</span>
  89    <span class="comment">// success_fraction = successes /trials.</span>
  90    <span class="comment">// Print out general info:</span>
  91    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span>
  92       <span class="string">"______________________________________________\n"</span>
  93       <span class="string">"2-Sided Confidence Limits For Success Fraction\n"</span>
  94       <span class="string">"______________________________________________\n\n"</span><span class="special">;</span>
  95    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">7</span><span class="special">);</span>
  96    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">40</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</span> <span class="special">&lt;&lt;</span> <span class="string">"Number of trials"</span> <span class="special">&lt;&lt;</span> <span class="string">" =  "</span> <span class="special">&lt;&lt;</span> <span class="identifier">trials</span> <span class="special">&lt;&lt;</span> <span class="string">"\n"</span><span class="special">;</span>
  97    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">40</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</span> <span class="special">&lt;&lt;</span> <span class="string">"Number of successes"</span> <span class="special">&lt;&lt;</span> <span class="string">" =  "</span> <span class="special">&lt;&lt;</span> <span class="identifier">successes</span> <span class="special">&lt;&lt;</span> <span class="string">"\n"</span><span class="special">;</span>
  98    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">40</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</span> <span class="special">&lt;&lt;</span> <span class="string">"Number of failures"</span> <span class="special">&lt;&lt;</span> <span class="string">" =  "</span> <span class="special">&lt;&lt;</span> <span class="identifier">trials</span> <span class="special">-</span> <span class="identifier">successes</span> <span class="special">&lt;&lt;</span> <span class="string">"\n"</span><span class="special">;</span>
  99    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">40</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</span> <span class="special">&lt;&lt;</span> <span class="string">"Observed frequency of occurrence"</span> <span class="special">&lt;&lt;</span> <span class="string">" =  "</span> <span class="special">&lt;&lt;</span> <span class="keyword">double</span><span class="special">(</span><span class="identifier">successes</span><span class="special">)</span> <span class="special">/</span> <span class="identifier">trials</span> <span class="special">&lt;&lt;</span> <span class="string">"\n"</span><span class="special">;</span>
 100
 101    <span class="comment">// Print table header:</span>
 102    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"\n\n"</span>
 103            <span class="string">"___________________________________________\n"</span>
 104            <span class="string">"Confidence        Lower          Upper\n"</span>
 105            <span class="string">" Value (%)        Limit          Limit\n"</span>
 106            <span class="string">"___________________________________________\n"</span><span class="special">;</span>
 107 </pre>
 108 <p>
 109             And now for the important part - the bounds themselves. For each value
 110             of <span class="emphasis"><em>alpha</em></span>, we call <code class="computeroutput"><span class="identifier">find_lower_bound_on_p</span></code>
 111             and <code class="computeroutput"><span class="identifier">find_upper_bound_on_p</span></code>
 112             to obtain lower and upper bounds respectively. Note that since we are
 113             calculating a two-sided interval, we must divide the value of alpha in
 114             two. Had we been calculating a single-sided interval, for example: <span class="emphasis"><em>"Calculate
 115             a lower bound so that we are P% sure that the true occurrence frequency
 116             is greater than some value"</em></span> then we would <span class="bold"><strong>not</strong></span>
 117             have divided by two.
 118           </p>
 119 <pre class="programlisting">   <span class="comment">// Now print out the upper and lower limits for the alpha table values.</span>
 120    <span class="keyword">for</span><span class="special">(</span><span class="keyword">unsigned</span> <span class="identifier">i</span> <span class="special">=</span> <span class="number">0</span><span class="special">;</span> <span class="identifier">i</span> <span class="special">&lt;</span> <span class="keyword">sizeof</span><span class="special">(</span><span class="identifier">alpha</span><span class="special">)/</span><span class="keyword">sizeof</span><span class="special">(</span><span class="identifier">alpha</span><span class="special">[</span><span class="number">0</span><span class="special">]);</span> <span class="special">++</span><span class="identifier">i</span><span class="special">)</span>
 121    <span class="special">{</span>
 122       <span class="comment">// Confidence value:</span>
 123       <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">fixed</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">3</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">10</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">right</span> <span class="special">&lt;&lt;</span> <span class="number">100</span> <span class="special">*</span> <span class="special">(</span><span class="number">1</span><span class="special">-</span><span class="identifier">alpha</span><span class="special">[</span><span class="identifier">i</span><span class="special">]);</span>
 124       <span class="comment">// Calculate bounds:</span>
 125       <span class="keyword">double</span> <span class="identifier">lower</span> <span class="special">=</span> <span class="identifier">negative_binomial</span><span class="special">::</span><span class="identifier">find_lower_bound_on_p</span><span class="special">(</span><span class="identifier">trials</span><span class="special">,</span> <span class="identifier">successes</span><span class="special">,</span> <span class="identifier">alpha</span><span class="special">[</span><span class="identifier">i</span><span class="special">]/</span><span class="number">2</span><span class="special">);</span>
 126       <span class="keyword">double</span> <span class="identifier">upper</span> <span class="special">=</span> <span class="identifier">negative_binomial</span><span class="special">::</span><span class="identifier">find_upper_bound_on_p</span><span class="special">(</span><span class="identifier">trials</span><span class="special">,</span> <span class="identifier">successes</span><span class="special">,</span> <span class="identifier">alpha</span><span class="special">[</span><span class="identifier">i</span><span class="special">]/</span><span class="number">2</span><span class="special">);</span>
 127       <span class="comment">// Print limits:</span>
 128       <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">fixed</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">5</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">15</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">right</span> <span class="special">&lt;&lt;</span> <span class="identifier">lower</span><span class="special">;</span>
 129       <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">fixed</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">5</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">15</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">right</span> <span class="special">&lt;&lt;</span> <span class="identifier">upper</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 130    <span class="special">}</span>
 131    <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 132 <span class="special">}</span> <span class="comment">// void confidence_limits_on_frequency(unsigned trials, unsigned successes)</span>
 133 </pre>
 134 <p>
 135             And then call confidence_limits_on_frequency with increasing numbers
 136             of trials, but always the same success fraction 0.1, or 1 in 10.
 137           </p>
 138 <pre class="programlisting"><span class="keyword">int</span> <span class="identifier">main</span><span class="special">()</span>
 139 <span class="special">{</span>
 140   <span class="identifier">confidence_limits_on_frequency</span><span class="special">(</span><span class="number">20</span><span class="special">,</span> <span class="number">2</span><span class="special">);</span> <span class="comment">// 20 trials, 2 successes, 2 in 20, = 1 in 10 = 0.1 success fraction.</span>
 141   <span class="identifier">confidence_limits_on_frequency</span><span class="special">(</span><span class="number">200</span><span class="special">,</span> <span class="number">20</span><span class="special">);</span> <span class="comment">// More trials, but same 0.1 success fraction.</span>
 142   <span class="identifier">confidence_limits_on_frequency</span><span class="special">(</span><span class="number">2000</span><span class="special">,</span> <span class="number">200</span><span class="special">);</span> <span class="comment">// Many more trials, but same 0.1 success fraction.</span>
 143
 144   <span class="keyword">return</span> <span class="number">0</span><span class="special">;</span>
 145 <span class="special">}</span> <span class="comment">// int main()</span>
 146 </pre>
 147 <p>
 148             Let's see some sample output for a 1 in 10 success ratio, first for a
 149             mere 20 trials:
 150           </p>
 151 <pre class="programlisting">______________________________________________
 152 2-Sided Confidence Limits For Success Fraction
 153 ______________________________________________
 154 Number of trials                         =  20
 155 Number of successes                      =  2
 156 Number of failures                       =  18
 157 Observed frequency of occurrence         =  0.1
 158 ___________________________________________
 159 Confidence        Lower          Upper
 160  Value (%)        Limit          Limit
 161 ___________________________________________
 162     50.000        0.04812        0.13554
 163     75.000        0.03078        0.17727
 164     90.000        0.01807        0.22637
 165     95.000        0.01235        0.26028
 166     99.000        0.00530        0.33111
 167     99.900        0.00164        0.41802
 168     99.990        0.00051        0.49202
 169     99.999        0.00016        0.55574
 170 </pre>
 171 <p>
 172             As you can see, even at the 95% confidence level the bounds (0.012 to
 173             0.26) are really very wide, and very asymmetric about the observed value
 174             0.1.
 175           </p>
 176 <p>
 177             Compare that with the program output for a mass 2000 trials:
 178           </p>
 179 <pre class="programlisting">______________________________________________
 180 2-Sided Confidence Limits For Success Fraction
 181 ______________________________________________
 182 Number of trials                         =  2000
 183 Number of successes                      =  200
 184 Number of failures                       =  1800
 185 Observed frequency of occurrence         =  0.1
 186 ___________________________________________
 187 Confidence        Lower          Upper
 188  Value (%)        Limit          Limit
 189 ___________________________________________
 190     50.000        0.09536        0.10445
 191     75.000        0.09228        0.10776
 192     90.000        0.08916        0.11125
 193     95.000        0.08720        0.11352
 194     99.000        0.08344        0.11802
 195     99.900        0.07921        0.12336
 196     99.990        0.07577        0.12795
 197     99.999        0.07282        0.13206
 198 </pre>
 199 <p>
 200             Now even when the confidence level is very high, the limits (at 99.999%,
 201             0.07 to 0.13) are really quite close and nearly symmetric to the observed
 202             value of 0.1.
 203           </p>
 204 </div>
 205 <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
 206 <td align="left"></td>
 207 <td align="right"><div class="copyright-footer">Copyright &#169; 2006-2010, 2012-2014 Nikhar Agrawal,
 208       Anton Bikineev, Paul A. Bristow, Marco Guazzone, Christopher Kormanyos, Hubert
 209       Holin, Bruno Lalande, John Maddock, Johan R&#229;de, Gautam Sewani, Benjamin Sobotta,
 210       Thijs van den Berg, Daryle Walker and Xiaogang Zhang<p>
 211         Distributed under the Boost Software License, Version 1.0. (See accompanying
 212         file LICENSE_1_0.txt or copy at <a href="http://www.boost.org/LICENSE_1_0.txt" target="_top">http://www.boost.org/LICENSE_1_0.txt</a>)
 213       </p>
 214 </div></td>
 215 </tr></table>
 216 <hr>
 217 <div class="spirit-nav">
 218 <a accesskey="p" href="../neg_binom_eg.html"><img src="../../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../neg_binom_eg.html"><img src="../../../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../../../index.html"><img src="../../../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="neg_binom_size_eg.html"><img src="../../../../../../../../doc/src/images/next.png" alt="Next"></a>
 219 </div>
 220 </body>
 221 </html>