diff --git a/bilby/core/prior.py b/bilby/core/prior.py
index 2bac06180c989a62b96e22b2c8c94721f54d1440..3692fa899b92db300d1e16c0beb5b6e7daf1e609 100644
--- a/bilby/core/prior.py
+++ b/bilby/core/prior.py
@@ -423,7 +423,8 @@ class Prior(object):
maximum: float, optional
Maximum of the domain, default=np.inf
periodic_boundary: bool, optional
- Whether or not the boundary condition is periodic. Not available in all samplers.
+ Whether or not the boundary condition is periodic.
+ Currently implemented in cpnest, dynesty and pymultinest.
"""
self.name = name
self.latex_label = latex_label
@@ -542,8 +543,8 @@ class Prior(object):
-------
ValueError: If val is not between 0 and 1
"""
- val = np.atleast_1d(val)
- tests = (val < 0) + (val > 1)
+ valarray = np.atleast_1d(val)
+ tests = (valarray < 0) + (valarray > 1)
if np.any(tests):
raise ValueError("Number to be rescaled should be in [0, 1]")
@@ -699,7 +700,7 @@ class DeltaFunction(Prior):
-------
float: Rescaled probability, equivalent to peak
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return self.peak * val ** 0
def prob(self, val):
@@ -761,7 +762,7 @@ class PowerLaw(Prior):
-------
float: Rescaled probability
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
if self.alpha == -1:
return self.minimum * np.exp(val * np.log(self.maximum / self.minimum))
else:
@@ -834,7 +835,7 @@ class Uniform(Prior):
periodic_boundary=periodic_boundary)
def rescale(self, val):
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return self.minimum + val * (self.maximum - self.minimum)
def prob(self, val):
@@ -938,7 +939,7 @@ class SymmetricLogUniform(Prior):
-------
float: Rescaled probability
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
if val < 0.5:
return -self.maximum * np.exp(-2 * val * np.log(self.maximum / self.minimum))
elif val > 0.5:
@@ -1006,7 +1007,7 @@ class Cosine(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
norm = 1 / (np.sin(self.maximum) - np.sin(self.minimum))
return np.arcsin(val / norm + np.sin(self.minimum))
@@ -1054,7 +1055,7 @@ class Sine(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
norm = 1 / (np.cos(self.minimum) - np.cos(self.maximum))
return np.arccos(np.cos(self.minimum) - val / norm)
@@ -1102,7 +1103,7 @@ class Gaussian(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return self.mu + erfinv(2 * val - 1) * 2 ** 0.5 * self.sigma
def prob(self, val):
@@ -1195,7 +1196,7 @@ class TruncatedGaussian(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return erfinv(2 * val * self.normalisation + erf(
(self.minimum - self.mu) / 2 ** 0.5 / self.sigma)) * 2 ** 0.5 * self.sigma + self.mu
@@ -1324,7 +1325,7 @@ class LogNormal(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return scipy.stats.lognorm.ppf(val, self.sigma, scale=np.exp(self.mu))
def prob(self, val):
@@ -1397,7 +1398,7 @@ class Exponential(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
return scipy.stats.expon.ppf(val, scale=self.mu)
def prob(self, val):
@@ -1458,7 +1459,7 @@ class StudentT(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
# use scipy distribution percentage point function (ppf)
return scipy.stats.t.ppf(val, self.df, loc=self.mu, scale=self.scale)
@@ -1526,7 +1527,7 @@ class Beta(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
# use scipy distribution percentage point function (ppf)
return self._dist.ppf(val)
@@ -1645,7 +1646,7 @@ class Logistic(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
# use scipy distribution percentage point function (ppf)
return scipy.stats.logistic.ppf(val, loc=self.mu, scale=self.scale)
@@ -1702,7 +1703,7 @@ class Cauchy(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
# use scipy distribution percentage point function (ppf)
return scipy.stats.cauchy.ppf(val, loc=self.alpha, scale=self.beta)
@@ -1785,7 +1786,7 @@ class Gamma(Prior):
This maps to the inverse CDF. This has been analytically solved for this case.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
# use scipy distribution percentage point function (ppf)
return scipy.stats.gamma.ppf(val, self.k, loc=0., scale=self.theta)
@@ -1915,7 +1916,7 @@ class Interped(Prior):
This maps to the inverse CDF. This is done using interpolation.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
rescaled = self.inverse_cumulative_distribution(val)
if rescaled.shape == ():
rescaled = float(rescaled)
@@ -2079,7 +2080,7 @@ class FermiDirac(Prior):
.. [1] M. Pitkin, M. Isi, J. Veitch & G. Woan, `arXiv:1705.08978v1
<https:arxiv.org/abs/1705.08978v1>`_, 2017.
"""
- Prior.test_valid_for_rescaling(val)
+ self.test_valid_for_rescaling(val)
inv = (-np.exp(-1. * self.r) + (1. + np.exp(self.r))**-val +
np.exp(-1. * self.r) * (1. + np.exp(self.r))**-val)
diff --git a/bilby/core/sampler/base_sampler.py b/bilby/core/sampler/base_sampler.py
index daa204293290a26dd2017ac32b52621a80a63a26..d1c0274b9630db21e1326ea374a1596ca50dc37d 100644
--- a/bilby/core/sampler/base_sampler.py
+++ b/bilby/core/sampler/base_sampler.py
@@ -494,7 +494,8 @@ class NestedSampler(Sampler):
class MCMCSampler(Sampler):
- nwalkers_equiv_kwargs = ['nwalker', 'nwalkers', 'draws']
+ nwalkers_equiv_kwargs = ['nwalker', 'nwalkers', 'draws', 'Niter']
+ nburn_equiv_kwargs = ['burn', 'nburn']
def print_nburn_logging_info(self):
""" Prints logging info as to how nburn was calculated """
diff --git a/bilby/core/sampler/dynesty.py b/bilby/core/sampler/dynesty.py
index e2787895f0cf220896b5e6c7aa232f046b33058e..94b9f6e689ad4573d4c7fcb5a80b5ca504619f16 100644
--- a/bilby/core/sampler/dynesty.py
+++ b/bilby/core/sampler/dynesty.py
@@ -71,7 +71,7 @@ class Dynesty(NestedSampler):
If true, resume run from checkpoint (if available)
"""
default_kwargs = dict(bound='multi', sample='rwalk',
- verbose=True,
+ verbose=True, periodic=None,
check_point_delta_t=600, nlive=500,
first_update=None,
npdim=None, rstate=None, queue_size=None, pool=None,
@@ -95,6 +95,7 @@ class Dynesty(NestedSampler):
self.n_check_point = n_check_point
self.check_point = check_point
self.resume = resume
+ self._apply_dynesty_boundaries()
if self.n_check_point is None:
# If the log_likelihood_eval_time is not calculable then
# check_point is set to False.
@@ -173,6 +174,15 @@ class Dynesty(NestedSampler):
sys.stderr.write(print_str)
sys.stderr.flush()
+ def _apply_dynesty_boundaries(self):
+ if self.kwargs['periodic'] is None:
+ logger.debug("Setting periodic boundaries for keys:")
+ self.kwargs['periodic'] = []
+ for ii, key in enumerate(self.search_parameter_keys):
+ if self.priors[key].periodic_boundary:
+ self.kwargs['periodic'].append(ii)
+ logger.debug(" {}".format(key))
+
def run_sampler(self):
import dynesty
self.sampler = dynesty.NestedSampler(
@@ -397,3 +407,23 @@ class Dynesty(NestedSampler):
self.result.log_evidence = np.nan
self.result.log_evidence_err = np.nan
return self.result
+
+ def prior_transform(self, theta):
+ """ Prior transform method that is passed into the external sampler.
+ cube we map this back to [0, 1].
+
+ Parameters
+ ----------
+ theta: list
+ List of sampled values on a unit interval
+
+ Returns
+ -------
+ list: Properly rescaled sampled values
+
+ Notes
+ -----
+ Since dynesty allows periodic parameters to wander outside the unit
+ """
+ theta = np.mod(theta, 1)
+ return self.priors.rescale(self._search_parameter_keys, theta)
diff --git a/bilby/core/sampler/ptmcmc.py b/bilby/core/sampler/ptmcmc.py
index bdb13962f75d5ece4977d622651e6d5d69bc6255..d8c93d7270b07b3ab824e12e064130dbfa19fe51 100644
--- a/bilby/core/sampler/ptmcmc.py
+++ b/bilby/core/sampler/ptmcmc.py
@@ -78,7 +78,7 @@ class PTMCMCSampler(MCMCSampler):
def _translate_kwargs(self, kwargs):
if 'Niter' not in kwargs:
- for equiv in self.nsteps_equiv_kwargs:
+ for equiv in self.nwalkers_equiv_kwargs:
if equiv in kwargs:
kwargs['Niter'] = kwargs.pop(equiv)
if 'burn' not in kwargs:
diff --git a/bilby/core/sampler/pymultinest.py b/bilby/core/sampler/pymultinest.py
index 7472e900e61a051dcea65eb5b1a791ea5e11bf16..0533e3bece39084c8f0a29c92e3a7914f216673f 100644
--- a/bilby/core/sampler/pymultinest.py
+++ b/bilby/core/sampler/pymultinest.py
@@ -46,6 +46,14 @@ class Pymultinest(NestedSampler):
context=0, write_output=True, log_zero=-1e100,
max_iter=0, init_MPI=False, dump_callback=None)
+ def __init__(self, likelihood, priors, outdir='outdir', label='label', use_ratio=False, plot=False,
+ skip_import_verification=False, **kwargs):
+ NestedSampler.__init__(self, likelihood=likelihood, priors=priors, outdir=outdir, label=label,
+ use_ratio=use_ratio, plot=plot,
+ skip_import_verification=skip_import_verification,
+ **kwargs)
+ self._apply_multinest_boundaries()
+
def _translate_kwargs(self, kwargs):
if 'n_live_points' not in kwargs:
for equiv in self.npoints_equiv_kwargs:
@@ -75,6 +83,15 @@ class Pymultinest(NestedSampler):
self.kwargs['outputfiles_basename'])
NestedSampler._verify_kwargs_against_default_kwargs(self)
+ def _apply_multinest_boundaries(self):
+ if self.kwargs['wrapped_params'] is None:
+ self.kwargs['wrapped_params'] = []
+ for param, value in self.priors.items():
+ if value.periodic_boundary:
+ self.kwargs['wrapped_params'].append(1)
+ else:
+ self.kwargs['wrapped_params'].append(0)
+
def run_sampler(self):
import pymultinest
self._verify_kwargs_against_default_kwargs()
diff --git a/test/sampler_test.py b/test/sampler_test.py
index 314605470f75034b3b6eba3384e92a7d66f81a64..7c1024e1ce6d3140e6a7d52f67d5dc3e37b63db6 100644
--- a/test/sampler_test.py
+++ b/test/sampler_test.py
@@ -132,7 +132,9 @@ class TestDynesty(unittest.TestCase):
def setUp(self):
self.likelihood = MagicMock()
- self.priors = dict()
+ self.priors = bilby.core.prior.PriorDict()
+ self.priors['a'] = bilby.core.prior.Prior(periodic_boundary=True)
+ self.priors['b'] = bilby.core.prior.Prior(periodic_boundary=False)
self.sampler = bilby.core.sampler.Dynesty(self.likelihood, self.priors,
outdir='outdir', label='label',
use_ratio=False, plot=False,
@@ -144,7 +146,7 @@ class TestDynesty(unittest.TestCase):
del self.sampler
def test_default_kwargs(self):
- expected = dict(bound='multi', sample='rwalk', verbose=True,
+ expected = dict(bound='multi', sample='rwalk', periodic=None, verbose=True,
check_point_delta_t=600, nlive=500, first_update=None,
npdim=None, rstate=None, queue_size=None, pool=None,
use_pool=None, live_points=None, logl_args=None, logl_kwargs=None,
@@ -152,12 +154,18 @@ class TestDynesty(unittest.TestCase):
enlarge=None, bootstrap=None, vol_dec=0.5, vol_check=2.0,
facc=0.5, slices=5, dlogz=0.1, maxiter=None, maxcall=None,
logl_max=np.inf, add_live=True, print_progress=True, save_bounds=True,
- walks=0, update_interval=300, print_func='func')
+ walks=10, update_interval=300, print_func='func')
self.sampler.kwargs['print_func'] = 'func' # set this manually as this is not testable otherwise
+ self.assertListEqual([0], self.sampler.kwargs['periodic']) # Check this separately
+ self.sampler.kwargs['periodic'] = None # The dict comparison can't handle lists
+ for key in self.sampler.kwargs.keys():
+ print(key)
+ print(expected[key])
+ print(self.sampler.kwargs[key])
self.assertDictEqual(expected, self.sampler.kwargs)
def test_translate_kwargs(self):
- expected = dict(bound='multi', sample='rwalk', verbose=True,
+ expected = dict(bound='multi', sample='rwalk', periodic=[0], verbose=True,
check_point_delta_t=600, nlive=250, first_update=None,
npdim=None, rstate=None, queue_size=None, pool=None,
use_pool=None, live_points=None, logl_args=None, logl_kwargs=None,
@@ -165,7 +173,7 @@ class TestDynesty(unittest.TestCase):
enlarge=None, bootstrap=None, vol_dec=0.5, vol_check=2.0,
facc=0.5, slices=5, dlogz=0.1, maxiter=None, maxcall=None,
logl_max=np.inf, add_live=True, print_progress=True, save_bounds=True,
- walks=0, update_interval=300, print_func='func')
+ walks=10, update_interval=300, print_func='func')
for equiv in bilby.core.sampler.base_sampler.NestedSampler.npoints_equiv_kwargs:
new_kwargs = self.sampler.kwargs.copy()
@@ -385,7 +393,9 @@ class TestPymultinest(unittest.TestCase):
def setUp(self):
self.likelihood = MagicMock()
- self.priors = dict()
+ self.priors = bilby.core.prior.PriorDict()
+ self.priors['a'] = bilby.core.prior.Prior(periodic_boundary=True)
+ self.priors['b'] = bilby.core.prior.Prior(periodic_boundary=False)
self.sampler = bilby.core.sampler.Pymultinest(self.likelihood, self.priors,
outdir='outdir', label='label',
use_ratio=False, plot=False,
@@ -408,6 +418,8 @@ class TestPymultinest(unittest.TestCase):
max_modes=100, mode_tolerance=-1e90, seed=-1,
context=0, write_output=True, log_zero=-1e100,
max_iter=0, init_MPI=False, dump_callback=None)
+ self.assertListEqual([1, 0], self.sampler.kwargs['wrapped_params']) # Check this separately
+ self.sampler.kwargs['wrapped_params'] = None # The dict comparison can't handle lists
self.assertDictEqual(expected, self.sampler.kwargs)
def test_translate_kwargs(self):
@@ -425,6 +437,7 @@ class TestPymultinest(unittest.TestCase):
for equiv in bilby.core.sampler.base_sampler.NestedSampler.npoints_equiv_kwargs:
new_kwargs = self.sampler.kwargs.copy()
del new_kwargs['n_live_points']
+ new_kwargs['wrapped_params'] = None # The dict comparison can't handle lists
new_kwargs[equiv] = 123
self.sampler.kwargs = new_kwargs
self.assertDictEqual(expected, self.sampler.kwargs)
@@ -444,8 +457,9 @@ class TestRunningSamplers(unittest.TestCase):
self.likelihood = bilby.likelihood.GaussianLikelihood(
self.x, self.y, self.model, self.sigma)
- self.priors = dict(
- m=bilby.core.prior.Uniform(0, 5), c=bilby.core.prior.Uniform(-2, 2))
+ self.priors = bilby.core.prior.PriorDict()
+ self.priors['m'] = bilby.core.prior.Uniform(0, 5, periodic_boundary=False)
+ self.priors['c'] = bilby.core.prior.Uniform(-2, 2, periodic_boundary=False)
bilby.core.utils.check_directory_exists_and_if_not_mkdir('outdir')
def tearDown(self):