apunkt/GPy
1
0
Fork 0
mirror of https://github.com/SheffieldML/GPy.git synced 2026-06-08 15:05:15 +02:00
Code Issues Projects Releases Packages Wiki Activity Actions

format on save

Browse source
This commit is contained in:
Martin Bubel 2023-10-10 19:51:04 +02:00
parent de670f8bcb
commit 5fde9d2edd
1 changed files with 32 additions and 31 deletions
Download patch file Download diff file Expand all files Collapse all files

63
GPy/testing/rv_transformation_tests.py
View file

@ -4,7 +4,6 @@ Test if hyperparameters in models are properly transformed.
"""
import unittest
import numpy as np
import scipy.stats as st
import GPy
@ -14,25 +13,25 @@ class TestModel(GPy.core.Model):
"""
A simple GPy model with one parameter.
"""
def __init__(self, theta=1.):
super(TestModel, self).__init__('test_model')
theta = GPy.core.Param('theta', theta)
def __init__(self, theta=1.0):
super(TestModel, self).__init__("test_model")
theta = GPy.core.Param("theta", theta)
self.link_parameter(theta)
def log_likelihood(self):
return 0.
return 0.0
class RVTransformationTestCase(unittest.TestCase):
def _test_trans(self, trans):
m = TestModel()
prior = GPy.priors.LogGaussian(.5, 0.1)
prior = GPy.priors.LogGaussian(0.5, 0.1)
m.theta.set_prior(prior)
m.theta.unconstrain()
m.theta.constrain(trans)
# The PDF of the transformed variables
p_phi = lambda phi : np.exp(-m._objective_grads(phi)[0])
p_phi = lambda phi: np.exp(-m._objective_grads(phi)[0])
# To the empirical PDF of:
theta_s = prior.rvs(1e5)
phi_s = trans.finv(theta_s)
@ -43,23 +42,25 @@ class RVTransformationTestCase(unittest.TestCase):
# The transformed PDF of phi should be this:
pdf_phi = np.array([p_phi(p) for p in phi])
# UNCOMMENT TO SEE GRAPHICAL COMPARISON
#import matplotlib.pyplot as plt
#fig, ax = plt.subplots()
#ax.hist(phi_s, normed=True, bins=100, alpha=0.25, label='Histogram')
#ax.plot(phi, kde(phi), '--', linewidth=2, label='Kernel Density Estimation')
#ax.plot(phi, pdf_phi, ':', linewidth=2, label='Transformed PDF')
#ax.set_xlabel(r'transformed $\theta$', fontsize=16)
#ax.set_ylabel('PDF', fontsize=16)
#plt.legend(loc='best')
#plt.show(block=True)
# import matplotlib.pyplot as plt
# fig, ax = plt.subplots()
# ax.hist(phi_s, normed=True, bins=100, alpha=0.25, label='Histogram')
# ax.plot(phi, kde(phi), '--', linewidth=2, label='Kernel Density Estimation')
# ax.plot(phi, pdf_phi, ':', linewidth=2, label='Transformed PDF')
# ax.set_xlabel(r'transformed $\theta$', fontsize=16)
# ax.set_ylabel('PDF', fontsize=16)
# plt.legend(loc='best')
# plt.show(block=True)
# END OF PLOT
# The following test cannot be very accurate
self.assertTrue(np.linalg.norm(pdf_phi - kde(phi)) / np.linalg.norm(kde(phi)) <= 1e-1)
self.assertTrue(
np.linalg.norm(pdf_phi - kde(phi)) / np.linalg.norm(kde(phi)) <= 1e-1
)
def _test_grad(self, trans):
np.random.seed(1234)
m = TestModel(np.random.uniform(.5, 1.5, 20))
prior = GPy.priors.LogGaussian(.5, 0.1)
m = TestModel(np.random.uniform(0.5, 1.5, 20))
prior = GPy.priors.LogGaussian(0.5, 0.1)
m.theta.set_prior(prior)
m.theta.constrain(trans)
m.randomize()
@ -70,26 +71,26 @@ class RVTransformationTestCase(unittest.TestCase):
self._test_trans(GPy.constraints.Logexp())
@unittest.skip("Gradient not checking right, @jameshensman what is going on here?")
def test_Logexp_grad(self):
def test_Logexp_grad(self):
self._test_grad(GPy.constraints.Logexp())
def test_Exponent(self):
self._test_trans(GPy.constraints.Exponent())
@unittest.skip("Gradient not checking right, @jameshensman what is going on here?")
def test_Exponent_grad(self):
self._test_grad(GPy.constraints.Exponent())
if __name__ == '__main__':
if __name__ == "__main__":
unittest.main()
quit()
m = TestModel()
prior = GPy.priors.LogGaussian(0., .9)
prior = GPy.priors.LogGaussian(0.0, 0.9)
m.theta.set_prior(prior)
# The following should return the PDF in terms of the transformed quantities
p_phi = lambda phi : np.exp(-m._objective_grads(phi)[0])
p_phi = lambda phi: np.exp(-m._objective_grads(phi)[0])
# Let's look at the transformation phi = log(exp(theta - 1))
trans = GPy.constraints.Exponent()
@ -103,14 +104,14 @@ if __name__ == '__main__':
# Transform it to the new variables
phi_s = trans.finv(theta_s)
# And draw their histogram
ax.hist(phi_s, normed=True, bins=100, alpha=0.25, label='Empirical')
ax.hist(phi_s, normed=True, bins=100, alpha=0.25, label="Empirical")
# This is to be compared to the PDF of the model expressed in terms of these new
# variables
ax.plot(phi, [p_phi(p) for p in phi], label='Transformed PDF', linewidth=2)
ax.plot(phi, [p_phi(p) for p in phi], label="Transformed PDF", linewidth=2)
ax.set_xlim(-3, 10)
ax.set_xlabel(r'transformed $\theta$', fontsize=16)
ax.set_ylabel('PDF', fontsize=16)
plt.legend(loc='best')
ax.set_xlabel(r"transformed $\theta$", fontsize=16)
ax.set_ylabel("PDF", fontsize=16)
plt.legend(loc="best")
# Now let's test the gradients
m.checkgrad(verbose=True)
# And show the plot