<!DOCTYPE html> <html xmlns="http://www.w3.org/1999/xhtml"> <head> <meta charset="utf-8" /> <meta http-equiv="Content-Type" content="text/html; charset=utf-8" /> <meta name="generator" content="pandoc" /> <meta name="author" content="Lei Sun" /> <meta name="date" content="2017-05-16" /> <title>Fitting N\left(0, \sigma^2\right) with Gaussian Derivatives</title> <script src="site_libs/jquery-1.11.3/jquery.min.js"></script> <meta name="viewport" content="width=device-width, initial-scale=1" /> <link href="site_libs/bootstrap-3.3.5/css/readable.min.css" rel="stylesheet" /> <script src="site_libs/bootstrap-3.3.5/js/bootstrap.min.js"></script> <script src="site_libs/bootstrap-3.3.5/shim/html5shiv.min.js"></script> <script src="site_libs/bootstrap-3.3.5/shim/respond.min.js"></script> <script src="site_libs/jqueryui-1.11.4/jquery-ui.min.js"></script> <link href="site_libs/tocify-1.9.1/jquery.tocify.css" rel="stylesheet" /> <script src="site_libs/tocify-1.9.1/jquery.tocify.js"></script> <script src="site_libs/navigation-1.1/tabsets.js"></script> <link href="site_libs/highlightjs-9.12.0/textmate.css" rel="stylesheet" /> <script src="site_libs/highlightjs-9.12.0/highlight.js"></script> <link href="site_libs/font-awesome-4.5.0/css/font-awesome.min.css" rel="stylesheet" /> <style type="text/css">code{white-space: pre;}</style> <style type="text/css"> pre:not([class]) { background-color: white; } </style> <script type="text/javascript"> if (window.hljs) { hljs.configure({languages: []}); hljs.initHighlightingOnLoad(); if (document.readyState && document.readyState === "complete") { window.setTimeout(function() { hljs.initHighlighting(); }, 0); } } </script> <style type="text/css"> h1 { font-size: 34px; } h1.title { font-size: 38px; } h2 { font-size: 30px; } h3 { font-size: 24px; } h4 { font-size: 18px; } h5 { font-size: 16px; } h6 { font-size: 12px; } .table th:not([align]) { text-align: left; } </style> </head> <body> <style type = "text/css"> .main-container { max-width: 940px; margin-left: auto; margin-right: auto; } code { color: inherit; background-color: rgba(0, 0, 0, 0.04); } img { max-width:100%; height: auto; } .tabbed-pane { padding-top: 12px; } button.code-folding-btn:focus { outline: none; } </style> <style type="text/css"> /* padding for bootstrap navbar */ body { padding-top: 51px; padding-bottom: 40px; } /* offset scroll position for anchor links (for fixed navbar) */ .section h1 { padding-top: 56px; margin-top: -56px; } .section h2 { padding-top: 56px; margin-top: -56px; } .section h3 { padding-top: 56px; margin-top: -56px; } .section h4 { padding-top: 56px; margin-top: -56px; } .section h5 { padding-top: 56px; margin-top: -56px; } .section h6 { padding-top: 56px; margin-top: -56px; } </style> <script> // manage active state of menu based on current page $(document).ready(function () { // active menu anchor href = window.location.pathname href = href.substr(href.lastIndexOf('/') + 1) if (href === "") href = "index.html"; var menuAnchor = $('a[href="' + href + '"]'); // mark it active menuAnchor.parent().addClass('active'); // if it's got a parent navbar menu mark it active as well menuAnchor.closest('li.dropdown').addClass('active'); }); </script> <div class="container-fluid main-container"> <!-- tabsets --> <script> $(document).ready(function () { window.buildTabsets("TOC"); }); </script> <!-- code folding --> <script> $(document).ready(function () { // move toc-ignore selectors from section div to header $('div.section.toc-ignore') .removeClass('toc-ignore') .children('h1,h2,h3,h4,h5').addClass('toc-ignore'); // establish options var options = { selectors: "h1,h2,h3", theme: "bootstrap3", context: '.toc-content', hashGenerator: function (text) { return text.replace(/[.\\/?&!#<>]/g, '').replace(/\s/g, '_').toLowerCase(); }, ignoreSelector: ".toc-ignore", scrollTo: 0 }; options.showAndHide = true; options.smoothScroll = true; // tocify var toc = $("#TOC").tocify(options).data("toc-tocify"); }); </script> <style type="text/css"> #TOC { margin: 25px 0px 20px 0px; } @media (max-width: 768px) { #TOC { position: relative; width: 100%; } } .toc-content { padding-left: 30px; padding-right: 40px; } div.main-container { max-width: 1200px; } div.tocify { width: 20%; max-width: 260px; max-height: 85%; } @media (min-width: 768px) and (max-width: 991px) { div.tocify { width: 25%; } } @media (max-width: 767px) { div.tocify { width: 100%; max-width: none; } } .tocify ul, .tocify li { line-height: 20px; } .tocify-subheader .tocify-item { font-size: 0.90em; padding-left: 25px; text-indent: 0; } .tocify .list-group-item { border-radius: 0px; } </style> <!-- setup 3col/9col grid for toc_float and main content --> <div class="row-fluid"> <div class="col-xs-12 col-sm-4 col-md-3"> <div id="TOC" class="tocify"> </div> </div> <div class="toc-content col-xs-12 col-sm-8 col-md-9"> <div class="navbar navbar-default navbar-fixed-top" role="navigation"> <div class="container"> <div class="navbar-header"> <button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar"> <span class="icon-bar"></span> <span class="icon-bar"></span> <span class="icon-bar"></span> </button> <a class="navbar-brand" href="index.html">truncash</a> </div> <div id="navbar" class="navbar-collapse collapse"> <ul class="nav navbar-nav"> <li> <a href="about.html">About</a> </li> <li> <a href="license.html">License</a> </li> </ul> <ul class="nav navbar-nav navbar-right"> <li> <a href="https://github.com/LSun/truncash"> <span class="fa fa-github"></span> </a> </li> </ul> </div><!--/.nav-collapse --> </div><!--/.container --> </div><!--/.navbar --> <!-- Add a small amount of space between sections. --> <style type="text/css"> div.section { padding-top: 12px; } </style> <div class="fluid-row" id="header"> <h1 class="title toc-ignore">Fitting <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> with Gaussian Derivatives</h1> <h4 class="author"><em>Lei Sun</em></h4> <h4 class="date"><em>2017-05-16</em></h4> </div> <p><strong>Last updated:</strong> 2018-05-15</p> <strong>workflowr checks:</strong> <small>(Click a bullet for more information)</small> <ul> <li> <details> <p><summary> <strong style="color:blue;">✔</strong> <strong>R Markdown file:</strong> up-to-date </summary></p> <p>Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.</p> </details> </li> <li> <details> <p><summary> <strong style="color:blue;">✔</strong> <strong>Environment:</strong> empty </summary></p> <p>Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.</p> </details> </li> <li> <details> <p><summary> <strong style="color:blue;">✔</strong> <strong>Seed:</strong> <code>set.seed(12345)</code> </summary></p> <p>The command <code>set.seed(12345)</code> was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.</p> </details> </li> <li> <details> <p><summary> <strong style="color:blue;">✔</strong> <strong>Session information:</strong> recorded </summary></p> <p>Great job! Recording the operating system, R version, and package versions is critical for reproducibility.</p> </details> </li> <li> <details> <p><summary> <strong style="color:blue;">✔</strong> <strong>Repository version:</strong> <a href="https://github.com/LSun/truncash/tree/388e65e06000e313c170a82f3ed57346f6024897" target="_blank">388e65e</a> </summary></p> Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility. The version displayed above was the version of the Git repository at the time these results were generated. <br><br> Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use <code>wflow_publish</code> or <code>wflow_git_commit</code>). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated: <pre><code> Ignored files: Ignored: .DS_Store Ignored: .Rhistory Ignored: .Rproj.user/ Ignored: analysis/.DS_Store Ignored: analysis/BH_robustness_cache/ Ignored: analysis/FDR_Null_cache/ Ignored: analysis/FDR_null_betahat_cache/ Ignored: analysis/Rmosek_cache/ Ignored: analysis/StepDown_cache/ Ignored: analysis/alternative2_cache/ Ignored: analysis/alternative_cache/ Ignored: analysis/ash_gd_cache/ Ignored: analysis/average_cor_gtex_2_cache/ Ignored: analysis/average_cor_gtex_cache/ Ignored: analysis/brca_cache/ Ignored: analysis/cash_deconv_cache/ Ignored: analysis/cash_fdr_1_cache/ Ignored: analysis/cash_fdr_2_cache/ Ignored: analysis/cash_fdr_3_cache/ Ignored: analysis/cash_fdr_4_cache/ Ignored: analysis/cash_fdr_5_cache/ Ignored: analysis/cash_fdr_6_cache/ Ignored: analysis/cash_plots_cache/ Ignored: analysis/cash_sim_1_cache/ Ignored: analysis/cash_sim_2_cache/ Ignored: analysis/cash_sim_3_cache/ Ignored: analysis/cash_sim_4_cache/ Ignored: analysis/cash_sim_5_cache/ Ignored: analysis/cash_sim_6_cache/ Ignored: analysis/cash_sim_7_cache/ Ignored: analysis/correlated_z_2_cache/ Ignored: analysis/correlated_z_3_cache/ Ignored: analysis/correlated_z_cache/ Ignored: analysis/create_null_cache/ Ignored: analysis/cutoff_null_cache/ Ignored: analysis/design_matrix_2_cache/ Ignored: analysis/design_matrix_cache/ Ignored: analysis/diagnostic_ash_cache/ Ignored: analysis/diagnostic_correlated_z_2_cache/ Ignored: analysis/diagnostic_correlated_z_3_cache/ Ignored: analysis/diagnostic_correlated_z_cache/ Ignored: analysis/diagnostic_plot_2_cache/ Ignored: analysis/diagnostic_plot_cache/ Ignored: analysis/efron_leukemia_cache/ Ignored: analysis/fitting_normal_cache/ Ignored: analysis/gaussian_derivatives_2_cache/ Ignored: analysis/gaussian_derivatives_3_cache/ Ignored: analysis/gaussian_derivatives_4_cache/ Ignored: analysis/gaussian_derivatives_5_cache/ Ignored: analysis/gaussian_derivatives_cache/ Ignored: analysis/gd-ash_cache/ Ignored: analysis/gd_delta_cache/ Ignored: analysis/gd_lik_2_cache/ Ignored: analysis/gd_lik_cache/ Ignored: analysis/gd_w_cache/ Ignored: analysis/knockoff_10_cache/ Ignored: analysis/knockoff_2_cache/ Ignored: analysis/knockoff_3_cache/ Ignored: analysis/knockoff_4_cache/ Ignored: analysis/knockoff_5_cache/ Ignored: analysis/knockoff_6_cache/ Ignored: analysis/knockoff_7_cache/ Ignored: analysis/knockoff_8_cache/ Ignored: analysis/knockoff_9_cache/ Ignored: analysis/knockoff_cache/ Ignored: analysis/knockoff_var_cache/ Ignored: analysis/marginal_z_alternative_cache/ Ignored: analysis/marginal_z_cache/ Ignored: analysis/mosek_reg_2_cache/ Ignored: analysis/mosek_reg_4_cache/ Ignored: analysis/mosek_reg_5_cache/ Ignored: analysis/mosek_reg_6_cache/ Ignored: analysis/mosek_reg_cache/ Ignored: analysis/pihat0_null_cache/ Ignored: analysis/plot_diagnostic_cache/ Ignored: analysis/poster_obayes17_cache/ Ignored: analysis/real_data_simulation_2_cache/ Ignored: analysis/real_data_simulation_3_cache/ Ignored: analysis/real_data_simulation_4_cache/ Ignored: analysis/real_data_simulation_5_cache/ Ignored: analysis/real_data_simulation_cache/ Ignored: analysis/rmosek_primal_dual_2_cache/ Ignored: analysis/rmosek_primal_dual_cache/ Ignored: analysis/seqgendiff_cache/ Ignored: analysis/simulated_correlated_null_2_cache/ Ignored: analysis/simulated_correlated_null_3_cache/ Ignored: analysis/simulated_correlated_null_cache/ Ignored: analysis/simulation_real_se_2_cache/ Ignored: analysis/simulation_real_se_cache/ Ignored: analysis/smemo_2_cache/ Ignored: data/LSI/ Ignored: docs/.DS_Store Ignored: docs/figure/.DS_Store Ignored: output/fig/ </code></pre> Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes. </details> </li> </ul> <details> <summary> <small><strong>Expand here to see past versions:</strong></small> </summary> <ul> <table style="border-collapse:separate; border-spacing:5px;"> <thead> <tr> <th style="text-align:left;"> File </th> <th style="text-align:left;"> Version </th> <th style="text-align:left;"> Author </th> <th style="text-align:left;"> Date </th> <th style="text-align:left;"> Message </th> </tr> </thead> <tbody> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/e05bc836b3c74dc6ebca415afb5938675d6c5436/docs/fitting_normal.html" target="_blank">e05bc83</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2018-05-12 </td> <td style="text-align:left;"> Update to 1.0 </td> </tr> <tr> <td style="text-align:left;"> rmd </td> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/cc0ab8379469bc3726f1508cd81e4ecd6fef1b1a/analysis/fitting_normal.rmd" target="_blank">cc0ab83</a> </td> <td style="text-align:left;"> Lei Sun </td> <td style="text-align:left;"> 2018-05-11 </td> <td style="text-align:left;"> update </td> </tr> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/0f36d998db26444c5dd01502ea1af7fbd1129b22/docs/fitting_normal.html" target="_blank">0f36d99</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-12-21 </td> <td style="text-align:left;"> Build site. </td> </tr> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/853a484bfacf347e109f6c8fb3ffaab5f4d6cc02/docs/fitting_normal.html" target="_blank">853a484</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-11-07 </td> <td style="text-align:left;"> Build site. </td> </tr> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/fa2c24ead93d0109e943745c837976a7b163d97c/docs/fitting_normal.html" target="_blank">fa2c24e</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-11-06 </td> <td style="text-align:left;"> transfer </td> </tr> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/87193bad023b2f8facb2b6f6ac5f06330bfe83d5/docs/fitting_normal.html" target="_blank">87193ba</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> <td style="text-align:left;"> fitting normal </td> </tr> <tr> <td style="text-align:left;"> rmd </td> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/055f9f0516cae31b3e345f95c18561c9ab952486/analysis/fitting_normal.rmd" target="_blank">055f9f0</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> <td style="text-align:left;"> websites </td> </tr> <tr> <td style="text-align:left;"> html </td> <td style="text-align:left;"> <a href="https://cdn.rawgit.com/LSun/truncash/a5300252ff521211d8a6c6a80a7cbb58b4d89438/docs/fitting_normal.html" target="_blank">a530025</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> <td style="text-align:left;"> mean 0 normal </td> </tr> <tr> <td style="text-align:left;"> rmd </td> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/4dab9a497efbd26184ea6fbde4b0fdf684486086/analysis/fitting_normal.rmd" target="_blank">4dab9a4</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> <td style="text-align:left;"> mean 0 normal </td> </tr> </tbody> </table> </ul> </details> <hr /> <div id="introduction" class="section level2"> <h2>Introduction</h2> <p>We know the empirical distribution of a number of correlated null <span class="math inline">\(N\left(0, 1\right)\)</span> <span class="math inline">\(z\)</span> scores can be approximated by Gaussian derivatives at least theoretically. Now we want to know if the <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> density can also be approximated by Gaussian derivatives. The question is important because it is related to the identifiability between the correlated null and true signals using the tool of Gaussian derivatives.</p> <p>The same question was investigated for <a href="alternative.html">small <span class="math inline">\(\sigma^2\)</span></a> and <a href="alternative2.html">large <span class="math inline">\(\sigma^2\)</span></a>. But these investigations was <em>empirical</em>, in the sense that we are fitting the Gaussian derivatives to a large number of independent samples simulated from <span class="math inline">\(N\left(0, \sigma^2\right)\)</span>. Now we are doing that <strong>theoretically</strong>.</p> </div> <div id="hermite-moments" class="section level2"> <h2>Hermite moments</h2> <p>If a probability density function (pdf) <span class="math inline">\(f\left(x\right)\)</span> can be approximated by normalized Gaussian derivatives, we should have</p> <p><span class="math display">\[ f\left(x\right) = w_0\varphi\left(x\right) + w_1\varphi^{(1)}\left(x\right) + w_2\frac{1}{\sqrt{2!}}\varphi^{(2)}\left(x\right) + \cdots + w_n\frac{1}{\sqrt{n!}}\varphi^{(n)}\left(x\right) + \cdots \ . \]</span> Using the fact that Hermite polynomials are orthonormal with respect to the Gaussian kernel and normalizing constants, we have</p> <p><span class="math display">\[ \int_{\mathbb{R}} \frac{1}{\sqrt{m!}} h_m\left(x\right) \frac{1}{\sqrt{n!}}\varphi^{(n)}\left(x\right)dx = \left(-1\right)^n\delta_{mn} \ . \]</span> Therefore, when using Gaussian derivatives to decompose a pdf <span class="math inline">\(f\)</span>, the coefficients <span class="math inline">\(w_n\)</span> can be obtained by</p> <p><span class="math display">\[ w_n = \left(-1\right)^n \int \frac{1}{\sqrt{n!}} h_n(x) f(x) dx = \frac{\left(-1\right)^n}{\sqrt{n!}} E_f\left[h_n\left(x\right)\right] \ . \]</span> Since <span class="math inline">\(h_n\)</span> is a degree <span class="math inline">\(n\)</span> polynomials, <span class="math inline">\(E_f\left[h_n\left(x\right)\right]\)</span> is essentially a linear combination of relevant moments up to order <span class="math inline">\(n\)</span>, denoted as <span class="math inline">\(E\left[x^i\right]\)</span>. For example, here are the first several <span class="math inline">\(E_f\left[h_n\left(x\right)\right]\)</span>.</p> <p><span class="math display">\[ \begin{array}{rclcl} E_f\left[h_0\left(x\right)\right] &=& E_f\left[1\right] &=& 1 \\ E_f\left[h_1\left(x\right)\right] &=& E_f\left[x\right] &=& E\left[x\right] \\ E_f\left[h_2\left(x\right)\right] &=& E_f\left[x^2 - 1\right] &=& E\left[x^2\right] - 1\\ E_f\left[h_3\left(x\right)\right] &=& E_f\left[x^3 - 3x\right] &=& E\left[x^3\right] - 3 E\left[x\right]\\ E_f\left[h_4\left(x\right)\right] &=& E_f\left[x^4 - 6x^2 + 3\right] &=& E\left[x^4\right] - 6 E\left[x^2\right] + 3\\ E_f\left[h_5\left(x\right)\right] &=& E_f\left[x^5 - 10x^3 + 15x\right] &=& E\left[x^5\right] - 10 E\left[x^3\right] + 15 E\left[x\right]\\ E_f\left[h_6\left(x\right)\right] &=& E_f\left[x^6 - 15x^4 + 45x^2 - 15\right] &=& E\left[x^6\right] - 15 E\left[x^4\right] + 45E\left[x^2\right] - 15\\ \end{array} \]</span> Some call these “Hermite moments.”</p> </div> <div id="when-f-nleft0-sigma2right-general-results" class="section level2"> <h2>When <span class="math inline">\(f = N\left(0, \sigma^2\right)\)</span>: General results</h2> <p>When <span class="math inline">\(f = N\left(0, \sigma^2\right)\)</span>, we can write all Hermite moments out analytically. After algebra, it turns out they can be written as</p> <p><span class="math display">\[ E_f\left[h_n\left(x\right)\right] = \begin{cases} 0 & n \text{ is odd}\\ \left(n - 1\right)!! \left(\sigma^2 - 1\right)^{n / 2} & n \text{ is even} \end{cases} \ , \]</span> where <span class="math inline">\(\left(n - 1\right)!!\)</span> is the <a href="https://en.wikipedia.org/wiki/Double_factorial">double factorial</a>.</p> <p>Hence, if we want to express the pdf of <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> in terms of normalized Gaussian derivatives, the coefficient</p> <p><span class="math display">\[ w_n = \begin{cases} 0 & n \text{ is odd}\\ \frac{\left(n - 1\right)!!}{\sqrt{n!}}\left(\sigma^2 - 1\right)^{n / 2} =\sqrt{\frac{\left(n - 1\right)!!}{n!!}}\left(\sigma^2 - 1\right)^{n / 2} & n \text{ is even} \end{cases} \ . \]</span> From now on we will only pay attention to the even-order coefficients <span class="math inline">\(w_n = \sqrt{\frac{\left(n - 1\right)!!}{n!!}}\left(\sigma^2 - 1\right)^{n / 2} := U_n\left(\sigma^2 - 1\right)^{n / 2}\)</span>, where the sequence <span class="math inline">\(U_n := \sqrt{\frac{\left(n - 1\right)!!}{n!!}}\)</span>. Note that this sequence <span class="math inline">\(U_n\)</span> is interesting. <a href="https://math.stackexchange.com/questions/584456/limit-of-a-fraction-of-double-factorials">It can be proved</a> that it’s going to zero as <span class="math inline">\(n\to\infty\)</span>. Actually it can be <a href="https://en.wikipedia.org/wiki/Double_factorial#Additional_identities">written in another way</a></p> <p><span class="math display">\[ \int_{0}^{\frac\pi2} \sin^nx dx = \int_{0}^{\frac\pi2} \cos^nx dx = \frac\pi2U_n^2 \ . \]</span> However, this sequence decays slowly. In particular, it doesn’t decay exponentially, as seen in the following plots.</p> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-3-1.png" width="672" style="display: block; margin: auto;" /></p> <details> <summary><em>Expand here to see past versions of unnamed-chunk-3-1.png:</em></summary> <table style = "border-collapse:separate; border-spacing:5px;"> <thead> <tr> <th style="text-align:left;"> Version </th> <th style="text-align:left;"> Author </th> <th style="text-align:left;"> Date </th> </tr> </thead> <tbody> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/0f36d998db26444c5dd01502ea1af7fbd1129b22/docs/figure/fitting_normal.rmd/unnamed-chunk-3-1.png" target="_blank">0f36d99</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-12-21 </td> </tr> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/a5300252ff521211d8a6c6a80a7cbb58b4d89438/docs/figure/fitting_normal.rmd/unnamed-chunk-3-1.png" target="_blank">a530025</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> </tr> </tbody> </table> </details> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-3-2.png" width="672" style="display: block; margin: auto;" /></p> <details> <summary><em>Expand here to see past versions of unnamed-chunk-3-2.png:</em></summary> <table style = "border-collapse:separate; border-spacing:5px;"> <thead> <tr> <th style="text-align:left;"> Version </th> <th style="text-align:left;"> Author </th> <th style="text-align:left;"> Date </th> </tr> </thead> <tbody> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/0f36d998db26444c5dd01502ea1af7fbd1129b22/docs/figure/fitting_normal.rmd/unnamed-chunk-3-2.png" target="_blank">0f36d99</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-12-21 </td> </tr> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/a5300252ff521211d8a6c6a80a7cbb58b4d89438/docs/figure/fitting_normal.rmd/unnamed-chunk-3-2.png" target="_blank">a530025</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> </tr> </tbody> </table> </details> <p><strong>Therefore, when <span class="math inline">\(\sigma^2 > 2\)</span>, <span class="math inline">\(w_n\)</span> is exploding!</strong> It means theoretially we cannot actually fit <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> with Gaussian derivatives when <span class="math inline">\(\sigma^2 > 2\)</span>.</p> <p>It has multiple implications. First, it shows why people usually say Gaussian derivatives can only fit a density that’s close enough to Gaussian, especially, a density whose variance is too “inflated” compared with the standard normal. Therefore, the fact that in many case, the correlated null can be satisfactorily fitted by Gaussian derivatives tells us the inflation caused by correlation is indeed peculiar.</p> <p>On the other hand, when <span class="math inline">\(0 < \sigma^2 < 2\)</span>, that is, <span class="math inline">\(\left|\sigma^2 - 1\right| < 1\)</span>, <span class="math inline">\(w_n\)</span> decays exponentially. It means that we are able to satisfactorily fit <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> with a limited number of Gaussian derivatives such that</p> <p><span class="math display">\[ f_N\left(x\right) = \varphi\left(x\right) + \sum\limits_{n = 1}^N w_n\frac{1}{\sqrt{n!}}\varphi^{(n)}\left(x\right) \approx N\left(0, \sigma^2\right) = \frac{1}{\sigma}\varphi\left(\frac{x}{\sigma}\right) \ . \]</span></p> <p>It also means that when the correlated null looks like <span class="math inline">\(N\left(0, \sigma^2\right)\)</span> with <span class="math inline">\(\sigma^2\in\left(1, 2\right]\)</span> (small inflation) or <span class="math inline">\(\sigma^2\in\left(0, 1\right)\)</span> (deflation), it’s hard to identify whether the inflation or deflation is caused by correlation or true effects.</p> </div> <div id="when-f-nleft0-sqrt22right-in-particular" class="section level2"> <h2>When <span class="math inline">\(f = N\left(0, \sqrt{2}^2\right)\)</span> in particular</h2> <p>That’s why in previous empirical investigations, we found that we <a href="alternative.html">could fit <span class="math inline">\(N\left(0, \sqrt{2}^2\right)\)</span> relatively well</a>, but <a href="alternative2.html">could not fit when <span class="math inline">\(\sigma^2 > 2\)</span></a>.</p> <p>Actually, <strong><span class="math inline">\(N\left(0, \sqrt{2}^2\right)\)</span> is a singularly interesting case</strong>. This case corresponds to when we have <span class="math inline">\(N\left(0, 1\right)\)</span> signal and <span class="math inline">\(N\left(0, 1\right)\)</span> independent noise, hence <span class="math inline">\(SNR = 0\)</span>. Theoretically the density curve can be fitted when the odd-order coefficients are <span class="math inline">\(0\)</span> and the even-order ones are <span class="math inline">\(U_n = \sqrt{\frac{\left(n - 1\right)!!}{n!!}}\)</span>.</p> <p>However, since <span class="math inline">\(U_n\)</span> is not decaying fast enough, the re-constructed curve using a limited number of Gaussian derivatives doesn’t look good. Meanwhile, if we truly have many random samples from <span class="math inline">\(N\left(0, \sqrt{2}^2\right)\)</span>, we can fit the data using Gaussian derivatives, and the estimated coefficients <span class="math inline">\(\hat w\)</span> usually give a better fit, as <a href="alternative.html">seen previously</a>.</p> <p>Here we plot the fitted curve using first <span class="math inline">\(N = 10, 50, 100\)</span> orders of Gaussian derivatives, compared with the true <span class="math inline">\(N\left(0, \sqrt{2}^2\right)\)</span> density and the curve estimated from random samples using <span class="math inline">\(L = 10\)</span> Gaussian derivatives.</p> <pre><code>Estimated Normalized w: 0 ~ 1 ; 1 ~ -0.0126 ; 2 ~ 0.72013 ; 3 ~ -0.01709 ; 4 ~ 0.6782 ; 5 ~ 0.01241 ; 6 ~ 0.62481 ; 7 ~ 0.03822 ; 8 ~ 0.4004 ; 9 ~ 0.00207 ; 10 ~ 0.11385 ; </code></pre> <pre><code>Theoretical Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ 0.70711 ; 3 ~ 0 ; 4 ~ 0.61237 ; 5 ~ 0 ; 6 ~ 0.55902 ; 7 ~ 0 ; 8 ~ 0.52291 ; 9 ~ 0 ; 10 ~ 0.49608 ; ...</code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-4-1.png" width="672" style="display: block; margin: auto;" /></p> <details> <summary><em>Expand here to see past versions of unnamed-chunk-4-1.png:</em></summary> <table style = "border-collapse:separate; border-spacing:5px;"> <thead> <tr> <th style="text-align:left;"> Version </th> <th style="text-align:left;"> Author </th> <th style="text-align:left;"> Date </th> </tr> </thead> <tbody> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/0f36d998db26444c5dd01502ea1af7fbd1129b22/docs/figure/fitting_normal.rmd/unnamed-chunk-4-1.png" target="_blank">0f36d99</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-12-21 </td> </tr> <tr> <td style="text-align:left;"> <a href="https://github.com/LSun/truncash/blob/a5300252ff521211d8a6c6a80a7cbb58b4d89438/docs/figure/fitting_normal.rmd/unnamed-chunk-4-1.png" target="_blank">a530025</a> </td> <td style="text-align:left;"> LSun </td> <td style="text-align:left;"> 2017-05-17 </td> </tr> </tbody> </table> </details> </div> <div id="sigma2-2-multiple-examples" class="section level2"> <h2><span class="math inline">\(0 < \sigma^2 < 2\)</span>: Multiple examples</h2> <div id="sigma2-0.2-n-40" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 0.2\)</span>, <span class="math inline">\(N = 40\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ -0.56569 ; 3 ~ 0 ; 4 ~ 0.39192 ; 5 ~ 0 ; 6 ~ -0.28622 ; 7 ~ 0 ; 8 ~ 0.21418 ; 9 ~ 0 ; 10 ~ -0.16255 ; 11 ~ 0 ; 12 ~ 0.12451 ; 13 ~ 0 ; 14 ~ -0.09598 ; 15 ~ 0 ; 16 ~ 0.07435 ; 17 ~ 0 ; 18 ~ -0.0578 ; 19 ~ 0 ; 20 ~ 0.04507 ; 21 ~ 0 ; 22 ~ -0.03523 ; 23 ~ 0 ; 24 ~ 0.02759 ; 25 ~ 0 ; 26 ~ -0.02164 ; 27 ~ 0 ; 28 ~ 0.017 ; 29 ~ 0 ; 30 ~ -0.01337 ; 31 ~ 0 ; 32 ~ 0.01053 ; 33 ~ 0 ; 34 ~ -0.0083 ; 35 ~ 0 ; 36 ~ 0.00655 ; 37 ~ 0 ; 38 ~ -0.00517 ; 39 ~ 0 ; 40 ~ 0.00408 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-5-1.png" width="672" style="display: block; margin: auto;" /></p> </div> <div id="sigma2-0.5-n-14" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 0.5\)</span>, <span class="math inline">\(N = 14\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ -0.35355 ; 3 ~ 0 ; 4 ~ 0.15309 ; 5 ~ 0 ; 6 ~ -0.06988 ; 7 ~ 0 ; 8 ~ 0.03268 ; 9 ~ 0 ; 10 ~ -0.0155 ; 11 ~ 0 ; 12 ~ 0.00742 ; 13 ~ 0 ; 14 ~ -0.00358 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-6-1.png" width="672" style="display: block; margin: auto;" /></p> </div> <div id="sigma2-0.8-n-6" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 0.8\)</span>, <span class="math inline">\(N = 6\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ -0.14142 ; 3 ~ 0 ; 4 ~ 0.02449 ; 5 ~ 0 ; 6 ~ -0.00447 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-7-1.png" width="672" style="display: block; margin: auto;" /></p> </div> <div id="sigma2-1.2-n-6" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 1.2\)</span>, <span class="math inline">\(N = 6\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ 0.14142 ; 3 ~ 0 ; 4 ~ 0.02449 ; 5 ~ 0 ; 6 ~ 0.00447 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-8-1.png" width="672" style="display: block; margin: auto;" /></p> </div> <div id="sigma2-1.5-n-14" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 1.5\)</span>, <span class="math inline">\(N = 14\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ 0.35355 ; 3 ~ 0 ; 4 ~ 0.15309 ; 5 ~ 0 ; 6 ~ 0.06988 ; 7 ~ 0 ; 8 ~ 0.03268 ; 9 ~ 0 ; 10 ~ 0.0155 ; 11 ~ 0 ; 12 ~ 0.00742 ; 13 ~ 0 ; 14 ~ 0.00358 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-9-1.png" width="672" style="display: block; margin: auto;" /></p> </div> <div id="sigma2-1.8-n-40" class="section level3"> <h3><span class="math inline">\(\sigma^2 = 1.8\)</span>, <span class="math inline">\(N = 40\)</span></h3> <pre><code>Normalized w: 0 ~ 1 ; 1 ~ 0 ; 2 ~ 0.56569 ; 3 ~ 0 ; 4 ~ 0.39192 ; 5 ~ 0 ; 6 ~ 0.28622 ; 7 ~ 0 ; 8 ~ 0.21418 ; 9 ~ 0 ; 10 ~ 0.16255 ; 11 ~ 0 ; 12 ~ 0.12451 ; 13 ~ 0 ; 14 ~ 0.09598 ; 15 ~ 0 ; 16 ~ 0.07435 ; 17 ~ 0 ; 18 ~ 0.0578 ; 19 ~ 0 ; 20 ~ 0.04507 ; 21 ~ 0 ; 22 ~ 0.03523 ; 23 ~ 0 ; 24 ~ 0.02759 ; 25 ~ 0 ; 26 ~ 0.02164 ; 27 ~ 0 ; 28 ~ 0.017 ; 29 ~ 0 ; 30 ~ 0.01337 ; 31 ~ 0 ; 32 ~ 0.01053 ; 33 ~ 0 ; 34 ~ 0.0083 ; 35 ~ 0 ; 36 ~ 0.00655 ; 37 ~ 0 ; 38 ~ 0.00517 ; 39 ~ 0 ; 40 ~ 0.00408 ; </code></pre> <p><img src="figure/fitting_normal.rmd/unnamed-chunk-10-1.png" width="672" style="display: block; margin: auto;" /></p> </div> </div> <div id="another-special-case-sigma2-0" class="section level2"> <h2>Another special case: <span class="math inline">\(\sigma^2 = 0\)</span></h2> <p>When <span class="math inline">\(\sigma^2 = 0\)</span>, <span class="math inline">\(f = N\left(0, 0\right) = \delta_0\)</span>. In this case,</p> <p><span class="math display">\[ w_n = \begin{cases} 0 & n \text{ is odd}\\ \sqrt{\frac{\left(n - 1\right)!!}{n!!}}\left(-1\right)^{n / 2} & n \text{ is even} \end{cases} \]</span> is equivalent to what we’ve <a href="gaussian_derivatives_4.html#extreme_case:_(rho_%7Bij%7D_equiv_1)">obtained previously</a></p> <p><span class="math display">\[ w_n = \frac{1}{\sqrt{n!}}h_n\left(0\right) \ . \]</span> The performance can be seen <a href="gd_delta.html#(z_=_0)">here</a>.</p> </div> <div id="session-information" class="section level2"> <h2>Session information</h2> <pre class="r"><code>sessionInfo()</code></pre> <pre><code>R version 3.4.3 (2017-11-30) Platform: x86_64-apple-darwin15.6.0 (64-bit) Running under: macOS High Sierra 10.13.4 Matrix products: default BLAS: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRblas.0.dylib LAPACK: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRlapack.dylib locale: [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] Rmosek_8.0.69 CVXR_0.95 REBayes_1.2 Matrix_1.2-12 [5] SQUAREM_2017.10-1 EQL_1.0-0 ttutils_1.0-1 PolynomF_1.0-1 loaded via a namespace (and not attached): [1] Rcpp_0.12.16 knitr_1.20 whisker_0.3-2 [4] magrittr_1.5 workflowr_1.0.1 bit_1.1-12 [7] lattice_0.20-35 R6_2.2.2 stringr_1.3.0 [10] tools_3.4.3 grid_3.4.3 R.oo_1.21.0 [13] git2r_0.21.0 scs_1.1-1 htmltools_0.3.6 [16] bit64_0.9-7 yaml_2.1.18 rprojroot_1.3-2 [19] digest_0.6.15 gmp_0.5-13.1 ECOSolveR_0.4 [22] R.utils_2.6.0 evaluate_0.10.1 rmarkdown_1.9 [25] stringi_1.1.6 Rmpfr_0.6-1 compiler_3.4.3 [28] backports_1.1.2 R.methodsS3_1.7.1</code></pre> </div> <!-- Adjust MathJax settings so that all math formulae are shown using TeX fonts only; see http://docs.mathjax.org/en/latest/configuration.html. This will make the presentation more consistent at the cost of the webpage sometimes taking slightly longer to load. Note that this only works because the footer is added to webpages before the MathJax javascript. --> <script type="text/x-mathjax-config"> MathJax.Hub.Config({ "HTML-CSS": { availableFonts: ["TeX"] } }); </script> <hr> <p> This reproducible <a href="http://rmarkdown.rstudio.com">R Markdown</a> analysis was created with <a href="https://github.com/jdblischak/workflowr">workflowr</a> 1.0.1 </p> <hr> </div> </div> </div> <script> // add bootstrap table styles to pandoc tables function bootstrapStylePandocTables() { $('tr.header').parent('thead').parent('table').addClass('table table-condensed'); } $(document).ready(function () { bootstrapStylePandocTables(); }); </script> <!-- dynamically load mathjax for compatibility with self-contained --> <script> (function () { var script = document.createElement("script"); script.type = "text/javascript"; script.src = "https://mathjax.rstudio.com/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"; document.getElementsByTagName("head")[0].appendChild(script); })(); </script> </body> </html>