-
Notifications
You must be signed in to change notification settings - Fork 18
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
Implementation of Robust Adaptive Metropolis #106
Conversation
Maybe you could have a look at this @devmotion ?:) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Code looks mostly good, I just had a comments on documentation and testing. Can't comment on the maths in adapt
since I haven't read the paper/understood the theory.
Added proper testing @mhauru ; thanks for the review:) |
src/RobustAdaptiveMetropolis.jl
Outdated
d = LogDensityProblems.dimension(f) | ||
|
||
# Initial parameter state. | ||
x = initial_params === nothing ? rand(rng, d) : initial_params |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Maybe ensure consistent types here as well?
x = initial_params === nothing ? rand(rng, d) : initial_params | |
x = initial_params === nothing ? rand(rng, eltype(sampler.γ), d) : initial_params |
By the way, rand(rng, d)
doesn't seem a good choice in general? The algorithm requires that you start with a point in the support of the target distribution and it's not clear if the target density is zero for this point. I wonder if it requires something like https://github.com/TuringLang/EllipticalSliceSampling.jl/blob/3296ae3566d329207875216837e65eeec3b809b2/src/interface.jl#L20-L29 in EllipticalSliceSampling.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
But this is using a RWMH as the main kernel, so IMO we're already assuming unconstrained support for this to be valid
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
And more generally, happy to deal with better initialisation, buuuuut prefer to do this in a separate PR as I'm imagining this will require some discussion
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
so IMO we're already assuming unconstrained support
But why prefer rand
over randn
in that case?
But I agree, probably this question (rand
/randn
/dedicated API) should be addressed separately.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
But why prefer rand over randn in that case?
Just did it quickly because I have some vague memory that it's generally preferred to do initialisation in a box near 0 for most of the linking transformations (believe this is the moticvation behind SampleFromUniform
in DPPL, though I think it technically initialses from a cube centered on 0?).
But it's w/e to me here; we need a better way in general for this, so I'll just change it to randn
👍
Co-authored-by: David Widmann <[email protected]>
Did most of what you suggested @devmotion; see if you're okay with it now:) I don't have much time these days, so I would strongly suggest that either we merge it as is or someone else takes over this PR 👍 |
Co-authored-by: David Widmann <[email protected]>
Co-authored-by: David Widmann <[email protected]>
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
"Request changes" because the quadruple backtick seems worth fixing, happy to merge once that is done. I left a couple of other comments too, but they are optional to attend to.
src/RobustAdaptiveMetropolis.jl
Outdated
|
||
julia> norm(cov(Array(chain)) - [1.0 0.5; 0.5 1.0]) < 0.2 | ||
true | ||
``` | ||
```` |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
That looks like one backtick too many. PS. How do I make a GitHub suggestion for this?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
````suggestion
```
````
if you have N backticks surround them with N+1 😄
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
On second thought, I just pushed a fix for this myself.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks @penelopeysm. :) I should have guessed that the universe would have been organised in such a way that the fix to quadruple backticks would be more quadruple backticks.
src/RobustAdaptiveMetropolis.jl
Outdated
initial_params=nothing, | ||
kwargs... | ||
sampler::RobustAdaptiveMetropolis; | ||
initial_params = nothing, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Has JuliaFormatter been run? Just wondering because I thought the style guide was to not have whitespace around kwarg defaults. Not important though, can be left as is.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Fixed the backticks myself, the rest are optional.
Co-authored-by: Markus Hauru <[email protected]>
Thanks @mhauru and @devmotion ! |
This PR introduces a simple implementation of Robust Adaptive Metropolis (RAM) that I have lying around.
This makes use of the more recent
AbstractMCMC.step_warmup
, which is meant to be used for exactly adaptation scenarios like this. The number of adaptation steps is therefore specified bynum_warmup
passed tosample
.Note that there have been a few previous attempts at adding adaptive MH samplers to AdvancdeMH.jl, but they've all "died out" as the complexity of adding such functionality to the existing "interface" is somewhat complicated. I would suggest we merge this as is (as I don't have enough time to push this to 100%), and then we later PRs can deal with either extending this to be more general and / or integrate it with the existing "interface".
Related: #91 #39 #57