joshspeagle / dynesty

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

1.0.2 issue with rejecting invalid multiellipsoid region #205

Closed segasai closed 3 years ago

segasai commented 4 years ago

Hi,

while running the latest master version (with the 'rejecting multiellipsoid region' change), I'm getting this exception. I didn't have it before (with the previously released version). The error seems to occur exactly after the warnign about rejecting the multiellipsoid region.

eff(%):  0.282 | loglstar: 10019.867 < 10026.423 < 10028.004 | logz: 9949.172 +/-  0.545 | stop:  1.314]
/home/skoposov/pyenv38/lib/python3.8/site-packages/dynesty/bounding.py:596: UserWarning: Rejecting invalid MultiEllipsoid region
  warnings.warn('Rejecting invalid MultiEllipsoid region')
59036it [8:18:40,  1.97it/s, batch: 8 | bound: 787 | nc: 332 | ncall: 20754337 | eff(%):  0.282 | loglstar: 10019.867 < 10026.430 < 10028.004 | logz: 9949.172 +/-  0.545 | stop:  1.314]
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-2-4fd4a960409f> in <module>
----> 1 R=qq.sampler()

~/science/dipper/qq.py in sampler()
    359                                             sample='slice'
    360                                             )
--> 361     dsampler.run_nested(dlogz_init=0.01,
    362                         nlive_init=500,
    363                         nlive_batch=500,

~/pyenv38/lib/python3.8/site-packages/dynesty/dynamicsampler.py in run_nested(self, nlive_init, maxiter_init, maxcall_init, dlogz_init, logl_max_init, n_effective_init, nlive_batch, wt_function, wt_kwargs, maxiter_batch, maxcall_batch, maxiter, maxcall, maxbatch, n_effective, stop_function, stop_kwargs, use_stop, save_bounds, print_progress, print_func, live_points)
   1658                     # Compute our sampling bounds using the provided
   1659                     # weight function.
-> 1660                     passback = self.add_batch(nlive=nlive_batch,
   1661                                               wt_function=wt_function,
   1662                                               wt_kwargs=wt_kwargs,

~/pyenv38/lib/python3.8/site-packages/dynesty/dynamicsampler.py in add_batch(self, nlive, wt_function, wt_kwargs, maxiter, maxcall, logl_bounds, save_bounds, print_progress, print_func, stop_val)
   1767                 if logl_bounds is None:
   1768                     logl_bounds = wt_function(res, wt_kwargs)
-> 1769                 for results in self.sample_batch(nlive_new=nlive,
   1770                                                  logl_bounds=logl_bounds,
   1771                                                  maxiter=maxiter,

~/pyenv38/lib/python3.8/site-packages/dynesty/dynamicsampler.py in sample_batch(self, nlive_new, update_interval, logl_bounds, maxiter, maxcall, save_bounds)
   1207         # sample past the original bounds "for free".
   1208         for i in range(1):
-> 1209             for it, results in enumerate(self.sampler.sample(dlogz=0.,
   1210                                          logl_max=logl_max,
   1211                                          maxiter=maxiter-nlive_new-1,

~/pyenv38/lib/python3.8/site-packages/dynesty/sampler.py in sample(self, maxiter, maxcall, dlogz, logl_max, n_effective, add_live, save_bounds, save_samples)
    780             # `logl > loglstar` using the bounding distribution and sampling
    781             # method from our sampler.
--> 782             u, v, logl, nc = self._new_point(loglstar_new, logvol)
    783             ncall += nc
    784             self.ncall += nc

~/pyenv38/lib/python3.8/site-packages/dynesty/sampler.py in _new_point(self, loglstar, logvol)
    378         while True:
    379             # Get the next point from the queue
--> 380             u, v, logl, nc, blob = self._get_point_value(loglstar)
    381             ncall += nc
    382 

~/pyenv38/lib/python3.8/site-packages/dynesty/sampler.py in _get_point_value(self, loglstar)
    362         # If the queue is empty, refill it.
    363         if self.nqueue <= 0:
--> 364             self._fill_queue(loglstar)
    365 
    366         # Grab the earliest entry.

~/pyenv38/lib/python3.8/site-packages/dynesty/sampler.py in _fill_queue(self, loglstar)
    331             if self._beyond_unit_bound(loglstar):
    332                 # Propose points using the provided sampling/bounding options.
--> 333                 point, axes = self.propose_point()
    334                 evolve_point = self.evolve_point
    335             else:

~/pyenv38/lib/python3.8/site-packages/dynesty/nestedsamplers.py in propose_live(self)
    664 
    665         # Pick a random ellipsoid that encompasses `u`.
--> 666         ell_idx = ell_idxs[self.rstate.randint(nidx)]
    667 
    668         # Choose axes.

mtrand.pyx in numpy.random.mtrand.RandomState.randint()

_bounded_integers.pyx in numpy.random._bounded_integers._rand_int64()

ValueError: low >= high

I'm not sure I'll be able to provide a backtrace on this one, as it was a very long run.

joshspeagle commented 4 years ago

Not a problem. I’ll look into this and see if I can find out what’s going on here. (These bounding distribution edge cases will be the end of me...)

segasai commented 4 years ago

Thanks. I am wondering honestly if the test suite is lacking challenging distributions. I.e. 100d Rosenbrock or funnel, or some almost degenerate high-D valley. I think that would help to catch those things. And if the running time is the issue, presumably individual tests can be scheduled as separate runs on travis, so each one is shorter than 1hr

joshspeagle commented 4 years ago

You're right -- I don't have those tests in there primarily because of the long(er) run times. But given some of these problems, it looks like I should just bite the bullet and throw in a challenging high-D problem.

segasai commented 4 years ago

(not to be critical), but personally, for several challenging problems I had so far in the last couple of years, I had to always switch back to pymultinest, because of some issues with sampling/boundaries/efficiency in dynesty. So It would be good to have a hard distribution testing suite to ensure that everything works like it should.

joshspeagle commented 4 years ago

No offense taken -- dynesty was the first public package that I made when I was learning how to code properly, and that's led to problems that have come back to bite me ever since! Hopefully, with future software I will have better tests in place to catch these earlier. Sorry that the package has caused you frustration.

jvines commented 3 years ago

Just because I don't want to open yet another issue with this same error. I could provide the data and run setup if that would help, as well as a dump of the sampler after it crashes..

joshspeagle commented 3 years ago

That'd be super. Thanks!

jvines commented 3 years ago

OK. This took long enough but I sent you an email with a script and instructions on how to replicate the errors I'm getting... Hope it helps you!