Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add more examples to docs #10

Open
1 of 3 tasks
sethaxen opened this issue Oct 31, 2021 · 4 comments
Open
1 of 3 tasks

Add more examples to docs #10

sethaxen opened this issue Oct 31, 2021 · 4 comments

Comments

@sethaxen
Copy link
Member

sethaxen commented Oct 31, 2021

We should add at least the following examples to the docs:

@joshualeond
Copy link

Exciting to see this work being done in Julia! Do you think it would make sense to include an example of how you could use the results of this package to initialize MCMC in Turing/Soss/Stan/etc?

@sethaxen
Copy link
Member Author

sethaxen commented May 4, 2022

Exciting to see this work being done in Julia! Do you think it would make sense to include an example of how you could use the results of this package to initialize MCMC in Turing/Soss/Stan/etc?

Sorry @joshualeond for the late reply! Yes, this is a great idea. We now have a Turing example, and a Soss one is planned. A Stan example would require one of the Stan interface packages exposing a log density function and its gradient, and I haven't found this functionality in any of the packages.

@burtonjosh
Copy link
Contributor

I'm not sure where you think these should go in the docs, but I've got a donut example, using $\text{logp\_donut}(x) = (\text{norm}(x)-3)^2/0.05$, which gives results like this:
image
and a simple multimodal example using a mixture of two multivariate normal distributions:
MixtureModel([MvNormal(zeros(2),0.1*ones(2)),MvNormal(ones(2),0.1*ones(2))]),
which works nicely
image
I found a more difficult looking multimodal distribution in this paper (see section 4.1) defined by the following hierarchical model
$$ \theta \sim p(\theta) = \mathcal{N}(0,\sigma_p^2\mathbb{I}_2),$$

$$ p(y\mid\theta) = \sum_{i=1}^2\frac{1}{2}\mathcal{N}(P(\theta_i),\sigma_l^2), $$

where $\sigma_l = \sigma_p = 1/4$ and $P(x) = (0.6-x)(-0.6-x)$, but I'm not sure how to figure out the log density. If you think this is worth trying and could help me with deriving the log density, I'd be happy to test on this distribution as well.

@sethaxen
Copy link
Member Author

sethaxen commented May 24, 2022

@burtonjosh these look very nice, and I'd be happy to include them as examples! I think a new docs section "More examples" with a subpage for each would be a good way to go. Those pages could be pretty simple, with a couple sentences explaining why the example is interesting (e.g. "The donut is a challenging distribution that is impossible to fit well with single-path Pathfinder."), the formula of the model, definition of logp, running Pathfinder, and viewing draws over a marginal contour map (if possible) or verses exact draws with HMC.

The logp of the multimodal distribution in the paper you referenced could be implemented something like:

using Distributions, LinearAlgebra

struct MultimodalModel{T,Y,P1,P2}
    σₚ::T
    σₗ::T
    y::Y
    p1::P1
    p2::P2
end

function (m::MultimodalModel)(θ)
    (; σₚ, σₗ, y, p1, p2) = m
    lp = logpdf(MvNormal(σₚ^2 * I(2)), θ)
    d = MixtureModel(Normal.([p1(θ[1]), p2(θ[2])], σₗ))
    lp += loglikelihood(d, y)
    return lp
end

Another multimodal distribution that might be worth looking into is a decaying periodic function in 2 dimensions. e.g. something like

using Plots
logp(θ) = prod(x -> sin(5x), θ) - sum(abs2, θ)/2
contour(-3:0.01:3, -3:0.01:3, exp  logp  Base.vect)

tmp

The first option is better for demonstrating performance with well-separated modes. The second option is better for demonstrating performance when the distribution is very bumpy, the modes are not well-separated, and each mode is non-normal.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants