/Src/Dependencies/Boost/libs/math/doc/sf_and_dist/html/math_toolkit/dist/stat_tut/weg/normal_example/normal_misc.html

http://hadesmem.googlecode.com/ · HTML · 616 lines · 579 code · 37 blank · 0 comment · 0 complexity · e819a13daba7bffa7ee9f1b34c6555c1 MD5 · raw file

Large files are truncated click here to view the full file

  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="../../../../../../../../../../doc/src/boostbook.css" type="text/css">
  6. <meta name="generator" content="DocBook XSL Stylesheets V1.74.0">
  7. <link rel="home" href="../../../../../index.html" title="Math Toolkit">
  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="../nccs_eg.html" title="Non Central Chi Squared 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="../nccs_eg.html"><img src="../../../../../../../../../../doc/src/images/next.png" alt="Next"></a>
  24. </div>
  25. <div class="section" lang="en">
  26. <div class="titlepage"><div><div><h6 class="title">
  27. <a name="math_toolkit.dist.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. </h6></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. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.traditional_tables"></a><h5>
  35. <a name="id1125992"></a>
  36. <a class="link" href="normal_misc.html#math_toolkit.dist.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
  41. some std output of course).
  42. </p>
  43. <p>
  44. </p>
  45. <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
  46. </span> <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.
  47. </span>
  48. <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iostream</span><span class="special">&gt;</span>
  49. <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>
  50. <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">iomanip</span><span class="special">&gt;</span>
  51. <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>
  52. <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">limits</span><span class="special">&gt;</span>
  53. <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special">;</span>
  54. <span class="keyword">int</span> <span class="identifier">main</span><span class="special">()</span>
  55. <span class="special">{</span>
  56. <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Example: Normal distribution, Miscellaneous Applications."</span><span class="special">;</span>
  57. <span class="keyword">try</span>
  58. <span class="special">{</span>
  59. <span class="special">{</span> <span class="comment">// Traditional tables and values.</span></pre>
  60. <p>
  61. </p>
  62. <p>
  63. Let's start by printing some traditional tables.
  64. </p>
  65. <p>
  66. </p>
  67. <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
  68. </span><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.
  69. </span><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.
  70. </span><span class="comment">// but std::numeric_limits&lt;double&gt;::max_digits10; on new Standard Libraries gives
  71. </span><span class="comment">// 17, the maximum number of digits that can possibly be significant.
  72. </span><span class="comment">// std::numeric_limits&lt;double&gt;::digits10; == 15 is number of guaranteed digits,
  73. </span><span class="comment">// the other two digits being 'noisy'.
  74. </span>
  75. <span class="comment">// Construct a standard normal distribution s
  76. </span> <span class="identifier">normal</span> <span class="identifier">s</span><span class="special">;</span> <span class="comment">// (default mean = zero, and standard deviation = unity)
  77. </span> <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>
  78. <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></pre>
  79. <p>
  80. </p>
  81. <p>
  82. First the probability distribution function (pdf).
  83. </p>
  84. <p>
  85. </p>
  86. <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>
  87. <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>
  88. <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>
  89. <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>
  90. <span class="special">{</span>
  91. <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>
  92. <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>
  93. <span class="special">}</span>
  94. <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></pre>
  95. <p>
  96. </p>
  97. <p>
  98. And the area under the normal curve from -&#8734; up to z, the cumulative distribution
  99. function (cdf).
  100. </p>
  101. <p>
  102. </p>
  103. <pre class="programlisting"><span class="comment">// For a standard normal distribution
  104. </span><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>
  105. <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>
  106. <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>
  107. <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>
  108. <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>
  109. <span class="special">{</span>
  110. <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>
  111. <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>
  112. <span class="special">}</span>
  113. <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></pre>
  114. <p>
  115. </p>
  116. <p>
  117. And all this you can do with a nanoscopic amount of work compared to
  118. the team of <span class="bold"><strong>human computers</strong></span> toiling
  119. with Milton Abramovitz and Irene Stegen at the US National Bureau of
  120. Standards (now <a href="http://www.nist.gov" target="_top">NIST</a>). Starting
  121. in 1938, their "Handbook of Mathematical Functions with Formulas,
  122. Graphs and Mathematical Tables", was eventually published in 1964,
  123. and has been reprinted numerous times since. (A major replacement is
  124. planned at <a href="http://dlmf.nist.gov" target="_top">Digital Library of Mathematical
  125. Functions</a>).
  126. </p>
  127. <p>
  128. Pretty-printing a traditional 2-dimensional table is left as an exercise
  129. for the student, but why bother now that the Math Toolkit lets you
  130. write
  131. </p>
  132. <p>
  133. </p>
  134. <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>
  135. <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></pre>
  136. <p>
  137. </p>
  138. <p>
  139. Correspondingly, we can obtain the traditional 'critical' values for
  140. significance levels. For the 95% confidence level, the significance
  141. level usually called alpha, is 0.05 = 1 - 0.95 (for a one-sided test),
  142. so we can write
  143. </p>
  144. <p>
  145. </p>
  146. <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>
  147. <span class="comment">// 95% of area has a z below 1.64485</span></pre>
  148. <p>
  149. </p>
  150. <p>
  151. and a two-sided test (a comparison between two levels, rather than
  152. a one-sided test)
  153. </p>
  154. <p>
  155. </p>
  156. <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>
  157. <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>
  158. <span class="comment">// 95% of area has a z between 1.95996 and -1.95996</span></pre>
  159. <p>
  160. </p>
  161. <p>
  162. First, define a table of significance levels: these are the probabilities
  163. that the true occurrence frequency lies outside the calculated interval.
  164. </p>
  165. <p>
  166. It is convenient to have an alpha level for the probability that z
  167. lies outside just one standard deviation. This will not be some nice
  168. neat number like 0.05, but we can easily calculate it,
  169. </p>
  170. <p>
  171. </p>
  172. <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
  173. </span><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></pre>
  174. <p>
  175. </p>
  176. <p>
  177. and place in our array of favorite alpha values.
  178. </p>
  179. <p>
  180. </p>
  181. <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.
  182. </span> <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></pre>
  183. <p>
  184. </p>
  185. <p>
  186. Confidence value as % is (1 - alpha) * 100 (so alpha 0.05 == 95% confidence)
  187. that the true occurrence frequency lies <span class="bold"><strong>inside</strong></span>
  188. the calculated interval.
  189. </p>
  190. <p>
  191. </p>
  192. <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>
  193. <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>
  194. <span class="keyword">for</span> <span class="special">(</span><span class="keyword">int</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>
  195. <span class="special">{</span>
  196. <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>
  197. <span class="comment">// Use quantile(complement(s, alpha[i]/2)) to avoid potential loss of accuracy from quantile(s, 1 - alpha[i]/2)
  198. </span><span class="special">}</span>
  199. <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">endl</span><span class="special">;</span></pre>
  200. <p>
  201. </p>
  202. <p>
  203. Notice the distinction between one-sided (also called one-tailed) where
  204. we are using a &gt; <span class="bold"><strong>or</strong></span> &lt; test (and
  205. not both) and considering the area of the tail (integral) from z up
  206. to +&#8734;, and a two-sided test where we are using two &gt; <span class="bold"><strong>and</strong></span>
  207. &lt; tests, and thus considering two tails, from -&#8734; up to z low and z
  208. high up to +&#8734;.
  209. </p>
  210. <p>
  211. So the 2-sided values alpha[i] are calculated using alpha[i]/2.
  212. </p>
  213. <p>
  214. If we consider a simple example of alpha = 0.05, then for a two-sided
  215. test, the lower tail area from -&#8734; up to -1.96 is 0.025 (alpha/2) and
  216. the upper tail area from +z up to +1.96 is also 0.025 (alpha/2), and
  217. the area between -1.96 up to 12.96 is alpha = 0.95. and the sum of
  218. the two tails is 0.025 + 0.025 = 0.05,
  219. </p>
  220. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.standard_deviations_either_side_of_the_mean"></a><h5>
  221. <a name="id1129710"></a>
  222. <a class="link" href="normal_misc.html#math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.standard_deviations_either_side_of_the_mean">Standard
  223. deviations either side of the Mean</a>
  224. </h5>
  225. <p>
  226. Armed with the cumulative distribution function, we can easily calculate
  227. the easy to remember proportion of values that lie within 1, 2 and
  228. 3 standard deviations from the mean.
  229. </p>
  230. <p>
  231. </p>
  232. <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>
  233. <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>
  234. <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
  235. </span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"cdf(complement(s, s.standard_deviation())) = "</span>
  236. <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>
  237. <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>
  238. <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>
  239. <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>
  240. <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>
  241. <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>
  242. <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></pre>
  243. <p>
  244. </p>
  245. <p>
  246. To a useful precision, the 1, 2 &amp; 3 percentages are 68, 95 and
  247. 99.7, and these are worth memorising as useful 'rules of thumb', as,
  248. for example, in <a href="http://en.wikipedia.org/wiki/Standard_deviation" target="_top">standard
  249. deviation</a>:
  250. </p>
  251. <pre class="programlisting">Fraction 1 standard deviation within either side of mean is 0.683
  252. Fraction 2 standard deviations within either side of mean is 0.954
  253. Fraction 3 standard deviations within either side of mean is 0.997
  254. </pre>
  255. <p>
  256. We could of course get some really accurate values for these <a href="http://en.wikipedia.org/wiki/Confidence_interval" target="_top">confidence intervals</a>
  257. by using cout.precision(15);
  258. </p>
  259. <pre class="programlisting">Fraction 1 standard deviation within either side of mean is 0.682689492137086
  260. Fraction 2 standard deviations within either side of mean is 0.954499736103642
  261. Fraction 3 standard deviations within either side of mean is 0.997300203936740
  262. </pre>
  263. <p>
  264. But before you get too excited about this impressive precision, don't
  265. forget that the <span class="bold"><strong>confidence intervals of the standard
  266. deviation</strong></span> are surprisingly wide, especially if you have
  267. estimated the standard deviation from only a few measurements.
  268. </p>
  269. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.some_simple_examples"></a><h5>
  270. <a name="id1130231"></a>
  271. <a class="link" href="normal_misc.html#math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.some_simple_examples">Some
  272. simple examples</a>
  273. </h5>
  274. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.life_of_light_bulbs"></a><h5>
  275. <a name="id1130244"></a>
  276. <a class="link" href="normal_misc.html#math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.life_of_light_bulbs">Life
  277. of light bulbs</a>
  278. </h5>
  279. <p>
  280. Examples from K. Krishnamoorthy, Handbook of Statistical Distributions
  281. with Applications, ISBN 1 58488 635 8, page 125... implemented using
  282. the Math Toolkit library.
  283. </p>
  284. <p>
  285. A few very simple examples are shown here:
  286. </p>
  287. <p>
  288. </p>
  289. <pre class="programlisting"><span class="comment">// K. Krishnamoorthy, Handbook of Statistical Distributions with Applications,
  290. </span> <span class="comment">// ISBN 1 58488 635 8, page 125, example 10.3.5</span></pre>
  291. <p>
  292. </p>
  293. <p>
  294. Mean lifespan of 100 W bulbs is 1100 h with standard deviation of 100
  295. h. Assuming, perhaps with little evidence and much faith, that the
  296. distribution is normal, we construct a normal distribution called
  297. <span class="emphasis"><em>bulbs</em></span> with these values:
  298. </p>
  299. <p>
  300. </p>
  301. <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>
  302. <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>
  303. <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>
  304. <span class="keyword">double</span> <span class="identifier">expected_life</span> <span class="special">=</span> <span class="number">1000.</span><span class="special">;</span></pre>
  305. <p>
  306. </p>
  307. <p>
  308. The we can use the Cumulative distribution function to predict fractions
  309. (or percentages, if * 100) that will last various lifetimes.
  310. </p>
  311. <p>
  312. </p>
  313. <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)
  314. </span> <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>
  315. <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)
  316. </span> <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>
  317. <span class="keyword">double</span> <span class="identifier">min_life</span> <span class="special">=</span> <span class="number">900</span><span class="special">;</span>
  318. <span class="keyword">double</span> <span class="identifier">max_life</span> <span class="special">=</span> <span class="number">1200</span><span class="special">;</span>
  319. <span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Fraction of bulbs that will last between "</span>
  320. <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>
  321. <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)
  322. </span> <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></pre>
  323. <p>
  324. </p>
  325. <div class="note"><table border="0" summary="Note">
  326. <tr>
  327. <td rowspan="2" align="center" valign="top" width="25"><img alt="[Note]" src="../../../../../../../../../../doc/src/images/note.png"></td>
  328. <th align="left">Note</th>
  329. </tr>
  330. <tr><td align="left" valign="top"><p>
  331. Real-life failures are often very ab-normal, with a significant number
  332. that 'dead-on-arrival' or suffer failure very early in their life:
  333. the lifetime of the survivors of 'early mortality' may be well described
  334. by the normal distribution.
  335. </p></td></tr>
  336. </table></div>
  337. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.how_many_onions_"></a><h5>
  338. <a name="id1130752"></a>
  339. <a class="link" href="normal_misc.html#math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.how_many_onions_">How
  340. many onions?</a>
  341. </h5>
  342. <p>
  343. Weekly demand for 5 lb sacks of onions at a store is normally distributed
  344. with mean 140 sacks and standard deviation 10.
  345. </p>
  346. <p>
  347. </p>
  348. <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.
  349. </span><span class="keyword">double</span> <span class="identifier">standard_deviation</span> <span class="special">=</span> <span class="number">10</span><span class="special">;</span>
  350. <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>
  351. <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.
  352. </span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Percentage of weeks overstocked "</span>
  353. <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)
  354. </span><span class="comment">// Percentage of weeks overstocked 97.7</span></pre>
  355. <p>
  356. </p>
  357. <p>
  358. So there will be lots of mouldy onions! So we should be able to say
  359. what stock level will meet demand 95% of the weeks.
  360. </p>
  361. <p>
  362. </p>
  363. <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>
  364. <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></pre>
  365. <p>
  366. </p>
  367. <p>
  368. And it is easy to estimate how to meet 80% of demand, and waste even
  369. less.
  370. </p>
  371. <p>
  372. </p>
  373. <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>
  374. <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>
  375. </pre>
  376. <p>
  377. </p>
  378. <a name="math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.packing_beef"></a><h5>
  379. <a name="id1131883"></a>
  380. <a class="link" href="normal_misc.html#math_toolkit.dist.stat_tut.weg.normal_example.normal_misc.packing_beef">Packing
  381. beef</a>
  382. </h5>
  383. <p>
  384. A machine is set to pack 3 kg of ground beef per pack. Over a long
  385. period of time it is found that the average packed was 3 kg with a
  386. standard deviation of 0.1 kg. Assuming the packing is normally distributed,
  387. we can find the fraction (or %) of packages that weigh more than 3.1
  388. kg.
  389. </p>
  390. <p>
  391. </p>
  392. <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
  393. </span><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
  394. </span><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>
  395. <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
  396. </span><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>
  397. <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)
  398. </span>
  399. <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>
  400. <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>
  401. <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>
  402. <span class="comment">// fraction of packs &lt;= 2.9 with a mean of 3 is 0.841345
  403. </span><span class="comment">// This is 0.84 - more than the target 0.95
  404. </span><span class="comment">// Want 95% to be over this weight, so what should we set the mean weight to be?
  405. </span><span class="comment">// KK StatCalc says:
  406. </span><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>
  407. <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>
  408. <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>
  409. <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>
  410. <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>
  411. <span class="comment">// fraction of packs &gt;= 2.9 with a mean of 3.06449 is 0.950005
  412. </span><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
  413. </span><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>
  414. <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>
  415. <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>
  416. <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>
  417. <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>
  418. <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>
  419. <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><spa