libs/math/doc/html/math_toolkit/stat_tut/weg/normal_example/normal_misc.html

   1 <html>
   2 <head>
   3 <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
   4 <title>Some Miscellaneous Examples of the Normal (Gaussian) Distribution</title>
   5 <link rel="stylesheet" href="../../../../math.css" type="text/css">
   6 <meta name="generator" content="DocBook XSL Stylesheets V1.79.1">
   7 <link rel="home" href="../../../../index.html" title="Math Toolkit 2.11.0">
   8 <link rel="up" href="../normal_example.html" title="Normal Distribution Examples">
   9 <link rel="prev" href="../normal_example.html" title="Normal Distribution Examples">
  10 <link rel="next" href="../inverse_chi_squared_eg.html" title="Inverse Chi-Squared Distribution Bayes Example">
  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="../normal_example.html"><img src="../../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../normal_example.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="../inverse_chi_squared_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.normal_example.normal_misc"></a><a class="link" href="normal_misc.html" title="Some Miscellaneous Examples of the Normal (Gaussian) Distribution">Some
  28           Miscellaneous Examples of the Normal (Gaussian) Distribution</a>
  29 </h5></div></div></div>
  30 <p>
  31             The sample program <a href="../../../../../../example/normal_misc_examples.cpp" target="_top">normal_misc_examples.cpp</a>
  32             illustrates their use.
  33           </p>
  34 <h5>
  35 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h0"></a>
  36             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.traditional_tables"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.traditional_tables">Traditional
  37             Tables</a>
  38           </h5>
  39 <p>
  40             First we need some includes to access the normal distribution (and some
  41             std output of course).
  42           </p>
  43 <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">normal</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span> <span class="comment">// for normal_distribution</span>
  44   <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">normal</span><span class="special">;</span> <span class="comment">// typedef provides default type is double.</span>
  45
  46 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iostream</span><span class="special">&gt;</span>
  47   <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> <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">showpoint</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">noshowpoint</span><span class="special">;</span>
  48 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iomanip</span><span class="special">&gt;</span>
  49   <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">setprecision</span><span class="special">;</span>
  50 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">limits</span><span class="special">&gt;</span>
  51   <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special">;</span>
  52
  53 <span class="keyword">int</span> <span class="identifier">main</span><span class="special">()</span>
  54 <span class="special">{</span>
  55   <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Example: Normal distribution, Miscellaneous Applications."</span><span class="special">;</span>
  56
  57   <span class="keyword">try</span>
  58   <span class="special">{</span>
  59     <span class="special">{</span> <span class="comment">// Traditional tables and values.</span>
  60 </pre>
  61 <p>
  62             Let's start by printing some traditional tables.
  63           </p>
  64 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">step</span> <span class="special">=</span> <span class="number">1.</span><span class="special">;</span> <span class="comment">// in z</span>
  65 <span class="keyword">double</span> <span class="identifier">range</span> <span class="special">=</span> <span class="number">4</span><span class="special">;</span> <span class="comment">// min and max z = -range to +range.</span>
  66 <span class="keyword">int</span> <span class="identifier">precision</span> <span class="special">=</span> <span class="number">17</span><span class="special">;</span> <span class="comment">// traditional tables are only computed to much lower precision.</span>
  67 <span class="comment">// but std::numeric_limits&lt;double&gt;::max_digits10; on new Standard Libraries gives</span>
  68 <span class="comment">// 17, the maximum number of digits that can possibly be significant.</span>
  69 <span class="comment">// std::numeric_limits&lt;double&gt;::digits10; == 15 is number of guaranteed digits,</span>
  70 <span class="comment">// the other two digits being 'noisy'.</span>
  71
  72 <span class="comment">// Construct a standard normal distribution s</span>
  73   <span class="identifier">normal</span> <span class="identifier">s</span><span class="special">;</span> <span class="comment">// (default mean = zero, and standard deviation = unity)</span>
  74   <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Standard normal distribution, mean = "</span><span class="special">&lt;&lt;</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span>
  75     <span class="special">&lt;&lt;</span> <span class="string">", standard deviation = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  76 </pre>
  77 <p>
  78             First the probability distribution function (pdf).
  79           </p>
  80 <pre class="programlisting"><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Probability distribution function values"</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  81 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"  z "</span> <span class="string">"      pdf "</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  82 <span class="identifier">cout</span><span class="special">.</span><span class="identifier">precision</span><span class="special">(</span><span class="number">5</span><span class="special">);</span>
  83 <span class="keyword">for</span> <span class="special">(</span><span class="keyword">double</span> <span class="identifier">z</span> <span class="special">=</span> <span class="special">-</span><span class="identifier">range</span><span class="special">;</span> <span class="identifier">z</span> <span class="special">&lt;</span> <span class="identifier">range</span> <span class="special">+</span> <span class="identifier">step</span><span class="special">;</span> <span class="identifier">z</span> <span class="special">+=</span> <span class="identifier">step</span><span class="special">)</span>
  84 <span class="special">{</span>
  85   <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</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">6</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">z</span> <span class="special">&lt;&lt;</span> <span class="string">" "</span>
  86     <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="identifier">precision</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">12</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">pdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">z</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  87 <span class="special">}</span>
  88 <span class="identifier">cout</span><span class="special">.</span><span class="identifier">precision</span><span class="special">(</span><span class="number">6</span><span class="special">);</span> <span class="comment">// default</span>
  89 </pre>
  90 <p>
  91             And the area under the normal curve from -&#8734; up to z, the cumulative distribution
  92             function (cdf).
  93           </p>
  94 <pre class="programlisting"><span class="comment">// For a standard normal distribution</span>
  95 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Standard normal mean = "</span><span class="special">&lt;&lt;</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span>
  96   <span class="special">&lt;&lt;</span> <span class="string">", standard deviation = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  97 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Integral (area under the curve) from - infinity up to z "</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  98 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"  z "</span> <span class="string">"      cdf "</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
  99 <span class="keyword">for</span> <span class="special">(</span><span class="keyword">double</span> <span class="identifier">z</span> <span class="special">=</span> <span class="special">-</span><span class="identifier">range</span><span class="special">;</span> <span class="identifier">z</span> <span class="special">&lt;</span> <span class="identifier">range</span> <span class="special">+</span> <span class="identifier">step</span><span class="special">;</span> <span class="identifier">z</span> <span class="special">+=</span> <span class="identifier">step</span><span class="special">)</span>
 100 <span class="special">{</span>
 101   <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">left</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">6</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">z</span> <span class="special">&lt;&lt;</span> <span class="string">" "</span>
 102     <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="identifier">precision</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">setw</span><span class="special">(</span><span class="number">12</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">z</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 103 <span class="special">}</span>
 104 <span class="identifier">cout</span><span class="special">.</span><span class="identifier">precision</span><span class="special">(</span><span class="number">6</span><span class="special">);</span> <span class="comment">// default</span>
 105 </pre>
 106 <p>
 107             And all this you can do with a nanoscopic amount of work compared to
 108             the team of <span class="bold"><strong>human computers</strong></span> toiling
 109             with Milton Abramovitz and Irene Stegen at the US National Bureau of
 110             Standards (now <a href="http://www.nist.gov" target="_top">NIST</a>). Starting
 111             in 1938, their "Handbook of Mathematical Functions with Formulas,
 112             Graphs and Mathematical Tables", was eventually published in 1964,
 113             and has been reprinted numerous times since. (A major replacement is
 114             planned at <a href="http://dlmf.nist.gov" target="_top">Digital Library of Mathematical
 115             Functions</a>).
 116           </p>
 117 <p>
 118             Pretty-printing a traditional 2-dimensional table is left as an exercise
 119             for the student, but why bother now that the Math Toolkit lets you write
 120           </p>
 121 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">z</span> <span class="special">=</span> <span class="number">2.</span><span class="special">;</span>
 122 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Area for z = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">z</span> <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">z</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">// to get the area for z.</span>
 123 </pre>
 124 <p>
 125             Correspondingly, we can obtain the traditional 'critical' values for
 126             significance levels. For the 95% confidence level, the significance level
 127             usually called alpha, is 0.05 = 1 - 0.95 (for a one-sided test), so we
 128             can write
 129           </p>
 130 <pre class="programlisting">  <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"95% of area has a z below "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">0.95</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 131 <span class="comment">// 95% of area has a z below 1.64485</span>
 132 </pre>
 133 <p>
 134             and a two-sided test (a comparison between two levels, rather than a
 135             one-sided test)
 136           </p>
 137 <pre class="programlisting">  <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"95% of area has a z between "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">0.975</span><span class="special">)</span>
 138     <span class="special">&lt;&lt;</span> <span class="string">" and "</span> <span class="special">&lt;&lt;</span> <span class="special">-</span><span class="identifier">quantile</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">0.975</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 139 <span class="comment">// 95% of area has a z between 1.95996 and -1.95996</span>
 140 </pre>
 141 <p>
 142             First, define a table of significance levels: these are the probabilities
 143             that the true occurrence frequency lies outside the calculated interval.
 144           </p>
 145 <p>
 146             It is convenient to have an alpha level for the probability that z lies
 147             outside just one standard deviation. This will not be some nice neat
 148             number like 0.05, but we can easily calculate it,
 149           </p>
 150 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">alpha1</span> <span class="special">=</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="special">-</span><span class="number">1</span><span class="special">)</span> <span class="special">*</span> <span class="number">2</span><span class="special">;</span> <span class="comment">// 0.3173105078629142</span>
 151 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">17</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="string">"Significance level for z == 1 is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">alpha1</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 152 </pre>
 153 <p>
 154             and place in our array of favorite alpha values.
 155           </p>
 156 <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.3173105078629142</span><span class="special">,</span> <span class="comment">// z for 1 standard deviation.</span>
 157   <span class="number">0.20</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>
 158 </pre>
 159 <p>
 160             Confidence value as % is (1 - alpha) * 100 (so alpha 0.05 == 95% confidence)
 161             that the true occurrence frequency lies <span class="bold"><strong>inside</strong></span>
 162             the calculated interval.
 163           </p>
 164 <pre class="programlisting"><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"level of significance (alpha)"</span> <span class="special">&lt;&lt;</span> <span class="identifier">setprecision</span><span class="special">(</span><span class="number">4</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 165 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"2-sided       1 -sided          z(alpha) "</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 166 <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>
 167 <span class="special">{</span>
 168   <span class="identifier">cout</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">alpha</span><span class="special">[</span><span class="identifier">i</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">alpha</span><span class="special">[</span><span class="identifier">i</span><span class="special">]</span> <span class="special">/</span><span class="number">2</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">quantile</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">s</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> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 169   <span class="comment">// Use quantile(complement(s, alpha[i]/2)) to avoid potential loss of accuracy from quantile(s,  1 - alpha[i]/2)</span>
 170 <span class="special">}</span>
 171 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 172 </pre>
 173 <p>
 174             Notice the distinction between one-sided (also called one-tailed) where
 175             we are using a &gt; <span class="bold"><strong>or</strong></span> &lt; test (and
 176             not both) and considering the area of the tail (integral) from z up to
 177             +&#8734;, and a two-sided test where we are using two &gt; <span class="bold"><strong>and</strong></span>
 178             &lt; tests, and thus considering two tails, from -&#8734; up to z low and z high
 179             up to +&#8734;.
 180           </p>
 181 <p>
 182             So the 2-sided values alpha[i] are calculated using alpha[i]/2.
 183           </p>
 184 <p>
 185             If we consider a simple example of alpha = 0.05, then for a two-sided
 186             test, the lower tail area from -&#8734; up to -1.96 is 0.025 (alpha/2) and the
 187             upper tail area from +z up to +1.96 is also 0.025 (alpha/2), and the
 188             area between -1.96 up to 12.96 is alpha = 0.95. and the sum of the two
 189             tails is 0.025 + 0.025 = 0.05,
 190           </p>
 191 <h5>
 192 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h1"></a>
 193             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.standard_deviations_either_side_"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.standard_deviations_either_side_">Standard
 194             deviations either side of the Mean</a>
 195           </h5>
 196 <p>
 197             Armed with the cumulative distribution function, we can easily calculate
 198             the easy to remember proportion of values that lie within 1, 2 and 3
 199             standard deviations from the mean.
 200           </p>
 201 <pre class="programlisting"><span class="identifier">cout</span><span class="special">.</span><span class="identifier">precision</span><span class="special">(</span><span class="number">3</span><span class="special">);</span>
 202 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">showpoint</span> <span class="special">&lt;&lt;</span> <span class="string">"cdf(s, s.standard_deviation()) = "</span>
 203   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">())</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>  <span class="comment">// from -infinity to 1 sd</span>
 204 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"cdf(complement(s, s.standard_deviation())) = "</span>
 205   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 206 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction 1 standard deviation within either side of mean is "</span>
 207   <span class="special">&lt;&lt;</span> <span class="number">1</span> <span class="special">-</span>  <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()))</span> <span class="special">*</span> <span class="number">2</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 208 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction 2 standard deviations within either side of mean is "</span>
 209   <span class="special">&lt;&lt;</span> <span class="number">1</span> <span class="special">-</span>  <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">2</span> <span class="special">*</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()))</span> <span class="special">*</span> <span class="number">2</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 210 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction 3 standard deviations within either side of mean is "</span>
 211   <span class="special">&lt;&lt;</span> <span class="number">1</span> <span class="special">-</span>  <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">3</span> <span class="special">*</span> <span class="identifier">s</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()))</span> <span class="special">*</span> <span class="number">2</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 212 </pre>
 213 <p>
 214             To a useful precision, the 1, 2 &amp; 3 percentages are 68, 95 and 99.7,
 215             and these are worth memorising as useful 'rules of thumb', as, for example,
 216             in <a href="http://en.wikipedia.org/wiki/Standard_deviation" target="_top">standard
 217             deviation</a>:
 218           </p>
 219 <pre class="programlisting">Fraction 1 standard deviation within either side of mean is 0.683
 220 Fraction 2 standard deviations within either side of mean is 0.954
 221 Fraction 3 standard deviations within either side of mean is 0.997
 222 </pre>
 223 <p>
 224             We could of course get some really accurate values for these <a href="http://en.wikipedia.org/wiki/Confidence_interval" target="_top">confidence
 225             intervals</a> by using cout.precision(15);
 226           </p>
 227 <pre class="programlisting">Fraction 1 standard deviation within either side of mean is 0.682689492137086
 228 Fraction 2 standard deviations within either side of mean is 0.954499736103642
 229 Fraction 3 standard deviations within either side of mean is 0.997300203936740
 230 </pre>
 231 <p>
 232             But before you get too excited about this impressive precision, don't
 233             forget that the <span class="bold"><strong>confidence intervals of the standard
 234             deviation</strong></span> are surprisingly wide, especially if you have estimated
 235             the standard deviation from only a few measurements.
 236           </p>
 237 <h5>
 238 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h2"></a>
 239             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.some_simple_examples"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.some_simple_examples">Some
 240             simple examples</a>
 241           </h5>
 242 <h5>
 243 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h3"></a>
 244             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.life_of_light_bulbs"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.life_of_light_bulbs">Life
 245             of light bulbs</a>
 246           </h5>
 247 <p>
 248             Examples from K. Krishnamoorthy, Handbook of Statistical Distributions
 249             with Applications, ISBN 1 58488 635 8, page 125... implemented using
 250             the Math Toolkit library.
 251           </p>
 252 <p>
 253             A few very simple examples are shown here:
 254           </p>
 255 <pre class="programlisting"><span class="comment">// K. Krishnamoorthy, Handbook of Statistical Distributions with Applications,</span>
 256  <span class="comment">// ISBN 1 58488 635 8, page 125, example 10.3.5</span>
 257 </pre>
 258 <p>
 259             Mean lifespan of 100 W bulbs is 1100 h with standard deviation of 100
 260             h. Assuming, perhaps with little evidence and much faith, that the distribution
 261             is normal, we construct a normal distribution called <span class="emphasis"><em>bulbs</em></span>
 262             with these values:
 263           </p>
 264 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">mean_life</span> <span class="special">=</span> <span class="number">1100.</span><span class="special">;</span>
 265 <span class="keyword">double</span> <span class="identifier">life_standard_deviation</span> <span class="special">=</span> <span class="number">100.</span><span class="special">;</span>
 266 <span class="identifier">normal</span> <span class="identifier">bulbs</span><span class="special">(</span><span class="identifier">mean_life</span><span class="special">,</span> <span class="identifier">life_standard_deviation</span><span class="special">);</span>
 267 <span class="keyword">double</span> <span class="identifier">expected_life</span> <span class="special">=</span> <span class="number">1000.</span><span class="special">;</span>
 268 </pre>
 269 <p>
 270             The we can use the Cumulative distribution function to predict fractions
 271             (or percentages, if * 100) that will last various lifetimes.
 272           </p>
 273 <pre class="programlisting"><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction of bulbs that will last at best (&lt;=) "</span> <span class="comment">// P(X &lt;= 1000)</span>
 274   <span class="special">&lt;&lt;</span> <span class="identifier">expected_life</span> <span class="special">&lt;&lt;</span> <span class="string">" is "</span><span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bulbs</span><span class="special">,</span> <span class="identifier">expected_life</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 275 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction of bulbs that will last at least (&gt;) "</span> <span class="comment">// P(X &gt; 1000)</span>
 276   <span class="special">&lt;&lt;</span> <span class="identifier">expected_life</span> <span class="special">&lt;&lt;</span> <span class="string">" is "</span><span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">bulbs</span><span class="special">,</span> <span class="identifier">expected_life</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 277 <span class="keyword">double</span> <span class="identifier">min_life</span> <span class="special">=</span> <span class="number">900</span><span class="special">;</span>
 278 <span class="keyword">double</span> <span class="identifier">max_life</span> <span class="special">=</span> <span class="number">1200</span><span class="special">;</span>
 279 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction of bulbs that will last between "</span>
 280   <span class="special">&lt;&lt;</span> <span class="identifier">min_life</span> <span class="special">&lt;&lt;</span> <span class="string">" and "</span> <span class="special">&lt;&lt;</span> <span class="identifier">max_life</span> <span class="special">&lt;&lt;</span> <span class="string">" is "</span>
 281   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bulbs</span><span class="special">,</span> <span class="identifier">max_life</span><span class="special">)</span>  <span class="comment">// P(X &lt;= 1200)</span>
 282    <span class="special">-</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bulbs</span><span class="special">,</span> <span class="identifier">min_life</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">// P(X &lt;= 900)</span>
 283 </pre>
 284 <div class="note"><table border="0" summary="Note">
 285 <tr>
 286 <td rowspan="2" align="center" valign="top" width="25"><img alt="[Note]" src="../../../../../../../../doc/src/images/note.png"></td>
 287 <th align="left">Note</th>
 288 </tr>
 289 <tr><td align="left" valign="top"><p>
 290               Real-life failures are often very ab-normal, with a significant number
 291               that 'dead-on-arrival' or suffer failure very early in their life:
 292               the lifetime of the survivors of 'early mortality' may be well described
 293               by the normal distribution.
 294             </p></td></tr>
 295 </table></div>
 296 <h5>
 297 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h4"></a>
 298             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.how_many_onions"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.how_many_onions">How
 299             many onions?</a>
 300           </h5>
 301 <p>
 302             Weekly demand for 5 lb sacks of onions at a store is normally distributed
 303             with mean 140 sacks and standard deviation 10.
 304           </p>
 305 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">mean</span> <span class="special">=</span> <span class="number">140.</span><span class="special">;</span> <span class="comment">// sacks per week.</span>
 306 <span class="keyword">double</span> <span class="identifier">standard_deviation</span> <span class="special">=</span> <span class="number">10</span><span class="special">;</span>
 307 <span class="identifier">normal</span> <span class="identifier">sacks</span><span class="special">(</span><span class="identifier">mean</span><span class="special">,</span> <span class="identifier">standard_deviation</span><span class="special">);</span>
 308
 309 <span class="keyword">double</span> <span class="identifier">stock</span> <span class="special">=</span> <span class="number">160.</span><span class="special">;</span> <span class="comment">// per week.</span>
 310 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Percentage of weeks overstocked "</span>
 311   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">sacks</span><span class="special">,</span> <span class="identifier">stock</span><span class="special">)</span> <span class="special">*</span> <span class="number">100.</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">// P(X &lt;=160)</span>
 312 <span class="comment">// Percentage of weeks overstocked 97.7</span>
 313 </pre>
 314 <p>
 315             So there will be lots of mouldy onions! So we should be able to say what
 316             stock level will meet demand 95% of the weeks.
 317           </p>
 318 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">stock_95</span> <span class="special">=</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">sacks</span><span class="special">,</span> <span class="number">0.95</span><span class="special">);</span>
 319 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Store should stock "</span> <span class="special">&lt;&lt;</span> <span class="keyword">int</span><span class="special">(</span><span class="identifier">stock_95</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="string">" sacks to meet 95% of demands."</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 320 </pre>
 321 <p>
 322             And it is easy to estimate how to meet 80% of demand, and waste even
 323             less.
 324           </p>
 325 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">stock_80</span> <span class="special">=</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">sacks</span><span class="special">,</span> <span class="number">0.80</span><span class="special">);</span>
 326 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Store should stock "</span> <span class="special">&lt;&lt;</span> <span class="keyword">int</span><span class="special">(</span><span class="identifier">stock_80</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="string">" sacks to meet 8 out of 10 demands."</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 327 </pre>
 328 <h5>
 329 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h5"></a>
 330             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.packing_beef"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.packing_beef">Packing
 331             beef</a>
 332           </h5>
 333 <p>
 334             A machine is set to pack 3 kg of ground beef per pack. Over a long period
 335             of time it is found that the average packed was 3 kg with a standard
 336             deviation of 0.1 kg. Assuming the packing is normally distributed, we
 337             can find the fraction (or %) of packages that weigh more than 3.1 kg.
 338           </p>
 339 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">mean</span> <span class="special">=</span> <span class="number">3.</span><span class="special">;</span> <span class="comment">// kg</span>
 340 <span class="keyword">double</span> <span class="identifier">standard_deviation</span> <span class="special">=</span> <span class="number">0.1</span><span class="special">;</span> <span class="comment">// kg</span>
 341 <span class="identifier">normal</span> <span class="identifier">packs</span><span class="special">(</span><span class="identifier">mean</span><span class="special">,</span> <span class="identifier">standard_deviation</span><span class="special">);</span>
 342
 343 <span class="keyword">double</span> <span class="identifier">max_weight</span> <span class="special">=</span> <span class="number">3.1</span><span class="special">;</span> <span class="comment">// kg</span>
 344 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Percentage of packs &gt; "</span> <span class="special">&lt;&lt;</span> <span class="identifier">max_weight</span> <span class="special">&lt;&lt;</span> <span class="string">" is "</span>
 345 <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">packs</span><span class="special">,</span> <span class="identifier">max_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">// P(X &gt; 3.1)</span>
 346
 347 <span class="keyword">double</span> <span class="identifier">under_weight</span> <span class="special">=</span> <span class="number">2.9</span><span class="special">;</span>
 348 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span><span class="string">"fraction of packs &lt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span> <span class="special">&lt;&lt;</span> <span class="string">" with a mean of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">mean</span>
 349   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">packs</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 350 <span class="comment">// fraction of packs &lt;= 2.9 with a mean of 3 is 0.841345</span>
 351 <span class="comment">// This is 0.84 - more than the target 0.95</span>
 352 <span class="comment">// Want 95% to be over this weight, so what should we set the mean weight to be?</span>
 353 <span class="comment">// KK StatCalc says:</span>
 354 <span class="keyword">double</span> <span class="identifier">over_mean</span> <span class="special">=</span> <span class="number">3.0664</span><span class="special">;</span>
 355 <span class="identifier">normal</span> <span class="identifier">xpacks</span><span class="special">(</span><span class="identifier">over_mean</span><span class="special">,</span> <span class="identifier">standard_deviation</span><span class="special">);</span>
 356 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"fraction of packs &gt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span>
 357 <span class="special">&lt;&lt;</span> <span class="string">" with a mean of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">xpacks</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span>
 358   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">xpacks</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 359 <span class="comment">// fraction of packs &gt;= 2.9 with a mean of 3.06449 is 0.950005</span>
 360 <span class="keyword">double</span> <span class="identifier">under_fraction</span> <span class="special">=</span> <span class="number">0.05</span><span class="special">;</span>  <span class="comment">// so 95% are above the minimum weight mean - sd = 2.9</span>
 361 <span class="keyword">double</span> <span class="identifier">low_limit</span> <span class="special">=</span> <span class="identifier">standard_deviation</span><span class="special">;</span>
 362 <span class="keyword">double</span> <span class="identifier">offset</span> <span class="special">=</span> <span class="identifier">mean</span> <span class="special">-</span> <span class="identifier">low_limit</span> <span class="special">-</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">packs</span><span class="special">,</span> <span class="identifier">under_fraction</span><span class="special">);</span>
 363 <span class="keyword">double</span> <span class="identifier">nominal_mean</span> <span class="special">=</span> <span class="identifier">mean</span> <span class="special">+</span> <span class="identifier">offset</span><span class="special">;</span>
 364
 365 <span class="identifier">normal</span> <span class="identifier">nominal_packs</span><span class="special">(</span><span class="identifier">nominal_mean</span><span class="special">,</span> <span class="identifier">standard_deviation</span><span class="special">);</span>
 366 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Setting the packer to "</span> <span class="special">&lt;&lt;</span> <span class="identifier">nominal_mean</span> <span class="special">&lt;&lt;</span> <span class="string">" will mean that "</span>
 367   <span class="special">&lt;&lt;</span> <span class="string">"fraction of packs &gt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span>
 368   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">nominal_packs</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 369 </pre>
 370 <p>
 371             Setting the packer to 3.06449 will mean that fraction of packs &gt;=
 372             2.9 is 0.95.
 373           </p>
 374 <p>
 375             Setting the packer to 3.13263 will mean that fraction of packs &gt;=
 376             2.9 is 0.99, but will more than double the mean loss from 0.0644 to 0.133.
 377           </p>
 378 <p>
 379             Alternatively, we could invest in a better (more precise) packer with
 380             a lower standard deviation.
 381           </p>
 382 <p>
 383             To estimate how much better (how much smaller standard deviation) it
 384             would have to be, we need to get the 5% quantile to be located at the
 385             under_weight limit, 2.9
 386           </p>
 387 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">p</span> <span class="special">=</span> <span class="number">0.05</span><span class="special">;</span> <span class="comment">// wanted p th quantile.</span>
 388 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Quantile of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">p</span> <span class="special">&lt;&lt;</span> <span class="string">" = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">packs</span><span class="special">,</span> <span class="identifier">p</span><span class="special">)</span>
 389   <span class="special">&lt;&lt;</span> <span class="string">", mean = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">packs</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="string">", sd = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">packs</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">//</span>
 390 </pre>
 391 <p>
 392             Quantile of 0.05 = 2.83551, mean = 3, sd = 0.1
 393           </p>
 394 <p>
 395             With the current packer (mean = 3, sd = 0.1), the 5% quantile is at 2.8551
 396             kg, a little below our target of 2.9 kg. So we know that the standard
 397             deviation is going to have to be smaller.
 398           </p>
 399 <p>
 400             Let's start by guessing that it (now 0.1) needs to be halved, to a standard
 401             deviation of 0.05
 402           </p>
 403 <pre class="programlisting"><span class="identifier">normal</span> <span class="identifier">pack05</span><span class="special">(</span><span class="identifier">mean</span><span class="special">,</span> <span class="number">0.05</span><span class="special">);</span>
 404 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Quantile of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">p</span> <span class="special">&lt;&lt;</span> <span class="string">" = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">pack05</span><span class="special">,</span> <span class="identifier">p</span><span class="special">)</span>
 405   <span class="special">&lt;&lt;</span> <span class="string">", mean = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack05</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="string">", sd = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack05</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 406
 407 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span><span class="string">"Fraction of packs &gt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span> <span class="special">&lt;&lt;</span> <span class="string">" with a mean of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">mean</span>
 408   <span class="special">&lt;&lt;</span> <span class="string">" and standard deviation of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack05</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span>
 409   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">pack05</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 410 <span class="comment">//</span>
 411 </pre>
 412 <p>
 413             Fraction of packs &gt;= 2.9 with a mean of 3 and standard deviation of
 414             0.05 is 0.9772
 415           </p>
 416 <p>
 417             So 0.05 was quite a good guess, but we are a little over the 2.9 target,
 418             so the standard deviation could be a tiny bit more. So we could do some
 419             more guessing to get closer, say by increasing to 0.06
 420           </p>
 421 <pre class="programlisting"><span class="identifier">normal</span> <span class="identifier">pack06</span><span class="special">(</span><span class="identifier">mean</span><span class="special">,</span> <span class="number">0.06</span><span class="special">);</span>
 422 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Quantile of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">p</span> <span class="special">&lt;&lt;</span> <span class="string">" = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">pack06</span><span class="special">,</span> <span class="identifier">p</span><span class="special">)</span>
 423   <span class="special">&lt;&lt;</span> <span class="string">", mean = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack06</span><span class="special">.</span><span class="identifier">mean</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="string">", sd = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack06</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 424
 425 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span><span class="string">"Fraction of packs &gt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span> <span class="special">&lt;&lt;</span> <span class="string">" with a mean of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">mean</span>
 426   <span class="special">&lt;&lt;</span> <span class="string">" and standard deviation of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack06</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span>
 427   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">pack06</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 428 </pre>
 429 <p>
 430             Fraction of packs &gt;= 2.9 with a mean of 3 and standard deviation of
 431             0.06 is 0.9522
 432           </p>
 433 <p>
 434             Now we are getting really close, but to do the job properly, we could
 435             use root finding method, for example the tools provided, and used elsewhere,
 436             in the Math Toolkit, see <a class="link" href="../../../roots_noderiv.html" title="Root Finding Without Derivatives">root-finding
 437             without derivatives</a>.
 438           </p>
 439 <p>
 440             But in this normal distribution case, we could be even smarter and make
 441             a direct calculation.
 442           </p>
 443 <pre class="programlisting"><span class="identifier">normal</span> <span class="identifier">s</span><span class="special">;</span> <span class="comment">// For standard normal distribution,</span>
 444 <span class="keyword">double</span> <span class="identifier">sd</span> <span class="special">=</span> <span class="number">0.1</span><span class="special">;</span>
 445 <span class="keyword">double</span> <span class="identifier">x</span> <span class="special">=</span> <span class="number">2.9</span><span class="special">;</span> <span class="comment">// Our required limit.</span>
 446 <span class="comment">// then probability p = N((x - mean) / sd)</span>
 447 <span class="comment">// So if we want to find the standard deviation that would be required to meet this limit,</span>
 448 <span class="comment">// so that the p th quantile is located at x,</span>
 449 <span class="comment">// in this case the 0.95 (95%) quantile at 2.9 kg pack weight, when the mean is 3 kg.</span>
 450
 451 <span class="keyword">double</span> <span class="identifier">prob</span> <span class="special">=</span>  <span class="identifier">pdf</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="special">(</span><span class="identifier">x</span> <span class="special">-</span> <span class="identifier">mean</span><span class="special">)</span> <span class="special">/</span> <span class="identifier">sd</span><span class="special">);</span>
 452 <span class="keyword">double</span> <span class="identifier">qp</span> <span class="special">=</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">s</span><span class="special">,</span> <span class="number">0.95</span><span class="special">);</span>
 453 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"prob = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">prob</span> <span class="special">&lt;&lt;</span> <span class="string">", quantile(p) "</span> <span class="special">&lt;&lt;</span> <span class="identifier">qp</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span> <span class="comment">// p = 0.241971, quantile(p) 1.64485</span>
 454 <span class="comment">// Rearranging, we can directly calculate the required standard deviation:</span>
 455 <span class="keyword">double</span> <span class="identifier">sd95</span> <span class="special">=</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">abs</span><span class="special">((</span><span class="identifier">x</span> <span class="special">-</span> <span class="identifier">mean</span><span class="special">))</span> <span class="special">/</span> <span class="identifier">qp</span><span class="special">;</span>
 456
 457 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"If we want the "</span><span class="special">&lt;&lt;</span> <span class="identifier">p</span> <span class="special">&lt;&lt;</span> <span class="string">" th quantile to be located at "</span>
 458   <span class="special">&lt;&lt;</span> <span class="identifier">x</span> <span class="special">&lt;&lt;</span> <span class="string">", would need a standard deviation of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">sd95</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 459
 460 <span class="identifier">normal</span> <span class="identifier">pack95</span><span class="special">(</span><span class="identifier">mean</span><span class="special">,</span> <span class="identifier">sd95</span><span class="special">);</span>  <span class="comment">// Distribution of the 'ideal better' packer.</span>
 461 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span><span class="string">"Fraction of packs &gt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">under_weight</span> <span class="special">&lt;&lt;</span> <span class="string">" with a mean of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">mean</span>
 462   <span class="special">&lt;&lt;</span> <span class="string">" and standard deviation of "</span> <span class="special">&lt;&lt;</span> <span class="identifier">pack95</span><span class="special">.</span><span class="identifier">standard_deviation</span><span class="special">()</span>
 463   <span class="special">&lt;&lt;</span> <span class="string">" is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">pack95</span><span class="special">,</span> <span class="identifier">under_weight</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 464
 465 <span class="comment">// Fraction of packs &gt;= 2.9 with a mean of 3 and standard deviation of 0.0608 is 0.95</span>
 466 </pre>
 467 <p>
 468             Notice that these two deceptively simple questions (do we over-fill or
 469             measure better) are actually very common. The weight of beef might be
 470             replaced by a measurement of more or less anything. But the calculations
 471             rely on the accuracy of the standard deviation - something that is almost
 472             always less good than we might wish, especially if based on a few measurements.
 473           </p>
 474 <h5>
 475 <a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.h6"></a>
 476             <span class="phrase"><a name="math_toolkit.stat_tut.weg.normal_example.normal_misc.length_of_bolts"></a></span><a class="link" href="normal_misc.html#math_toolkit.stat_tut.weg.normal_example.normal_misc.length_of_bolts">Length
 477             of bolts</a>
 478           </h5>
 479 <p>
 480             A bolt is usable if between 3.9 and 4.1 long. From a large batch of bolts,
 481             a sample of 50 show a mean length of 3.95 with standard deviation 0.1.
 482             Assuming a normal distribution, what proportion is usable? The true sample
 483             mean is unknown, but we can use the sample mean and standard deviation
 484             to find approximate solutions.
 485           </p>
 486 <pre class="programlisting">    <span class="identifier">normal</span> <span class="identifier">bolts</span><span class="special">(</span><span class="number">3.95</span><span class="special">,</span> <span class="number">0.1</span><span class="special">);</span>
 487     <span class="keyword">double</span> <span class="identifier">top</span> <span class="special">=</span> <span class="number">4.1</span><span class="special">;</span>
 488     <span class="keyword">double</span> <span class="identifier">bottom</span> <span class="special">=</span> <span class="number">3.9</span><span class="special">;</span>
 489
 490 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction long enough [ P(X &lt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">top</span> <span class="special">&lt;&lt;</span> <span class="string">") ] is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="identifier">top</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 491 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction too short [ P(X &lt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">bottom</span> <span class="special">&lt;&lt;</span> <span class="string">") ] is "</span> <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="identifier">bottom</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 492 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction OK  -between "</span> <span class="special">&lt;&lt;</span> <span class="identifier">bottom</span> <span class="special">&lt;&lt;</span> <span class="string">" and "</span> <span class="special">&lt;&lt;</span> <span class="identifier">top</span>
 493   <span class="special">&lt;&lt;</span> <span class="string">"[ P(X &lt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">top</span>  <span class="special">&lt;&lt;</span> <span class="string">") - P(X&lt;= "</span> <span class="special">&lt;&lt;</span> <span class="identifier">bottom</span> <span class="special">&lt;&lt;</span> <span class="string">" ) ] is "</span>
 494   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="identifier">top</span><span class="special">)</span> <span class="special">-</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="identifier">bottom</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 495
 496 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction too long [ P(X &gt; "</span> <span class="special">&lt;&lt;</span> <span class="identifier">top</span> <span class="special">&lt;&lt;</span> <span class="string">") ] is "</span>
 497   <span class="special">&lt;&lt;</span> <span class="identifier">cdf</span><span class="special">(</span><span class="identifier">complement</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="identifier">top</span><span class="special">))</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 498
 499 <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"95% of bolts are shorter than "</span> <span class="special">&lt;&lt;</span> <span class="identifier">quantile</span><span class="special">(</span><span class="identifier">bolts</span><span class="special">,</span> <span class="number">0.95</span><span class="special">)</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span>
 500 </pre>
 501 </div>
 502 <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
 503 <td align="left"></td>
 504 <td align="right"><div class="copyright-footer">Copyright &#169; 2006-2019 Nikhar
 505       Agrawal, Anton Bikineev, Paul A. Bristow, Marco Guazzone, Christopher Kormanyos,
 506       Hubert Holin, Bruno Lalande, John Maddock, Jeremy Murphy, Matthew Pulver, Johan
 507       R&#229;de, Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg,
 508       Daryle Walker and Xiaogang Zhang<p>
 509         Distributed under the Boost Software License, Version 1.0. (See accompanying
 510         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>)
 511       </p>
 512 </div></td>
 513 </tr></table>
 514 <hr>
 515 <div class="spirit-nav">
 516 <a accesskey="p" href="../normal_example.html"><img src="../../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../normal_example.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="../inverse_chi_squared_eg.html"><img src="../../../../../../../../doc/src/images/next.png" alt="Next"></a>
 517 </div>
 518 </body>
 519 </html>