|
10 | 10 | components[3] (prior = 0.2000): Distributions.Normal{Float64}(μ=20.0, σ=1.5) |
11 | 11 | </code></pre><p>This is a simple 1-dimensional distribution, so let's visualize it:</p><pre><code class="language-julia hljs">using StatsPlots |
12 | 12 | figsize = (800, 400) |
13 | | -plot(target_distribution; components=false, label=nothing, size=figsize)</code></pre><img src="30120a1f.svg" alt="Example block output"/><p>We can convert a <code>Distribution</code> from <a href="https://github.com/JuliaStats/Distributions.jl">Distributions.jl</a> into something we can pass to <code>sample</code> for many different samplers by implementing the <a href="https://github.com/tpapp/LogDensityProblems.jl">LogDensityProblems.jl</a> interface:</p><pre><code class="language-julia hljs">using LogDensityProblems: LogDensityProblems |
| 13 | +plot(target_distribution; components=false, label=nothing, size=figsize)</code></pre><img src="4f602687.svg" alt="Example block output"/><p>We can convert a <code>Distribution</code> from <a href="https://github.com/JuliaStats/Distributions.jl">Distributions.jl</a> into something we can pass to <code>sample</code> for many different samplers by implementing the <a href="https://github.com/tpapp/LogDensityProblems.jl">LogDensityProblems.jl</a> interface:</p><pre><code class="language-julia hljs">using LogDensityProblems: LogDensityProblems |
14 | 14 |
|
15 | 15 | struct DistributionLogDensity{D} |
16 | 16 | d::D |
|
57 | 57 | </span><span class="sgr90"> Symbol Float64 Float64 Float64 Float64 Float64 |
58 | 58 |
|
59 | 59 | x -5.6359 -3.2669 -1.4454 2.5327 5.3961 |
60 | | -</span></span></code></pre><pre><code class="language-julia hljs">plot(chain; size=figsize)</code></pre><img src="38b0b6c1.svg" alt="Example block output"/><p>This doesn't look quite like what we're expecting.</p><pre><code class="language-julia hljs">plot(target_distribution; components=false, linewidth=2) |
| 60 | +</span></span></code></pre><pre><code class="language-julia hljs">plot(chain; size=figsize)</code></pre><img src="1791bcdd.svg" alt="Example block output"/><p>This doesn't look quite like what we're expecting.</p><pre><code class="language-julia hljs">plot(target_distribution; components=false, linewidth=2) |
61 | 61 | density!(chain) |
62 | | -plot!(size=figsize)</code></pre><img src="47635022.svg" alt="Example block output"/><p>Notice how <code>chain</code> has zero probability mass in the left-most component of the mixture!</p><p>Let's instead try to use a <em>tempered</em> version of <code>RWMH</code>. <em>But</em> before we do that, we need to make sure that AdvancedMH.jl is compatible with MCMCTempering.jl.</p><p>To do that we need to implement two methods. First we need to tell MCMCTempering how to extract the parameters, and potentially the log-probabilities, from a <code>AdvancedMH.Transition</code>:</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="MCMCTempering.getparams_and_logprob-getting-started" href="#MCMCTempering.getparams_and_logprob-getting-started"><code>MCMCTempering.getparams_and_logprob</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">getparams_and_logprob([model, ]state)</code></pre><p>Return a vector of parameters from the <code>state</code>.</p><p>See also: <a href="../api/#MCMCTempering.setparams_and_logprob!!"><code>setparams_and_logprob!!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/c7dde12b17a62981062215319e01373119f12342/src/abstractmcmc.jl#L26-L32">source</a></section></article><p>And similarly, we need a way to <em>update</em> the parameters and the log-probabilities of a <code>AdvancedMH.Transition</code>:</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="MCMCTempering.setparams_and_logprob!!-getting-started" href="#MCMCTempering.setparams_and_logprob!!-getting-started"><code>MCMCTempering.setparams_and_logprob!!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">setparams_and_logprob!!([model, ]state, params)</code></pre><p>Set the parameters in the state to <code>params</code>, possibly mutating if it makes sense.</p><p>See also: <a href="../api/#MCMCTempering.getparams_and_logprob"><code>getparams_and_logprob</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/c7dde12b17a62981062215319e01373119f12342/src/abstractmcmc.jl#L35-L41">source</a></section></article><p>Luckily, implementing these is quite easy:</p><pre><code class="language-julia hljs">using MCMCTempering |
| 62 | +plot!(size=figsize)</code></pre><img src="ba7467f1.svg" alt="Example block output"/><p>Notice how <code>chain</code> has zero probability mass in the left-most component of the mixture!</p><p>Let's instead try to use a <em>tempered</em> version of <code>RWMH</code>. <em>But</em> before we do that, we need to make sure that AdvancedMH.jl is compatible with MCMCTempering.jl.</p><p>To do that we need to implement two methods. First we need to tell MCMCTempering how to extract the parameters, and potentially the log-probabilities, from a <code>AdvancedMH.Transition</code>:</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="MCMCTempering.getparams_and_logprob-getting-started" href="#MCMCTempering.getparams_and_logprob-getting-started"><code>MCMCTempering.getparams_and_logprob</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">getparams_and_logprob([model, ]state)</code></pre><p>Return a vector of parameters from the <code>state</code>.</p><p>See also: <a href="../api/#MCMCTempering.setparams_and_logprob!!"><code>setparams_and_logprob!!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/d9bc809cdc1c99b86e4adf46b31adf5e46338920/src/abstractmcmc.jl#L26-L32">source</a></section></article><p>And similarly, we need a way to <em>update</em> the parameters and the log-probabilities of a <code>AdvancedMH.Transition</code>:</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="MCMCTempering.setparams_and_logprob!!-getting-started" href="#MCMCTempering.setparams_and_logprob!!-getting-started"><code>MCMCTempering.setparams_and_logprob!!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">setparams_and_logprob!!([model, ]state, params)</code></pre><p>Set the parameters in the state to <code>params</code>, possibly mutating if it makes sense.</p><p>See also: <a href="../api/#MCMCTempering.getparams_and_logprob"><code>getparams_and_logprob</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/d9bc809cdc1c99b86e4adf46b31adf5e46338920/src/abstractmcmc.jl#L35-L41">source</a></section></article><p>Luckily, implementing these is quite easy:</p><pre><code class="language-julia hljs">using MCMCTempering |
63 | 63 |
|
64 | 64 | MCMCTempering.getparams_and_logprob(transition::AdvancedMH.Transition) = transition.params, transition.lp |
65 | 65 | function MCMCTempering.setparams_and_logprob!!(transition::AdvancedMH.Transition, params, lp) |
|
136 | 136 | end |
137 | 137 | density!(chain_tempered_all[1], color="green", size=figsize) |
138 | 138 | plot!(size=figsize)</code></pre><p>Works like a charm!</p><p><em>But</em> we're recomputing both the logdensity and the gradient of the logdensity upon every <a href="../api/#MCMCTempering.setparams_and_logprob!!"><code>MCMCTempering.setparams_and_logprob!!</code></a> above! This seems wholly unnecessary in the tempering case, since</p><p class="math-container">\[\pi_{\beta_1}(x) = \pi(x)^{\beta_1} = \big( \pi(x)^{\beta_2} \big)^{\beta_1 / \beta_2} = \pi_{\beta_2}^{\beta_1 / \beta_2}\]</p><p>i.e. if <code>model</code> in the above is tempered with <span>$\beta_1$</span> and the <code>params</code> are coming from a model with <span>$\beta_2$</span>, we can could just compute it as</p><pre><code class="language-julia hljs">(β_1 / β_2) * logprob</code></pre><p>and similarly for the gradient! Luckily, it's possible to tell MCMCTempering that this should be done by overloading the <a href="#MCMCTempering.state_from"><code>MCMCTempering.state_from</code></a> method. In particular, we'll specify that when we're working with two models of type <a href="../api/#MCMCTempering.TemperedLogDensityProblem"><code>MCMCTempering.TemperedLogDensityProblem</code></a> and two states of type <code>AdvancedHMC.HMCState</code>, then we can just re-use scale the logdensity and gradient computation from the <a href="#MCMCTempering.state_from"><code>MCMCTempering.state_from</code></a> to get the quantities we want, thus avoiding unnecessary computations:</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="MCMCTempering.state_from" href="#MCMCTempering.state_from"><code>MCMCTempering.state_from</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">state_from(model_source, state_target, state_source) |
139 | | -state_from(model_source, model_target, state_target, state_source)</code></pre><p>Return a new state similar to <code>state_target</code> but updated from <code>state_source</code>, which could be a different type of state.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/c7dde12b17a62981062215319e01373119f12342/src/abstractmcmc.jl#L44-L50">source</a></section></article><pre><code class="language-julia hljs">using AbstractMCMC: AbstractMCMC |
| 139 | +state_from(model_source, model_target, state_target, state_source)</code></pre><p>Return a new state similar to <code>state_target</code> but updated from <code>state_source</code>, which could be a different type of state.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/TuringLang/MCMCTempering.jl/blob/d9bc809cdc1c99b86e4adf46b31adf5e46338920/src/abstractmcmc.jl#L44-L50">source</a></section></article><pre><code class="language-julia hljs">using AbstractMCMC: AbstractMCMC |
140 | 140 |
|
141 | 141 | function MCMCTempering.state_from( |
142 | 142 | # AdvancedHMC.jl works with `LogDensityModel`, and by default `AbstractMCMC` will wrap |
|
183 | 183 | density!(chain_tempered, color="green", alpha=inv(sqrt(length(chain_tempered_all)))) |
184 | 184 | end |
185 | 185 | density!(chain_tempered_all[1], color="green", size=figsize) |
186 | | -plot!(size=figsize)</code></pre></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../">« Home</a><a class="docs-footer-nextpage" href="../api/">API »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="auto">Automatic (OS)</option><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option><option value="catppuccin-latte">catppuccin-latte</option><option value="catppuccin-frappe">catppuccin-frappe</option><option value="catppuccin-macchiato">catppuccin-macchiato</option><option value="catppuccin-mocha">catppuccin-mocha</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> version 1.7.0 on <span class="colophon-date" title="Monday 7 October 2024 13:51">Monday 7 October 2024</span>. Using Julia version 1.10.5.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html> |
| 186 | +plot!(size=figsize)</code></pre></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../">« Home</a><a class="docs-footer-nextpage" href="../api/">API »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="auto">Automatic (OS)</option><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option><option value="catppuccin-latte">catppuccin-latte</option><option value="catppuccin-frappe">catppuccin-frappe</option><option value="catppuccin-macchiato">catppuccin-macchiato</option><option value="catppuccin-mocha">catppuccin-mocha</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> version 1.7.0 on <span class="colophon-date" title="Monday 7 October 2024 15:38">Monday 7 October 2024</span>. Using Julia version 1.10.5.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html> |
0 commit comments