joshspeagle / dynesty

Dynamic Nested Sampling package for computing Bayesian posteriors and evidences
https://dynesty.readthedocs.io/
MIT License
347 stars 76 forks source link

Performance versus nestle #170

Closed andrewfowlie closed 4 years ago

andrewfowlie commented 4 years ago

I was looking at the performance of a few nested sampling implementations. I expected nestle and dynesty to be very similar. Here is a comparison with the classic eggbox problem. I tried to make the nestle and dynesty settings the same for a fair test (multi-ellipsoid sampling with the same nlive etc).

import numpy as np
import time

import nestle
import dynesty

# Classic eggbox problem
tmax = 5. * np.pi

def loglike(x):
    t = tmax * x - 0.5 * tmax
    return (2. + np.cos(t[0]) * np.cos(t[1]))**5

def prior(x):
    return x

t0 = time.time()
res = nestle.sample(loglike, prior, 2, npoints=200, method='multi', update_interval=50, dlogz=0.5, enlarge=1.25)
dt = time.time() - t0
print("nestle. time = {}. like calls = {}. iterations = {}".format(dt, res.ncall, res.niter))

t0 = time.time()
sampler = dynesty.NestedSampler(loglike, prior, 2, nlive=200, bound='multi', sample="unif", update_interval=50, enlarge=1.25)
sampler.run_nested(print_progress=False, dlogz=0.5)
dt = time.time() - t0
res = sampler.results
print("dynesty. time = {}. like calls = {}. iterations = {}".format(dt, sum(res.ncall), res.niter))

I ran the above code and found (and similar on repeats)

time = 0.659329891205. like calls = 3607. iterations = 1506
time = 1.22893309593. like calls = 8692. iterations = 1576

i.e., dynesty appeared to take around twice and long and use about 2-3 times as many likelihood calls, in roughly the same number of iterations.

Are my settings equivalent? If so, what's going on here?

joshspeagle commented 4 years ago

I thought I had included a direct answer in the FAQ page of the documentation for this, but realized it was totally missing and only addressed deep in the bounding options section of the documentation and somewhat in the FAQ.

In brief, the reason for this slowdown has to do with two things. First, dynesty generally has a more conservative threshold for generating/splitting ellipsoids than nestle (where this option is set internally and can't be easily modified), leading it to generally be less efficient overall in order to be more robust. You can adjust these with the vol_dec and vol_check arguments, although I actually forget now what the nestle defaults actually are.

The second and more important reason has to do with the first_update argument. Unlike nestle or MultiNest, dynesty actually waits to construct the first bound until the sampling efficiency (averaged over the proposals) hits 10%. This helps to avoid "shredding" the live points into a ton of ellipsoids, which can happen when trying to approximate a uniform distribution of points in a cube with a bunch of overlapping ellipsoids. Disabling that option using the first_update argument should resolve most of the discrepancy, although I still expect nestle to slightly outperform dynesty on this toy problem.

I've actually added in additional text into both the FAQ and modified the bounding section of the documentation to emphasize these points better. I'll close this issue once those are pushed to the master branch.

Hope this helps! Thanks for bringing this to my attention.

andrewfowlie commented 4 years ago

Thanks Josh, that’s now very clear. PS I’ll post the full benchmarking code on GitHub gists at some point. I am making comparisons of many Nested Sampling implementations.