From f50fac810418a8d56e2ecf1729f6eca14f346f0e Mon Sep 17 00:00:00 2001 From: Ricardo Date: Mon, 3 Feb 2014 15:29:06 +0000 Subject: [PATCH 1/2] minor changes --- GPy/kern/constructors.py | 8 ++++---- GPy/kern/kern.py | 4 ++-- GPy/kern/parts/prod.py | 25 +++++++++++++------------ 3 files changed, 19 insertions(+), 18 deletions(-) diff --git a/GPy/kern/constructors.py b/GPy/kern/constructors.py index 11da5e9b..a362aff8 100644 --- a/GPy/kern/constructors.py +++ b/GPy/kern/constructors.py @@ -411,7 +411,7 @@ def periodic_Matern52(input_dim, variance=1., lengthscale=None, period=2 * np.pi part = parts.periodic_Matern52.PeriodicMatern52(input_dim, variance, lengthscale, period, n_freq, lower, upper) return kern(input_dim, [part]) -def prod(k1,k2,tensor=False): +def prod(k1,k2,**kwargs): """ Construct a product kernel over input_dim from two kernels over input_dim @@ -422,7 +422,7 @@ def prod(k1,k2,tensor=False): :rtype: kernel object """ - part = parts.prod.Prod(k1, k2, tensor) + part = parts.prod.Prod(k1, k2, **kwargs) return kern(part.input_dim, [part]) def symmetric(k): @@ -433,7 +433,7 @@ def symmetric(k): k_.parts = [symmetric.Symmetric(p) for p in k.parts] return k_ -def coregionalize(output_dim,rank=1, W=None, kappa=None): +def coregionalize(output_dim,rank=1, W=None, kappa=None,name='coregion'): """ Coregionlization matrix B, of the form: @@ -458,7 +458,7 @@ def coregionalize(output_dim,rank=1, W=None, kappa=None): :rtype: kernel object """ - p = parts.coregionalize.Coregionalize(output_dim,rank,W,kappa) + p = parts.coregionalize.Coregionalize(output_dim,rank,W,kappa,name) return kern(1,[p]) diff --git a/GPy/kern/kern.py b/GPy/kern/kern.py index de87ff14..6c24ab7d 100644 --- a/GPy/kern/kern.py +++ b/GPy/kern/kern.py @@ -170,7 +170,7 @@ class kern(Parameterized): """ return self.prod(other, tensor=True) - def prod(self, other, tensor=False): + def prod(self, other, tensor=False, **kwargs): """ Multiply two kernels (either on the same space, or on the tensor product of the input space). @@ -189,7 +189,7 @@ class kern(Parameterized): s1[sl1], s2[sl2] = [True], [True] slices += [s1 + s2] - newkernparts = [prod(k1, k2, tensor) for k1, k2 in itertools.product(K1._parameters_, K2._parameters_)] + newkernparts = [prod(k1, k2, tensor, **kwargs) for k1, k2 in itertools.product(K1._parameters_, K2._parameters_)] if tensor: newkern = kern(K1.input_dim + K2.input_dim, newkernparts, slices) diff --git a/GPy/kern/parts/prod.py b/GPy/kern/parts/prod.py index 2569c51c..3e624756 100644 --- a/GPy/kern/parts/prod.py +++ b/GPy/kern/parts/prod.py @@ -17,21 +17,21 @@ class Prod(Kernpart): :rtype: kernel object """ - def __init__(self,k1,k2,tensor=False): + def __init__(self,k1,k2,tensor=False,name="product"): if tensor: - super(Prod, self).__init__(k1.input_dim + k2.input_dim, '['+k1.name + '**' + k2.name +']') + super(Prod, self).__init__(k1.input_dim + k2.input_dim, name) else: assert k1.input_dim == k2.input_dim, "Error: The input spaces of the kernels to sum don't have the same dimension." - super(Prod, self).__init__(k1.input_dim, '['+k1.name + '*' + k2.name +']') + super(Prod, self).__init__(k1.input_dim, name) #self.num_params = k1.num_params + k2.num_params self.k1 = k1 self.k2 = k2 -# if tensor: -# self.slice1 = slice(0,self.k1.input_dim) -# self.slice2 = slice(self.k1.input_dim,self.k1.input_dim+self.k2.input_dim) -# else: -# self.slice1 = slice(0,self.input_dim) -# self.slice2 = slice(0,self.input_dim) + if tensor: + self.slice1 = slice(0,self.k1.input_dim) + self.slice2 = slice(self.k1.input_dim,self.k1.input_dim+self.k2.input_dim) + else: + self.slice1 = slice(0,self.input_dim) + self.slice2 = slice(0,self.input_dim) self._X, self._X2 = np.empty(shape=(2,1)) self.add_parameters(self.k1, self.k2) @@ -117,18 +117,19 @@ class Prod(Kernpart): def _K_computations(self,X,X2): if not (np.array_equal(X,self._X) and np.array_equal(X2,self._X2) and np.array_equal(self._params , self._get_params())): - self._X = X.copy() - self._params == self._get_params().copy() + #self._X = X.copy() + #self._params == self._get_params().copy() if X2 is None: self._X2 = None self._K1 = np.zeros((X.shape[0],X.shape[0])) self._K2 = np.zeros((X.shape[0],X.shape[0])) self.k1.K(X[:,self.slice1],None,self._K1) self.k2.K(X[:,self.slice2],None,self._K2) + #self.k1.K(X[:,self.k1.input_slices],None,self._K1) + #self.k2.K(X[:,self.k2_input_slices],None,self._K2) else: self._X2 = X2.copy() self._K1 = np.zeros((X.shape[0],X2.shape[0])) self._K2 = np.zeros((X.shape[0],X2.shape[0])) self.k1.K(X[:,self.slice1],X2[:,self.slice1],self._K1) self.k2.K(X[:,self.slice2],X2[:,self.slice2],self._K2) - From 0cf1e44c685e27051edb0bdc9cefe7e4ca7519fe Mon Sep 17 00:00:00 2001 From: Ricardo Date: Wed, 26 Feb 2014 11:38:46 +0000 Subject: [PATCH 2/2] some work on ep, and some messing with wher ethe derivatives are computed (in the model, not the inference --- GPy/core/gp.py | 2 +- GPy/inference/latent_function_inference/ep.py | 15 +++++++++++ .../exact_gaussian_inference.py | 3 +-- GPy/util/multioutput.py | 25 +++++++++++++++++-- 4 files changed, 40 insertions(+), 5 deletions(-) diff --git a/GPy/core/gp.py b/GPy/core/gp.py index 81985773..50f9d342 100644 --- a/GPy/core/gp.py +++ b/GPy/core/gp.py @@ -70,7 +70,7 @@ class GP(Model): def parameters_changed(self): self.posterior, self._log_marginal_likelihood, grad_dict = self.inference_method.inference(self.kern, self.X, self.likelihood, self.Y, Y_metadata=self.Y_metadata) - self._dL_dK = grad_dict['dL_dK'] + self.kern.update_gradients_full(grad_dict['dL_dK']) def log_likelihood(self): return self._log_marginal_likelihood diff --git a/GPy/inference/latent_function_inference/ep.py b/GPy/inference/latent_function_inference/ep.py index 87c08221..1904d48c 100644 --- a/GPy/inference/latent_function_inference/ep.py +++ b/GPy/inference/latent_function_inference/ep.py @@ -23,11 +23,26 @@ class EP(object): self.old_mutilde, self.old_vtilde = None, None def inference(self, kern, X, likelihood, Y, Y_metadata=None): + num_data, output_dim = X.shape + assert output_dim ==1, "ep in 1D only (for now!)" K = kern.K(X) mu_tilde, tau_tilde = self.expectation_propagation() + Wi, LW, LWi, W_logdet = pdinv(K + np.diag(1./tau_tilde) + + alpha, _ = dpotrs(LW, mu_tilde, lower=1) + + log_marginal = 0.5*(-num_data * log_2_pi - W_logdet - np.sum(alpha * mu_tilde)) + + dL_dK = 0.5 * (tdot(alpha[:,None]) - Wi) + + #TODO: what abot derivatives of the likelihood parameters? + + return Posterior(woodbury_inv=Wi, woodbury_vector=alpha, K=K), log_marginal, {'dL_dK':dL_dK} + + def expectation_propagation(self, K, Y, Y_metadata, likelihood) diff --git a/GPy/inference/latent_function_inference/exact_gaussian_inference.py b/GPy/inference/latent_function_inference/exact_gaussian_inference.py index 2c3bfa6a..922b52f4 100644 --- a/GPy/inference/latent_function_inference/exact_gaussian_inference.py +++ b/GPy/inference/latent_function_inference/exact_gaussian_inference.py @@ -49,8 +49,7 @@ class ExactGaussianInference(object): dL_dK = 0.5 * (tdot(alpha) - Y.shape[1] * Wi) - kern.update_gradients_full(dL_dK, X) - + #TODO: does this really live here? likelihood.update_gradients(np.diag(dL_dK)) return Posterior(woodbury_chol=LW, woodbury_vector=alpha, K=K), log_marginal, {'dL_dK':dL_dK} diff --git a/GPy/util/multioutput.py b/GPy/util/multioutput.py index a57593a7..eb4d8d08 100644 --- a/GPy/util/multioutput.py +++ b/GPy/util/multioutput.py @@ -2,6 +2,27 @@ import numpy as np import warnings from .. import kern +def build_XY(input_list,output_list=None,index=None): + num_outputs = len(input_list) + _s = [0] + [ _x.shape[0] for _x in input_list ] + _s = np.cumsum(_s) + slices = [slice(a,b) for a,b in zip(_s[:-1],_s[1:])] + if output_list is not None: + assert num_outputs == len(output_list) + Y = np.vstack(output_list) + else: + Y = None + + if index is not None: + assert len(index) == num_outputs + I = np.vstack( [j*np.ones((_x.shape[0],1)) for _x,j in zip(input_list,index)] ) + else: + I = np.vstack( [j*np.ones((_x.shape[0],1)) for _x,j in zip(input_list,range(num_outputs))] ) + + X = np.vstack(input_list) + X = np.hstack([X,I]) + return X,Y,slices + def build_lcm(input_dim, num_outputs, CK = [], NC = [], W_columns=1,W=None,kappa=None): #TODO build_icm or build_lcm """ @@ -25,9 +46,9 @@ def build_lcm(input_dim, num_outputs, CK = [], NC = [], W_columns=1,W=None,kappa k.input_dim = input_dim + 1 warnings.warn("kernel's input dimension overwritten to fit input_dim parameter.") - kernel = CK[0].prod(kern.coregionalize(num_outputs,W_columns,W,kappa),tensor=True) + kernel = CK[0].prod(kern.Coregionalize(num_outputs,W_columns,W,kappa),tensor=True) for k in CK[1:]: - k_coreg = kern.coregionalize(num_outputs,W_columns,W,kappa) + k_coreg = kern.Coregionalize(num_outputs,W_columns,W,kappa) kernel += k.prod(k_coreg,tensor=True) for k in NC: kernel += k