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

Merge branch 'devel' of github.com:SheffieldML/GPy into devel

Browse source
This commit is contained in:
Max Zwiessele 2014-11-04 12:40:27 +00:00
commit 54d852d55f
8 changed files with 190 additions and 52 deletions
Download patch file Download diff file Expand all files Collapse all files

68
GPy/core/gp.py
View file

@ -93,21 +93,51 @@ class GP(Model):
self.link_parameter(self.kern)
self.link_parameter(self.likelihood)
def set_XY(self, X=None, Y=None):
"""
Set the input / output of the model
:param X: input observations
:param Y: output observations
"""
self.update_model(False)
if Y is not None:
if self.normalizer is not None:
self.normalizer.scale_by(Y)
self.Y_normalized = ObsAr(self.normalizer.normalize(Y))
self.Y = Y
else:
self.Y = ObsAr(Y)
self.Y_normalized = self.Y
if X is not None:
if self.X in self.parameters:
# LVM models
from ..core.parameterization.variational import VariationalPosterior
if isinstance(self.X, VariationalPosterior):
assert isinstance(X, type(self.X)), "The given X must have the same type as the X in the model!"
self.unlink_parameter(self.X)
self.X = X
self.link_parameters(self.X)
else:
self.unlink_parameter(self.X)
from ..core import Param
self.X = Param('latent mean',X)
self.link_parameters(self.X)
else:
self.X = ObsAr(X)
self.update_model(True)
def set_X(self,X):
# TODO: it does not work with BGPLVM
if isinstance(X, ObsAr):
self.X = X
else:
self.X = ObsAr(X)
"""
Set the input of the model
"""
self.set_XY(X=X)
def set_Y(self,Y):
if self.normalizer is not None:
self.normalizer.scale_by(Y)
self.Y_normalized = ObsAr(self.normalizer.normalize(Y))
self.Y = Y
else:
self.Y = ObsAr(Y)
self.Y_normalized = self.Y
"""
Set the input of the model
"""
self.set_XY(Y=Y)
def parameters_changed(self):
self.posterior, self._log_marginal_likelihood, self.grad_dict = self.inference_method.inference(self.kern, self.X, self.likelihood, self.Y_normalized, self.Y_metadata)
@ -354,3 +384,17 @@ class GP(Model):
print "KeyboardInterrupt caught, calling on_optimization_end() to round things up"
self.inference_method.on_optimization_end()
raise
def infer_newX(self, Y_new, optimize=True, ):
"""
Infer the distribution of X for the new observed data *Y_new*.
:param Y_new: the new observed data for inference
:type Y_new: numpy.ndarray
:param optimize: whether to optimize the location of new X (True by default)
:type optimize: boolean
:return: a tuple containing the posterior estimation of X and the model that optimize X
:rtype: (GPy.core.parameterization.variational.VariationalPosterior or numpy.ndarray, GPy.core.Model)
"""
from ..inference.latent_function_inference.inferenceX import infer_newX
return infer_newX(self, Y_new, optimize=optimize)

6
GPy/core/parameterization/parameter_core.py
View file

@ -518,7 +518,7 @@ class Indexable(Nameable, Observable):
self.constrain_negative(warning)
elif prior.domain is _REAL:
rav_i = self._raveled_index()
assert all(all(c.domain is _REAL for c in con) for con in self.constraints.properties_for(rav_i))
assert all(all(False if c is __fixed__ else c.domain is _REAL for c in con) for con in self.constraints.properties_for(rav_i)), 'Domain of prior and constraint have to match, please unconstrain if you REALLY wish to use this prior'
def unset_priors(self, *priors):
"""
@ -824,7 +824,7 @@ class OptimizationHandlable(Indexable):
#===========================================================================
# Randomizeable
#===========================================================================
def randomize(self, rand_gen=np.random.normal, *args, **kwargs):
def randomize(self, rand_gen=None, *args, **kwargs):
"""
Randomize the model.
Make this draw from the prior if one exists, else draw from given random generator
@ -834,6 +834,8 @@ class OptimizationHandlable(Indexable):
:param float scale: scale parameter for random number generator
:param args, kwargs: will be passed through to random number generator
"""
if rand_gen is None:
rand_gen = np.random.normal
# first take care of all parameters (from N(0,1))
x = rand_gen(size=self._size_transformed(), *args, **kwargs)
updates = self.update_model()

69
GPy/inference/latent_function_inference/inferenceX.py
View file

@ -51,8 +51,18 @@ class InferenceX(Model):
self.kern.GPU(True)
from copy import deepcopy
self.posterior = deepcopy(model.posterior)
self.variational_prior = model.variational_prior.copy()
self.Z = model.Z.copy()
if hasattr(model, 'variational_prior'):
self.uncertain_input = True
self.variational_prior = model.variational_prior.copy()
else:
self.uncertain_input = False
if hasattr(model, 'inducing_inputs'):
self.sparse_gp = True
self.Z = model.Z.copy()
else:
self.sparse_gp = False
self.uncertain_input = False
self.Z = model.X.copy()
self.Y = Y
self.X = self._init_X(model, Y, init=init)
self.compute_dL()
@ -72,6 +82,8 @@ class InferenceX(Model):
dist = -2.*Y_new.dot(Y.T) + np.square(Y_new).sum(axis=1)[:,None]+ np.square(Y).sum(axis=1)[None,:]
elif init=='NCC':
dist = Y_new.dot(Y.T)
elif init=='rand':
dist = np.random.rand(Y_new.shape[0],Y.shape[0])
idx = dist.argmin(axis=1)
from ...models import SSGPLVM
@ -81,7 +93,11 @@ class InferenceX(Model):
if model.group_spike:
X.gamma.fix()
else:
X = variational.NormalPosterior(param_to_array(model.X.mean[idx]), param_to_array(model.X.variance[idx]))
if self.uncertain_input and self.sparse_gp:
X = variational.NormalPosterior(param_to_array(model.X.mean[idx]), param_to_array(model.X.variance[idx]))
else:
from ...core import Param
X = Param('latent mean',param_to_array(model.X[idx]).copy())
return X
@ -99,29 +115,42 @@ class InferenceX(Model):
else:
self.dL_dpsi2 = beta*(output_dim*self.posterior.woodbury_inv - np.einsum('md,od->mo',wv, wv))/2.
self.dL_dpsi1 = beta*np.dot(self.Y, wv.T)
self.dL_dpsi0 = -beta/2.* np.ones(self.Y.shape[0]) #self.dL_dpsi0[:] = 0
self.dL_dpsi0 = -beta/2.* np.ones(self.Y.shape[0])
def parameters_changed(self):
psi0 = self.kern.psi0(self.Z, self.X)
psi1 = self.kern.psi1(self.Z, self.X)
psi2 = self.kern.psi2(self.Z, self.X)
if self.uncertain_input:
psi0 = self.kern.psi0(self.Z, self.X)
psi1 = self.kern.psi1(self.Z, self.X)
psi2 = self.kern.psi2(self.Z, self.X)
else:
psi0 = self.kern.Kdiag(self.X)
psi1 = self.kern.K(self.X, self.Z)
psi2 = np.dot(psi1.T,psi1)
self._log_marginal_likelihood = (self.dL_dpsi2*psi2).sum()+(self.dL_dpsi1*psi1).sum()+(self.dL_dpsi0*psi0).sum()
X_grad = self.kern.gradients_qX_expectations(variational_posterior=self.X, Z=self.Z, dL_dpsi0=self.dL_dpsi0, dL_dpsi1=self.dL_dpsi1, dL_dpsi2=self.dL_dpsi2)
self.X.set_gradients(X_grad)
from ...core.parameterization.variational import SpikeAndSlabPrior
if isinstance(self.variational_prior, SpikeAndSlabPrior):
# Update Log-likelihood
KL_div = self.variational_prior.KL_divergence(self.X, N=self.Y.shape[0])
# update for the KL divergence
self.variational_prior.update_gradients_KL(self.X, N=self.Y.shape[0])
if self.uncertain_input:
X_grad = self.kern.gradients_qX_expectations(variational_posterior=self.X, Z=self.Z, dL_dpsi0=self.dL_dpsi0, dL_dpsi1=self.dL_dpsi1, dL_dpsi2=self.dL_dpsi2)
self.X.set_gradients(X_grad)
else:
# Update Log-likelihood
KL_div = self.variational_prior.KL_divergence(self.X)
# update for the KL divergence
self.variational_prior.update_gradients_KL(self.X)
self._log_marginal_likelihood += -KL_div
dL_dpsi1 = self.dL_dpsi1 + 2.*np.dot(psi1,self.dL_dpsi2)
X_grad = self.kern.gradients_X_diag(self.dL_dpsi0, self.X)
X_grad += self.kern.gradients_X(dL_dpsi1, self.X, self.Z)
self.X.gradient = X_grad
if self.uncertain_input:
from ...core.parameterization.variational import SpikeAndSlabPrior
if isinstance(self.variational_prior, SpikeAndSlabPrior):
# Update Log-likelihood
KL_div = self.variational_prior.KL_divergence(self.X, N=self.Y.shape[0])
# update for the KL divergence
self.variational_prior.update_gradients_KL(self.X, N=self.Y.shape[0])
else:
# Update Log-likelihood
KL_div = self.variational_prior.KL_divergence(self.X)
# update for the KL divergence
self.variational_prior.update_gradients_KL(self.X)
self._log_marginal_likelihood += -KL_div
def log_likelihood(self):
return self._log_marginal_likelihood

16
GPy/models/bayesian_gplvm.py
View file

@ -141,22 +141,6 @@ class BayesianGPLVM(SparseGP_MPI):
resolution, ax, marker, s,
fignum, plot_inducing, legend,
plot_limits, aspect, updates, predict_kwargs, imshow_kwargs)
def infer_newX(self, Y_new, optimize=True, ):
"""
Infer the distribution of X for the new observed data *Y_new*.
:param model: the GPy model used in inference
:type model: GPy.core.Model
:param Y_new: the new observed data for inference
:type Y_new: numpy.ndarray
:param optimize: whether to optimize the location of new X (True by default)
:type optimize: boolean
:return: a tuple containing the estimated posterior distribution of X and the model that optimize X
:rtype: (GPy.core.parameterization.variational.VariationalPosterior, GPy.core.Model)
"""
from ..inference.latent_function_inference.inferenceX import infer_newX
return infer_newX(self, Y_new, optimize=optimize)
def do_test_latents(self, Y):
"""

3
GPy/models/sparse_gplvm.py
View file

@ -26,7 +26,8 @@ class SparseGPLVM(SparseGPRegression):
def parameters_changed(self):
super(SparseGPLVM, self).parameters_changed()
self.X.gradient = self.kern.gradients_X(self.grad_dict['dL_dKnm'], self.X, self.Z)
self.X.gradient = self.kern.gradients_X_diag(self.grad_dict['dL_dKdiag'], self.X)
self.X.gradient += self.kern.gradients_X(self.grad_dict['dL_dKnm'], self.X, self.Z)
def plot_latent(self, labels=None, which_indices=None,
resolution=50, ax=None, marker='o', s=40,

11
GPy/testing/inference_tests.py
View file

@ -65,6 +65,17 @@ class InferenceXTestCase(unittest.TestCase):
self.assertTrue(np.allclose(m.X.mean, mi.X.mean))
self.assertTrue(np.allclose(m.X.variance, mi.X.variance))
def test_inferenceX_GPLVM(self):
Ys = self.genData()
m = GPy.models.GPLVM(Ys[0],3,kernel=GPy.kern.RBF(3,ARD=True))
x,mi = m.infer_newX(m.Y, optimize=False)
self.assertTrue(mi.checkgrad())
# m.optimize(max_iters=10000)
# x,mi = m.infer_newX(m.Y)
# self.assertTrue(np.allclose(m.X, x))
if __name__ == "__main__":

63
GPy/testing/prior_tests.py
View file

@ -45,6 +45,69 @@ class PriorTests(unittest.TestCase):
# should raise an assertionerror.
self.assertRaises(AssertionError, m.rbf.set_prior, gaussian)
def test_set_prior(self):
xmin, xmax = 1, 2.5*np.pi
b, C, SNR = 1, 0, 0.1
X = np.linspace(xmin, xmax, 500)
y = b*X + C + 1*np.sin(X)
y += 0.05*np.random.randn(len(X))
X, y = X[:, None], y[:, None]
m = GPy.models.GPRegression(X, y)
gaussian = GPy.priors.Gaussian(1, 1)
#m.rbf.set_prior(gaussian)
# setting a Gaussian prior on non-negative parameters
# should raise an assertionerror.
self.assertRaises(AssertionError, m.rbf.set_prior, gaussian)
def test_set_gaussian_for_reals(self):
xmin, xmax = 1, 2.5*np.pi
b, C, SNR = 1, 0, 0.1
X = np.linspace(xmin, xmax, 500)
y = b*X + C + 1*np.sin(X)
y += 0.05*np.random.randn(len(X))
X, y = X[:, None], y[:, None]
m = GPy.models.SparseGPRegression(X, y)
gaussian = GPy.priors.Gaussian(1, 1)
m.Z.set_prior(gaussian)
# setting a Gaussian prior on non-negative parameters
# should raise an assertionerror.
#self.assertRaises(AssertionError, m.Z.set_prior, gaussian)
def test_fixed_domain_check(self):
xmin, xmax = 1, 2.5*np.pi
b, C, SNR = 1, 0, 0.1
X = np.linspace(xmin, xmax, 500)
y = b*X + C + 1*np.sin(X)
y += 0.05*np.random.randn(len(X))
X, y = X[:, None], y[:, None]
m = GPy.models.GPRegression(X, y)
m.rbf.fix()
gaussian = GPy.priors.Gaussian(1, 1)
# setting a Gaussian prior on non-negative parameters
# should raise an assertionerror.
self.assertRaises(AssertionError, m.rbf.set_prior, gaussian)
def test_fixed_domain_check1(self):
xmin, xmax = 1, 2.5*np.pi
b, C, SNR = 1, 0, 0.1
X = np.linspace(xmin, xmax, 500)
y = b*X + C + 1*np.sin(X)
y += 0.05*np.random.randn(len(X))
X, y = X[:, None], y[:, None]
m = GPy.models.GPRegression(X, y)
m.kern.lengthscale.fix()
gaussian = GPy.priors.Gaussian(1, 1)
# setting a Gaussian prior on non-negative parameters
# should raise an assertionerror.
self.assertRaises(AssertionError, m.rbf.set_prior, gaussian)
if __name__ == "__main__":
print "Running unit tests, please be (very) patient..."

6
README.md
View file

@ -3,11 +3,11 @@ GPy
A Gaussian processes framework in Python.
* [GPy homepage](http://sheffieldml.github.io/GPy/)
* [User mailing list](https://lists.shef.ac.uk/sympa/subscribe/gpy-users)
* [Online documentation](https://gpy.readthedocs.org/en/latest/)
* [Unit tests (Travis-CI)](https://travis-ci.org/SheffieldML/GPy)
Continuous integration status: ![CI status](https://travis-ci.org/SheffieldML/GPy.png)
Citation
@ -20,6 +20,10 @@ Citation
year = {2012--2014}
}
Pronounciation
==============
We like to pronounce it 'Gee-pie'.
Getting started
===============
Installing with pip