wrap around dimensions

[philippeller]

We're currently using PyMultinest, but would like to switch to dynesty. The one important missing feature preventing us from doing so is wrap-around parameters. For cyclic dimensions like an azimuth angle it is necessary to be able to define that. Is this a feature that could be easily added?

[joshspeagle]

Happy to add in support for this. Could you provide an example?

[philippeller]

Definitely! We use nested sampling to explore a space in spherical coordinates (r, theta and phi). The azimuthal angle phi is cyclic in a sense that for example an angle phi + 2pi is the same as phi. Consequently, boundary conditions for drawing sampling points need to be able to wrap around.....meaning you could have an ellipse spanning for example from 6 (lower bound) to 0.5 radians (upper bound) crossing zero.

This is realized in Multinest by setting the pWrap flag for a dimension. I copied the following from the MN documentation:

Periodic Boundary Conditions:

In order to sample from parameters with periodic boundary conditions (or wraparound parameters), set pWrap[i],
where i is the index of the parameter to be wraparound, to a non-zero value. If pWrap[i] = 0, then the ith
parameter is not wraparound.

Let me know if that's enough information or not!

[philippeller]

I just checked in the MN implementation, and i think what happens (as far as my fortran reading abilities go ) that the sampler is simply allowed to walk outside the unit cube. And the paramter is reported back essentially as the modulo of the parmater to 1

[joshspeagle]

Now I'm a bit confused how this periodic behavior requires special treatment. Is the idea that you don't want to bound phi from 0 to 2pi and would prefer to sample completely unconstrained? Or that if you set the prior to be from [a, b] you then want to map back to [0, 2pi]? Because that seems as simple as just including an additional modulo operation within prior_transform.

[philippeller]

hey, the problem is not mapping the paramater (that, as you pointed out would be trivial to do in the prior) but if let's say your maximum of you function is close to the boundary, your nested sampling will need to include the points that map to the other side. So it starts being a problem once for example the hyper-ellipse starts constraining the parameter space. If you don't have that wrap around functionality you end up with two separate hyper ellispses at each boundary

[joshspeagle]

Ah, I get it now. Yes, that edge case would impede sampling by creating two modes when there really is just one. Let me think about how I want to implement this type of unit-cube shifting here and I'll try to get back to you by Monday.

[joshspeagle]

Okay, I have a local test version with wraparound working using a similar API as PyMultiNest. This functions similarly to the code you described above, where the live points are just allowed to walk off the unit cube and wrap around mod 1. The implementation currently modifies the _wrapper_function used to wrap the input prior_transform function, so the user doesn't have to change anything other than setting the appropriate flags.

On a sidenote, this turned out to be quite a bit more involved than I had expected because I had actually hard-coded unit cube checks everywhere and had to replace them with a more generalized check.

I'm finishing up some additional testing now and will try to push it by early next week.

[philippeller]

Hey, that's really cool! I'll be happy to test the code when available (but gonna be at a conference next week, so will take a bit longer)