Skip to content

Commit

Permalink
Added tests to Multinest
Browse files Browse the repository at this point in the history
addresses #282
  • Loading branch information
ben18785 committed Mar 12, 2021
1 parent 1996e4b commit a9eed9f
Show file tree
Hide file tree
Showing 3 changed files with 137 additions and 53 deletions.
33 changes: 33 additions & 0 deletions examples/sampling/nested-multinest.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -1306,6 +1306,39 @@
"pints.plot.pairwise(samples, kde=True, ref_parameters=real_parameters)\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"ename": "LinAlgError",
"evalue": "Singular matrix",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mLinAlgError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m<ipython-input-17-7f65263123d1>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0mpoints\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0marray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;36m0.5\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m0.5\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mcov\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcov\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtranspose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpoints\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mcov_inv\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlinalg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minv\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcov\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;32m<__array_function__ internals>\u001b[0m in \u001b[0;36minv\u001b[0;34m(*args, **kwargs)\u001b[0m\n",
"\u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.9/lib/python3.9/site-packages/numpy/linalg/linalg.py\u001b[0m in \u001b[0;36minv\u001b[0;34m(a)\u001b[0m\n\u001b[1;32m 543\u001b[0m \u001b[0msignature\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m'D->D'\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0misComplexType\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mt\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32melse\u001b[0m \u001b[0;34m'd->d'\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 544\u001b[0m \u001b[0mextobj\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mget_linalg_error_extobj\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0m_raise_linalgerror_singular\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 545\u001b[0;31m \u001b[0mainv\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_umath_linalg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minv\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msignature\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0msignature\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mextobj\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mextobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 546\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mwrap\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mainv\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mastype\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mresult_t\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcopy\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mFalse\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 547\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.9/lib/python3.9/site-packages/numpy/linalg/linalg.py\u001b[0m in \u001b[0;36m_raise_linalgerror_singular\u001b[0;34m(err, flag)\u001b[0m\n\u001b[1;32m 86\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 87\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m_raise_linalgerror_singular\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0merr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mflag\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 88\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mLinAlgError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"Singular matrix\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 89\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 90\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m_raise_linalgerror_nonposdef\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0merr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mflag\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;31mLinAlgError\u001b[0m: Singular matrix"
]
}
],
"source": [
"points = np.array([[0, 0], [0.5, 0.5], [1, 1]])\n",
"cov = np.cov(np.transpose(points))\n",
"cov_inv = np.linalg.inv(cov)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
Expand Down
1 change: 1 addition & 0 deletions pints/_nested/_multinest.py
Original file line number Diff line number Diff line change
Expand Up @@ -178,6 +178,7 @@ def __init__(self, log_prior):

self._convert_to_unit_cube = log_prior.convert_to_unit_cube
self._convert_from_unit_cube = log_prior.convert_from_unit_cube
self._ellipsoid_tree = None

def ask(self, n_points):
"""
Expand Down
156 changes: 103 additions & 53 deletions pints/tests/test_nested_multinest.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,58 +21,63 @@
unittest.TestCase.assertRaisesRegex = unittest.TestCase.assertRaisesRegexp


# class TestMultiNestSampler(unittest.TestCase):
# """
# Unit (not functional!) tests for :class:`MultinestSampler`.
# """
#
# @classmethod
# def setUpClass(cls):
# """ Prepare for the test. """
# # Create toy model
# model = pints.toy.LogisticModel()
# cls.real_parameters = [0.015, 500]
# times = np.linspace(0, 1000, 1000)
# values = model.simulate(cls.real_parameters, times)
#
# # Add noise
# np.random.seed(1)
# cls.noise = 10
# values += np.random.normal(0, cls.noise, values.shape)
# cls.real_parameters.append(cls.noise)
#
# # Create an object with links to the model and time series
# problem = pints.SingleOutputProblem(model, times, values)
#
# # Create a uniform prior over both the parameters and the new noise
# # variable
# cls.log_prior = pints.UniformLogPrior(
# [0.01, 400],
# [0.02, 600]
# )
#
# # Create a log-likelihood
# cls.log_likelihood = pints.GaussianKnownSigmaLogLikelihood(
# problem, cls.noise)
#
# def test_getters_and_setters(self):
# # tests various get() and set() methods.
# controller = pints.NestedController(self.log_likelihood,
# self.log_prior,
# method=pints.MultinestSampler)
# self.assertEqual(controller.sampler().f_s_threshold(), 1.1)
# controller.sampler().set_f_s_threshold(4)
# self.assertEqual(controller.sampler().f_s_threshold(), 4)
# self.assertRaises(ValueError, controller.sampler().set_f_s_threshold,
# 0.5)
#
# def test_runs(self):
# # tests that sampler runs
# sampler = pints.NestedController(self.log_likelihood, self.log_prior,
# method=pints.MultinestSampler)
# sampler.set_iterations(100)
# sampler.set_log_to_screen(False)
# sampler.run()
class TestMultiNestSampler(unittest.TestCase):
"""
Unit (not functional!) tests for :class:`MultinestSampler`.
"""

@classmethod
def setUpClass(cls):
""" Prepare for the test. """
# Create toy model
model = pints.toy.LogisticModel()
cls.real_parameters = [0.015, 500]
times = np.linspace(0, 1000, 1000)
values = model.simulate(cls.real_parameters, times)

# Add noise
np.random.seed(1)
cls.noise = 10
values += np.random.normal(0, cls.noise, values.shape)
cls.real_parameters.append(cls.noise)

# Create an object with links to the model and time series
problem = pints.SingleOutputProblem(model, times, values)

# Create a uniform prior over both the parameters and the new noise
# variable
cls.log_prior = pints.UniformLogPrior(
[0.01, 400],
[0.02, 600]
)

# Create a log-likelihood
cls.log_likelihood = pints.GaussianKnownSigmaLogLikelihood(
problem, cls.noise)

def test_getters_and_setters(self):
# tests various get() and set() methods.
controller = pints.NestedController(self.log_likelihood,
self.log_prior,
method=pints.MultinestSampler)
self.assertEqual(controller.sampler().f_s_threshold(), 1.1)
controller.sampler().set_f_s_threshold(4)
self.assertEqual(controller.sampler().f_s_threshold(), 4)
self.assertRaises(ValueError, controller.sampler().set_f_s_threshold,
0.5)

def test_runs(self):
# tests that sampler runs
controller = pints.NestedController(self.log_likelihood,
self.log_prior,
method=pints.MultinestSampler)
controller.set_iterations(450)
controller.set_log_to_screen(False)
controller.run()

# test getting ellipsoid tree post-run
et = controller.sampler().ellipsoid_tree()
self.assertTrue(et.n_leaf_ellipsoids() >= 1)


class TestEllipsoidTree(unittest.TestCase):
Expand Down Expand Up @@ -121,14 +126,28 @@ def test_getters(self):

# leaf nodes
self.assertTrue(tree.n_leaf_ellipsoids() >= 1)
leaves = tree.leaf_ellipsoids()
leaves = np.copy(tree.leaf_ellipsoids())
self.assertEqual(tree.n_leaf_ellipsoids(), len(leaves))
[self.assertTrue(isinstance(x, Ellipsoid)) for x in leaves]

# bounding ellipsoid
ellipsoid = tree.ellipsoid()
self.assertTrue(isinstance(ellipsoid, Ellipsoid))

# Check a given draw is within bounding ellipsoids
self.assertTrue(tree.count_within_leaf_ellipsoids(
self.draws[0]) >= 1)
self.assertTrue(tree.count_within_leaf_ellipsoids(
[-100, -100]) == 0)

# check f_s
self.assertTrue(tree.f_s() >= 1)

# updating triggers changes?
f_s_old = np.copy(tree.f_s())
tree.update_leaf_ellipsoids(1)
self.assertFalse(tree.f_s() == f_s_old)

def test_calculations(self):
# tests the calculations done

Expand All @@ -152,7 +171,38 @@ def test_calculations(self):
ellipsoid = Ellipsoid.minimum_volume_ellipsoid(draws)
self.assertRaises(ValueError, ellipsoidTree.vsk, ellipsoid)

# hk
ellipsoid = ellipsoidTree.ellipsoid()
d = Ellipsoid.mahalanobis_distance(self.draws[0],
ellipsoid.weight_matrix(),
ellipsoid.centroid())
V_S_k = 1
self.assertEqual(ellipsoid.volume() * d / V_S_k,
ellipsoidTree.h_k(self.draws[0], ellipsoid, V_S_k))

def test_sampling(self):
# tests sampling from leaf ellipsoids

# check all samples in leaf ellipsoids
tree = EllipsoidTree(self.draws, 100)
samples = tree.sample_leaf_ellipsoids(100)
for sample in samples:
self.assertTrue(tree.count_within_leaf_ellipsoids(
sample) >= 1)

def test_splitting(self):
# tests splitting and enlarging ellipsoids

# singular matrix error?
tree = EllipsoidTree(self.draws, 100)
points = np.array([[0, 0], [0.5, 0.5], [1, 1]])
assignments = np.array([0, 0, 0])
_, _, bool = tree.split_ellipsoids(points, assignments, 0)
self.assertFalse(bool)

# enlarging
ellipsoid = tree.ellipsoid()
tree.compare_enlarge(ellipsoid, 20)


if __name__ == '__main__':
Expand Down

0 comments on commit a9eed9f

Please sign in to comment.