You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
<h1>Implementation in <spanclass="math inline">\(\texttt{deep\_tensor}\)</span></h1>
486
486
<p>We will now use <spanclass="math inline">\(\texttt{deep\_tensor}\)</span> to construct a DIRT approximation to the posterior. To accelerate this process, we will use a reduced order model in place of the full model. Then, we will illustrate some debiasing techniques which use the DIRT approximation to the posterior, in combination with the full model, to accelerate the process of drawing exact posterior samples.</p>
<spanid="cb1-4"><ahref="#cb1-4" aria-hidden="true" tabindex="-1"></a><spanclass="im">import</span> deep_tensor <spanclass="im">as</span> dt</span></code><buttontitle="Copy to Clipboard" class="code-copy-button"><iclass="bi"></i></button></pre></div>
492
492
</div>
493
493
<p>We begin by defining the prior, (full) model and reduced order model.</p>
494
494
<p>The full model is implemented in <ahref="https://fenicsproject.org/download/archive/">FEniCS</a>, on a <spanclass="math inline">\(96 \times 32\)</span> grid, using piecwise linear basis functions. Timestepping is done using the backward Euler method. The reduced order model is constructed using the proper orthogonal decomposition <spanclass="citation" data-cites="Benner2015">(see, <em>e.g.</em>, <ahref="#ref-Benner2015" role="doc-biblioref">Benner, Gugercin, and Willcox 2015</a>)</span>.</p>
<spanid="cb2-3"><ahref="#cb2-3" aria-hidden="true" tabindex="-1"></a><spanclass="co"># Construct the prior, full model and reduced order model</span></span>
499
499
<spanid="cb2-4"><ahref="#cb2-4" aria-hidden="true" tabindex="-1"></a>prior, model, rom <spanclass="op">=</span> setup_heat_problem()</span></code><buttontitle="Copy to Clipboard" class="code-copy-button"><iclass="bi"></i></button></pre></div>
500
500
</div>
501
501
<p>Next, we will generate the true log-diffusion coefficient using a sample from the prior. The true log-diffusion coefficient is plotted in <ahref="#fig-ktrue" class="quarto-xref">Figure 1</a>.</p>
<spanid="cb3-2"><ahref="#cb3-2" aria-hidden="true" tabindex="-1"></a>logk_true <spanclass="op">=</span> prior.transform(xi_true)</span></code><buttontitle="Copy to Clipboard" class="code-copy-button"><iclass="bi"></i></button></pre></div>
505
505
</div>
@@ -528,7 +528,7 @@ <h1>Implementation in <span class="math inline">\(\texttt{deep\_tensor}\)</span>
528
528
</div>
529
529
</div>
530
530
<p>Next, we will solve the (full) model to obtain the modelled temperatures corresponding to the true diffusion coefficient, and use these to generate some synthetic data. <ahref="#fig-utrue" class="quarto-xref">Figure 2</a> shows the true temperature field at time <spanclass="math inline">\(T=10\)</span>, as well as the observation locations.</p>
<h2class="anchored" data-anchor-id="building-the-dirt-object">Building the DIRT Object</h2>
571
571
<p>Now we will build a DIRT object to approximate the posterior density of the log-diffusion coefficient for the reduced-order model. We begin by defining functions which return the potential associated with the likelihood and prior.</p>
<spanid="cb7-2"><ahref="#cb7-2" aria-hidden="true" tabindex="-1"></a><spanclass="co">"""Returns the negative log prior density evaluated a given set of </span></span>
@@ -599,27 +599,26 @@ <h2 class="anchored" data-anchor-id="building-the-dirt-object">Building the DIRT
599
599
<spanid="cb7-27"><ahref="#cb7-27" aria-hidden="true" tabindex="-1"></a><spanclass="cf">return</span> _negloglik(rom, xs)</span></code><buttontitle="Copy to Clipboard" class="code-copy-button"><iclass="bi"></i></button></pre></div>
600
600
</div>
601
601
<p>Next, we specify a preconditioner. Because the prior of the coefficients <spanclass="math inline">\(\{\xi^{(i)}\}_{i=1}^{d}\)</span> is the standard Gaussian, the mapping between a Gaussian reference and the prior is simply the identity mapping. This is an appropriate choice of preconditioner in the absence of any other information.</p>
<divclass="sourceCode cell-code" id="cb10"><preclass="sourceCode python code-with-copy"><codeclass="sourceCode python"><spanid="cb10-1"><ahref="#cb10-1" aria-hidden="true" tabindex="-1"></a><spanclass="co"># Reduce the initial and maximum tensor ranks to reduce the cost of </span></span>
613
-
<spanid="cb10-2"><ahref="#cb10-2" aria-hidden="true" tabindex="-1"></a><spanclass="co"># each layer</span></span>
<divclass="sourceCode cell-code" id="cb10"><preclass="sourceCode python code-with-copy"><codeclass="sourceCode python"><spanid="cb10-1"><ahref="#cb10-1" aria-hidden="true" tabindex="-1"></a><spanclass="co"># Reduce the initial and maximum tensor ranks to reduce the cost of each layer</span></span>
<p>First, we will illustrate how to use the DIRT density as part of an MCMC sampler. The simplest sampler, which we demonstrate here, is an independence sampler using the DIRT density as a proposal density.</p>
<divclass="sourceCode cell-code" id="cb11"><preclass="sourceCode python code-with-copy"><codeclass="sourceCode python"><spanid="cb11-1"><ahref="#cb11-1" aria-hidden="true" tabindex="-1"></a><spanclass="co"># Generate a set of samples from the DIRT density</span></span>
<spanid="cb6-6"><ahref="#cb6-6" aria-hidden="true" tabindex="-1"></a>)</span></code><buttontitle="Copy to Clipboard" class="code-copy-button"><iclass="bi"></i></button></pre></div>
855
855
<p>Evaluates the conditional inverse Rosenblatt transport.</p>
856
856
<p>Returns the conditional inverse Rosenblatt transport evaluated at a set of samples in the approximation domain.</p>
857
-
<p>The conditional inverse Rosenblatt transport takes the form</p>
858
-
<p><spanclass="math display">\[
859
-
X|Y = \mathcal{T}(\mathcal{R}_{k}(Y), R),
860
-
\]</span></p>
861
-
<p>where <spanclass="math inline">\(Y\)</span> is a <spanclass="math inline">\(k\)</span>-dimensional random variable, <spanclass="math inline">\(R\)</span> is a <spanclass="math inline">\(n-k\)</span>-dimensional reference random variable, <spanclass="math inline">\(\mathcal{R}(\,\cdot\,)\)</span> denotes the (full) Rosenblatt transport, <spanclass="math inline">\(\mathcal{T}(\,\cdot\,) := \mathcal{R}^{-1}(\,\cdot\,)\)</span>, denotes its inverse, and <spanclass="math inline">\(\mathcal{R}_{k}(\,\cdot\,)\)</span> denotes the Rosenblatt transport for the first (or last) <spanclass="math inline">\(k\)</span> variables.</p>
<p>A <spanclass="math inline">\(1 \times k\)</span> (if the same realisation of <spanclass="math inline">\(Y\)</span> is to be used for all samples in <code>rs</code>) or <spanclass="math inline">\(n \times k\)</span>matrix (if a different realisation of <spanclass="math inline">\(Y\)</span> is to be used for all samples in <code>rs</code>) containing samples from the approximation domain.</p>
862
+
<p>A matrix containing samples from the approximation domain. The matrix should have dimensions <spanclass="math inline">\(1 \times k\)</span> (if the same realisation of <spanclass="math inline">\(Y\)</span> is to be used for all samples in <code>rs</code>) or <spanclass="math inline">\(n \times k\)</span> (if a different realisation of <spanclass="math inline">\(Y\)</span> is to be used for each samples in <code>rs</code>).</p>
<p>The number of layers of the deep inverse Rosenblatt transport to push the samples forward under. If not specified, the samples will be pushed forward through all the layers.</p>
874
+
<p>The number of layers of the DIRT object to use when evaluating the CIRT. If not specified, all layers will be used.</p>
<p>An <spanclass="math inline">\(n\)</span>-dimensional vector containing the potential function of the approximation to the conditional density of <spanclass="math inline">\(X \textbar Y\)</span> evaluated at each sample in <code>rs</code>.</p>
887
+
<p>An <spanclass="math inline">\(n\)</span>-dimensional vector containing the potential function of the approximation to the conditional density of <spanclass="math inline">\(X \textbar Y\)</span> evaluated at each sample in <code>xs</code>.</p>
0 commit comments