signal_statistics.html 37 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336
  1. <html>
  2. <head>
  3. <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
  4. <title>Signal Statistics</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="../statistics.html" title="Chapter&#160;6.&#160;Statistics">
  9. <link rel="prev" href="bivariate_statistics.html" title="Bivariate Statistics">
  10. <link rel="next" href="anderson_darling.html" title="The Anderson-Darling Test">
  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="bivariate_statistics.html"><img src="../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../statistics.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="anderson_darling.html"><img src="../../../../../doc/src/images/next.png" alt="Next"></a>
  24. </div>
  25. <div class="section">
  26. <div class="titlepage"><div><div><h2 class="title" style="clear: both">
  27. <a name="math_toolkit.signal_statistics"></a><a class="link" href="signal_statistics.html" title="Signal Statistics">Signal Statistics</a>
  28. </h2></div></div></div>
  29. <h4>
  30. <a name="math_toolkit.signal_statistics.h0"></a>
  31. <span class="phrase"><a name="math_toolkit.signal_statistics.synopsis"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.synopsis">Synopsis</a>
  32. </h4>
  33. <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">statistics</span><span class="special">/</span><span class="identifier">signal_statistics</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
  34. <span class="keyword">namespace</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span> <span class="special">{</span>
  35. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  36. <span class="keyword">auto</span> <span class="identifier">absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">Container</span> <span class="special">&amp;</span> <span class="identifier">c</span><span class="special">);</span>
  37. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">ForwardIterator</span><span class="special">&gt;</span>
  38. <span class="keyword">auto</span> <span class="identifier">absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">ForwardIterator</span> <span class="identifier">first</span><span class="special">,</span> <span class="identifier">ForwardIterator</span> <span class="identifier">last</span><span class="special">);</span>
  39. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  40. <span class="keyword">auto</span> <span class="identifier">sample_absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">Container</span> <span class="special">&amp;</span> <span class="identifier">c</span><span class="special">);</span>
  41. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">ForwardIterator</span><span class="special">&gt;</span>
  42. <span class="keyword">auto</span> <span class="identifier">sample_absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">ForwardIterator</span> <span class="identifier">first</span><span class="special">,</span> <span class="identifier">ForwardIterator</span> <span class="identifier">last</span><span class="special">);</span>
  43. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  44. <span class="keyword">auto</span> <span class="identifier">hoyer_sparsity</span><span class="special">(</span><span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">c</span><span class="special">);</span>
  45. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">ForwardIterator</span><span class="special">&gt;</span>
  46. <span class="keyword">auto</span> <span class="identifier">hoyer_sparsity</span><span class="special">(</span><span class="identifier">ForwardIterator</span> <span class="identifier">first</span><span class="special">,</span> <span class="identifier">ForwardIterator</span> <span class="identifier">last</span><span class="special">);</span>
  47. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  48. <span class="keyword">auto</span> <span class="identifier">oracle_snr</span><span class="special">(</span><span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">signal</span><span class="special">,</span> <span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">noisy_signal</span><span class="special">);</span>
  49. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  50. <span class="keyword">auto</span> <span class="identifier">oracle_snr_db</span><span class="special">(</span><span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">signal</span><span class="special">,</span> <span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">noisy_signal</span><span class="special">);</span>
  51. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">ForwardIterator</span><span class="special">&gt;</span>
  52. <span class="keyword">auto</span> <span class="identifier">m2m4_snr_estimator</span><span class="special">(</span><span class="identifier">ForwardIterator</span> <span class="identifier">first</span><span class="special">,</span> <span class="identifier">ForwardIterator</span> <span class="identifier">last</span><span class="special">,</span> <span class="keyword">decltype</span><span class="special">(*</span><span class="identifier">first</span><span class="special">)</span> <span class="identifier">estimated_signal_kurtosis</span><span class="special">=</span><span class="number">1</span><span class="special">,</span> <span class="keyword">decltype</span><span class="special">(*</span><span class="identifier">first</span><span class="special">)</span> <span class="identifier">estimated_noise_kurtosis</span><span class="special">=</span><span class="number">3</span><span class="special">);</span>
  53. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  54. <span class="keyword">auto</span> <span class="identifier">m2m4_snr_estimator</span><span class="special">(</span><span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">noisy_signal</span><span class="special">,</span> <span class="keyword">typename</span> <span class="identifier">Container</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">estimated_signal_kurtosis</span><span class="special">=</span><span class="number">1</span><span class="special">,</span> <span class="keyword">typename</span> <span class="identifier">Container</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">estimate_noise_kurtosis</span><span class="special">=</span><span class="number">3</span><span class="special">);</span>
  55. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">ForwardIterator</span><span class="special">&gt;</span>
  56. <span class="keyword">auto</span> <span class="identifier">m2m4_snr_estimator_db</span><span class="special">(</span><span class="identifier">ForwardIterator</span> <span class="identifier">first</span><span class="special">,</span> <span class="identifier">ForwardIterator</span> <span class="identifier">last</span><span class="special">,</span> <span class="keyword">decltype</span><span class="special">(*</span><span class="identifier">first</span><span class="special">)</span> <span class="identifier">estimated_signal_kurtosis</span><span class="special">=</span><span class="number">1</span><span class="special">,</span> <span class="keyword">decltype</span><span class="special">(*</span><span class="identifier">first</span><span class="special">)</span> <span class="identifier">estimated_noise_kurtosis</span><span class="special">=</span><span class="number">3</span><span class="special">);</span>
  57. <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Container</span><span class="special">&gt;</span>
  58. <span class="keyword">auto</span> <span class="identifier">m2m4_snr_estimator_db</span><span class="special">(</span><span class="identifier">Container</span> <span class="keyword">const</span> <span class="special">&amp;</span> <span class="identifier">noisy_signal</span><span class="special">,</span><span class="keyword">typename</span> <span class="identifier">Container</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">estimated_signal_kurtosis</span><span class="special">=</span><span class="number">1</span><span class="special">,</span> <span class="keyword">typename</span> <span class="identifier">Container</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">estimate_noise_kurtosis</span><span class="special">=</span><span class="number">3</span><span class="special">);</span>
  59. <span class="special">}</span>
  60. </pre>
  61. <h4>
  62. <a name="math_toolkit.signal_statistics.h1"></a>
  63. <span class="phrase"><a name="math_toolkit.signal_statistics.description"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.description">Description</a>
  64. </h4>
  65. <p>
  66. The file <code class="computeroutput"><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">statistics</span><span class="special">/</span><span class="identifier">signal_statistics</span><span class="special">.</span><span class="identifier">hpp</span></code> is a
  67. set of facilities for computing quantities commonly used in signal analysis.
  68. </p>
  69. <p>
  70. Our examples use <code class="computeroutput"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span></code> to
  71. hold the data, but this not required. In general, you can store your data in
  72. an Eigen array, and Armadillo vector, <code class="computeroutput"><span class="identifier">std</span><span class="special">::</span><span class="identifier">array</span></code>,
  73. and for many of the routines, a <code class="computeroutput"><span class="identifier">std</span><span class="special">::</span><span class="identifier">forward_list</span></code>.
  74. These routines are usable in float, double, long double, and Boost.Multiprecision
  75. precision, as well as their complex extensions whenever the computation is
  76. well-defined.
  77. </p>
  78. <h4>
  79. <a name="math_toolkit.signal_statistics.h2"></a>
  80. <span class="phrase"><a name="math_toolkit.signal_statistics.absolute_gini_coefficient"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.absolute_gini_coefficient">Absolute
  81. Gini Coefficient</a>
  82. </h4>
  83. <p>
  84. The Gini coefficient, first used to measure wealth inequality, is also one
  85. of the best measures of the sparsity of an expansion in a basis. A sparse expansion
  86. has most of its norm concentrated in just a few coefficients, making the connection
  87. with wealth inequality obvious. See <a href="https://arxiv.org/pdf/0811.4706.pdf" target="_top">Hurley
  88. and Rickard</a> for details. However, for measuring sparsity, the phase
  89. of the numbers is irrelevant, so we provide the <code class="computeroutput"><span class="identifier">absolute_gini_coefficient</span></code>:
  90. </p>
  91. <pre class="programlisting"><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">statistics</span><span class="special">::</span><span class="identifier">sample_absolute_gini_coefficient</span><span class="special">;</span>
  92. <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">statistics</span><span class="special">::</span><span class="identifier">absolute_gini_coefficient</span><span class="special">;</span>
  93. <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;&gt;</span> <span class="identifier">v</span><span class="special">{{</span><span class="number">0</span><span class="special">,</span><span class="number">1</span><span class="special">},</span> <span class="special">{</span><span class="number">0</span><span class="special">,</span><span class="number">0</span><span class="special">},</span> <span class="special">{</span><span class="number">0</span><span class="special">,</span><span class="number">0</span><span class="special">},</span> <span class="special">{</span><span class="number">0</span><span class="special">,</span><span class="number">0</span><span class="special">}};</span>
  94. <span class="keyword">double</span> <span class="identifier">abs_gini</span> <span class="special">=</span> <span class="identifier">sample_absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">v</span><span class="special">);</span>
  95. <span class="comment">// now abs_gini = 1; maximally unequal</span>
  96. <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;&gt;</span> <span class="identifier">w</span><span class="special">{{</span><span class="number">0</span><span class="special">,</span><span class="number">1</span><span class="special">},</span> <span class="special">{</span><span class="number">1</span><span class="special">,</span><span class="number">0</span><span class="special">},</span> <span class="special">{</span><span class="number">0</span><span class="special">,-</span><span class="number">1</span><span class="special">},</span> <span class="special">{-</span><span class="number">1</span><span class="special">,</span><span class="number">0</span><span class="special">}};</span>
  97. <span class="identifier">abs_gini</span> <span class="special">=</span> <span class="identifier">absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">w</span><span class="special">);</span>
  98. <span class="comment">// now abs_gini = 0; every element of the vector has equal magnitude</span>
  99. <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">u</span><span class="special">{-</span><span class="number">1</span><span class="special">,</span> <span class="number">1</span><span class="special">,</span> <span class="special">-</span><span class="number">1</span><span class="special">};</span>
  100. <span class="identifier">abs_gini</span> <span class="special">=</span> <span class="identifier">absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">u</span><span class="special">);</span>
  101. <span class="comment">// now abs_gini = 0</span>
  102. <span class="comment">// Alternative call useful for computing over subset of the input:</span>
  103. <span class="identifier">abs_gini</span> <span class="special">=</span> <span class="identifier">absolute_gini_coefficient</span><span class="special">(</span><span class="identifier">u</span><span class="special">.</span><span class="identifier">begin</span><span class="special">(),</span> <span class="identifier">u</span><span class="special">.</span><span class="identifier">begin</span><span class="special">()</span> <span class="special">+</span> <span class="number">1</span><span class="special">);</span>
  104. </pre>
  105. <p>
  106. The sample Gini coefficient returns unity for a vector which has only one nonzero
  107. coefficient. The population Gini coefficient of a vector with one non-zero
  108. element is dependent on the length of the input.
  109. </p>
  110. <p>
  111. The sample Gini coefficient lacks one desirable property of the population
  112. Gini coefficient, namely that "cloning" a vector has the same Gini
  113. coefficient; though cloning holds to very high accuracy with the sample Gini
  114. coefficient and can easily be recovered by a rescaling.
  115. </p>
  116. <p>
  117. If sorting the input data is too much expense for a sparsity measure (is it
  118. going to be perfect anyway?), consider calculating the Hoyer sparsity instead.
  119. </p>
  120. <h4>
  121. <a name="math_toolkit.signal_statistics.h3"></a>
  122. <span class="phrase"><a name="math_toolkit.signal_statistics.hoyer_sparsity"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.hoyer_sparsity">Hoyer
  123. Sparsity</a>
  124. </h4>
  125. <p>
  126. The Hoyer sparsity measures a normalized ratio of the &#8467;<sup>1</sup> and &#8467;<sup>2</sup> norms.
  127. As the name suggests, it is used to measure the sparsity of an expansion in
  128. some basis.
  129. </p>
  130. <p>
  131. The Hoyer sparsity computes (&#8730;<span class="emphasis"><em>N</em></span> - &#8467;<sup>1</sup>(v)/&#8467;<sup>2</sup>(v))/(&#8730;N
  132. -1). For details, see <a href="http://www.jmlr.org/papers/volume5/hoyer04a/hoyer04a.pdf" target="_top">Hoyer</a>
  133. as well as <a href="https://arxiv.org/pdf/0811.4706.pdf" target="_top">Hurley and Rickard</a>.
  134. </p>
  135. <p>
  136. A few special cases will serve to clarify the intended use: If <span class="emphasis"><em>v</em></span>
  137. has only one nonzero coefficient, the Hoyer sparsity attains its maxima of
  138. 1. If the coefficients of <span class="emphasis"><em>v</em></span> all have the same magnitude,
  139. then the Hoyer sparsity attains its minima of zero. If the elements of <span class="emphasis"><em>v</em></span>
  140. are uniformly distributed on an interval [0, <span class="emphasis"><em>b</em></span>], then
  141. the Hoyer sparsity is approximately 0.133.
  142. </p>
  143. <p>
  144. Usage:
  145. </p>
  146. <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;</span> <span class="identifier">v</span><span class="special">{</span><span class="number">1</span><span class="special">,</span><span class="number">0</span><span class="special">,</span><span class="number">0</span><span class="special">};</span>
  147. <span class="identifier">Real</span> <span class="identifier">hs</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">hoyer_sparsity</span><span class="special">(</span><span class="identifier">v</span><span class="special">);</span>
  148. <span class="comment">// hs = 1</span>
  149. <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;</span> <span class="identifier">v</span><span class="special">{</span><span class="number">1</span><span class="special">,-</span><span class="number">1</span><span class="special">,</span><span class="number">1</span><span class="special">};</span>
  150. <span class="identifier">Real</span> <span class="identifier">hs</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">hoyer_sparsity</span><span class="special">(</span><span class="identifier">v</span><span class="special">.</span><span class="identifier">begin</span><span class="special">(),</span> <span class="identifier">v</span><span class="special">.</span><span class="identifier">end</span><span class="special">());</span>
  151. <span class="comment">// hs = 0</span>
  152. </pre>
  153. <p>
  154. The container must be forward iterable and the contents are not modified. Accepts
  155. real, complex, and integer inputs. If the input is an integral type, the output
  156. is a double precision float.
  157. </p>
  158. <h4>
  159. <a name="math_toolkit.signal_statistics.h4"></a>
  160. <span class="phrase"><a name="math_toolkit.signal_statistics.oracle_signal_to_noise_ratio"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.oracle_signal_to_noise_ratio">Oracle
  161. Signal-to-noise ratio</a>
  162. </h4>
  163. <p>
  164. The function <code class="computeroutput"><span class="identifier">oracle_snr</span></code> computes
  165. the ratio &#8214; <span class="emphasis"><em>s</em></span> &#8214;<sub>2</sub><sup>2</sup> / &#8214; <span class="emphasis"><em>s</em></span>
  166. - <span class="emphasis"><em>x</em></span> &#8214;<sub>2</sub><sup>2</sup>, where <span class="emphasis"><em>s</em></span> is signal
  167. and <span class="emphasis"><em>x</em></span> is a noisy signal. The function <code class="computeroutput"><span class="identifier">oracle_snr_db</span></code>
  168. computes 10<code class="computeroutput"><span class="identifier">log</span></code><sub>10</sub>(&#8214;
  169. <span class="emphasis"><em>s</em></span> &#8214;<sup>2</sup> / &#8214; <span class="emphasis"><em>s</em></span> - <span class="emphasis"><em>x</em></span>
  170. &#8214;<sup>2</sup>). The functions are so named because in general, one does not know
  171. how to decompose a real signal <span class="emphasis"><em>x</em></span> into <span class="emphasis"><em>s</em></span>
  172. + <span class="emphasis"><em>w</em></span> and as such <span class="emphasis"><em>s</em></span> is regarded as
  173. oracle information. Hence this function is mainly useful for unit testing other
  174. SNR estimators.
  175. </p>
  176. <p>
  177. Usage:
  178. </p>
  179. <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">signal</span><span class="special">(</span><span class="number">500</span><span class="special">,</span> <span class="number">3.2</span><span class="special">);</span>
  180. <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">noisy_signal</span><span class="special">(</span><span class="number">500</span><span class="special">);</span>
  181. <span class="comment">// fill 'noisy_signal' signal + noise</span>
  182. <span class="keyword">double</span> <span class="identifier">snr_db</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">oracle_snr_db</span><span class="special">(</span><span class="identifier">signal</span><span class="special">,</span> <span class="identifier">noisy_signal</span><span class="special">);</span>
  183. <span class="keyword">double</span> <span class="identifier">snr</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">oracle_snr</span><span class="special">(</span><span class="identifier">signal</span><span class="special">,</span> <span class="identifier">noisy_signal</span><span class="special">);</span>
  184. </pre>
  185. <p>
  186. The input can be real, complex, or integral. Integral inputs produce double
  187. precision floating point outputs. The input data is not modified and must satisfy
  188. the requirements of a <code class="computeroutput"><span class="identifier">RandomAccessContainer</span></code>.
  189. </p>
  190. <h4>
  191. <a name="math_toolkit.signal_statistics.h5"></a>
  192. <span class="phrase"><a name="math_toolkit.signal_statistics.m_sub_2_m_sub_4_snr_estimation"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.m_sub_2_m_sub_4_snr_estimation"><span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> SNR
  193. Estimation</a>
  194. </h4>
  195. <p>
  196. Estimates the SNR of a noisy signal via the <span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> method.
  197. See <a href="https://doi.org/10.1109/26.871393" target="_top">Pauluzzi and N.C. Beaulieu</a>
  198. and <a href="https://doi.org/10.1109/ISIT.1994.394869" target="_top">Matzner and Englberger</a>
  199. for details.
  200. </p>
  201. <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">noisy_signal</span><span class="special">(</span><span class="number">512</span><span class="special">);</span>
  202. <span class="comment">// fill noisy_signal with data contaminated by Gaussian white noise:</span>
  203. <span class="keyword">double</span> <span class="identifier">est_snr_db</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">m2m4_snr_estimator_db</span><span class="special">(</span><span class="identifier">noisy_signal</span><span class="special">);</span>
  204. </pre>
  205. <p>
  206. The <span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> SNR estimator is an "in-service"
  207. estimator, meaning that the estimate is made using the noisy, data-bearing
  208. signal, and does not require a background estimate. This estimator has been
  209. found to be work best between roughly -3 and 15db, tending to overestimate
  210. the noise below -3db, and underestimate the noise above 15db. See <a href="https://www.mdpi.com/2078-2489/8/3/75/pdf" target="_top">Xue
  211. et al</a> for details.
  212. </p>
  213. <p>
  214. The <span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> SNR estimator, by default,
  215. assumes that the kurtosis of the signal is 1 and the kurtosis of the noise
  216. is 3, the latter corresponding to Gaussian noise. These parameters, however,
  217. can be overridden:
  218. </p>
  219. <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">noisy_signal</span><span class="special">(</span><span class="number">512</span><span class="special">);</span>
  220. <span class="comment">// fill noisy_signal with the data:</span>
  221. <span class="keyword">double</span> <span class="identifier">signal_kurtosis</span> <span class="special">=</span> <span class="number">1.5</span><span class="special">;</span>
  222. <span class="comment">// Noise is assumed to follow Laplace distribution, which has kurtosis of 6:</span>
  223. <span class="keyword">double</span> <span class="identifier">noise_kurtosis</span> <span class="special">=</span> <span class="number">6</span><span class="special">;</span>
  224. <span class="keyword">double</span> <span class="identifier">est_snr</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">statistics</span><span class="special">::</span><span class="identifier">m2m4_snr_estimator_db</span><span class="special">(</span><span class="identifier">noisy_signal</span><span class="special">,</span> <span class="identifier">signal_kurtosis</span><span class="special">,</span> <span class="identifier">noise_kurtosis</span><span class="special">);</span>
  225. </pre>
  226. <p>
  227. Now, technically the method is a "blind SNR estimator", meaning that
  228. the no <span class="emphasis"><em>a-priori</em></span> information about the signal is required
  229. to use the method. However, the performance of the method is <span class="emphasis"><em>vastly</em></span>
  230. better if you can come up with a better estimate of the signal and noise kurtosis.
  231. How can we do this? Suppose we know that the SNR is much greater than 1. Then
  232. we can estimate the signal kurtosis simply by using the noisy signal kurtosis.
  233. If the SNR is much less than one, this method breaks down as the noisy signal
  234. kurtosis will tend to the noise kurtosis-though in this limit we have an excellent
  235. estimator of the noise kurtosis! In addition, if you have a model of what your
  236. signal should look like, you can precompute the signal kurtosis. For example,
  237. sinusoids have a kurtosis of 1.5. See <a href="http://www.jcomputers.us/vol8/jcp0808-21.pdf" target="_top">here</a>
  238. for a study which uses estimates of this sort to improve the performance of
  239. the <span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> estimator.
  240. </p>
  241. <p>
  242. <span class="emphasis"><em>Nota bene</em></span>: The traditional definition of SNR is <span class="emphasis"><em>not</em></span>
  243. mean invariant. By this we mean that if a constant is added to every sample
  244. of a signal, the SNR is changed. For example, adding DC bias to a signal changes
  245. its SNR. For most use cases, this is really not what you intend; for example
  246. a signal consisting of zeros plus Gaussian noise has an SNR of zero, whereas
  247. a signal with a constant DC bias and random Gaussian noise might have a very
  248. large SNR.
  249. </p>
  250. <p>
  251. The <span class="emphasis"><em>M</em></span><sub>2</sub><span class="emphasis"><em>M</em></span><sub>4</sub> SNR estimator is computed
  252. from mean-invariant quantities, and hence it should really be compared to the
  253. mean-invariant SNR.
  254. </p>
  255. <p>
  256. <span class="emphasis"><em>Nota bene</em></span>: This computation requires the solution of a
  257. system of quadratic equations involving the noise kurtosis, the signal kurtosis,
  258. and the second and fourth moments of the data. There is no guarantee that a
  259. solution of this system exists for all value of these parameters, in fact nonexistence
  260. can easily be demonstrated for certain data. If there is no solution to the
  261. system, then failure is communicated by returning NaNs. This happens distressingly
  262. often; if a user is aware of any blind SNR estimators which do not suffer from
  263. this drawback, please open a github ticket and let us know.
  264. </p>
  265. <p>
  266. The author has not managed to fully characterize the conditions under which
  267. a real solution with <span class="emphasis"><em>S &gt; 0</em></span> and <span class="emphasis"><em>N &gt;0</em></span>
  268. exists. However, a very intuitive example demonstrates why nonexistence can
  269. occur. Suppose the signal and noise kurtosis are equal. Then the method has
  270. no way to distinguish between the signal and the noise, and the solution is
  271. non-unique.
  272. </p>
  273. <h4>
  274. <a name="math_toolkit.signal_statistics.h6"></a>
  275. <span class="phrase"><a name="math_toolkit.signal_statistics.references"></a></span><a class="link" href="signal_statistics.html#math_toolkit.signal_statistics.references">References</a>
  276. </h4>
  277. <div class="itemizedlist"><ul class="itemizedlist" style="list-style-type: disc; ">
  278. <li class="listitem">
  279. Mallat, Stephane. <span class="emphasis"><em>A wavelet tour of signal processing: the sparse
  280. way.</em></span> Academic press, 2008.
  281. </li>
  282. <li class="listitem">
  283. Hurley, Niall, and Scott Rickard. <span class="emphasis"><em>Comparing measures of sparsity.</em></span>
  284. IEEE Transactions on Information Theory 55.10 (2009): 4723-4741.
  285. </li>
  286. <li class="listitem">
  287. Jensen, Arne, and Anders la Cour-Harbo. <span class="emphasis"><em>Ripples in mathematics:
  288. the discrete wavelet transform.</em></span> Springer Science &amp; Business
  289. Media, 2001.
  290. </li>
  291. <li class="listitem">
  292. D. R. Pauluzzi and N. C. Beaulieu, <span class="emphasis"><em>A comparison of SNR estimation
  293. techniques for the AWGN channel,</em></span> IEEE Trans. Communications,
  294. Vol. 48, No. 10, pp. 1681-1691, 2000.
  295. </li>
  296. <li class="listitem">
  297. Hoyer, Patrik O. <span class="emphasis"><em>Non-negative matrix factorization with sparseness
  298. constraints.</em></span>, Journal of machine learning research 5.Nov (2004):
  299. 1457-1469.
  300. </li>
  301. </ul></div>
  302. </div>
  303. <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
  304. <td align="left"></td>
  305. <td align="right"><div class="copyright-footer">Copyright &#169; 2006-2019 Nikhar
  306. Agrawal, Anton Bikineev, Paul A. Bristow, Marco Guazzone, Christopher Kormanyos,
  307. Hubert Holin, Bruno Lalande, John Maddock, Jeremy Murphy, Matthew Pulver, Johan
  308. R&#229;de, Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg,
  309. Daryle Walker and Xiaogang Zhang<p>
  310. Distributed under the Boost Software License, Version 1.0. (See accompanying
  311. 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>)
  312. </p>
  313. </div></td>
  314. </tr></table>
  315. <hr>
  316. <div class="spirit-nav">
  317. <a accesskey="p" href="bivariate_statistics.html"><img src="../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../statistics.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="anderson_darling.html"><img src="../../../../../doc/src/images/next.png" alt="Next"></a>
  318. </div>
  319. </body>
  320. </html>