From c2568d2e9f81bbed537ef6a2e24076a97a55256c Mon Sep 17 00:00:00 2001 From: Zhenwen Dai Date: Thu, 29 May 2014 13:21:14 +0100 Subject: [PATCH 1/5] [splitkern] bug fix --- GPy/kern/_src/splitKern.py | 28 +++++++++++++++------------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/GPy/kern/_src/splitKern.py b/GPy/kern/_src/splitKern.py index 0aed7127..dc7bef8b 100644 --- a/GPy/kern/_src/splitKern.py +++ b/GPy/kern/_src/splitKern.py @@ -6,6 +6,8 @@ import numpy as np from kern import Kern,CombinationKernel from .independent_outputs import index_to_slices import itertools +from GPy.kern import Linear,RBF + class SplitKern(CombinationKernel): """ @@ -45,9 +47,9 @@ class SplitKern(CombinationKernel): # diagonal blocks [[target.__setitem__((s,s2), self.kern.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[i], slices2[i])] for i in xrange(min(len(slices),len(slices)))] if len(slices)>1: - [target.__setitem__((s,s2), self.kern.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[1], slices2[0])] + [target.__setitem__((s,s2), self.kern_cross.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[1], slices2[0])] if len(slices2)>1: - [target.__setitem__((s,s2), self.kern.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[0], slices2[1])] + [target.__setitem__((s,s2), self.kern_cross.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[0], slices2[1])] return target def Kdiag(self,X): @@ -60,7 +62,7 @@ class SplitKern(CombinationKernel): def collate_grads(dL, X, X2, cross=False): if cross: self.kern_cross.update_gradients_full(dL,X,X2) - target[:] += self.kern_cross.gradient + target[:] += self.kern_cross.kern.gradient else: self.kern.update_gradients_full(dL,X,X2) target[:] += self.kern.gradient @@ -102,22 +104,22 @@ class SplitKern_cross(Kern): def update_gradients_full(self, dL_dK, X, X2=None): if X2 is None: X2 = X - + k1 = self.kern.K(X,self.Xp) k2 = self.kern.K(self.Xp,X2) k3 = self.kern.K(self.Xp,self.Xp) - dL_dk1 = np.einsum('ij,j->i',dL_dK,k2.flat)/k3.flat - dL_dk2 = np.einsum('ij,i->j',dL_dK,k1.flat)/k3.flat - dL_dk3 = np.einsum('ij,ij->',dL_dK,-np.dot(k1,k2)/(k3.flat*k3.flat)) - + dL_dk1 = np.einsum('ij,j->i',dL_dK,k2[0])/k3[0,0] + dL_dk2 = np.einsum('ij,i->j',dL_dK,k1[:,0])/k3[0,0] + dL_dk3 = np.einsum('ij,ij->',dL_dK,-np.dot(k1,k2)/(k3[0,0]*k3[0,0])) + self.kern.update_gradients_full(dL_dk1[:,None],X,self.Xp) - grad1 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dk2[None,:],self.Xp,X) - grad2 = self.kern.gradient.copy() + grad = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dk2[None,:],self.Xp,X2) + grad += self.kern.gradient.copy() self.kern.update_gradients_full(np.array([[dL_dk3]]),self.Xp,self.Xp) - grad3 = self.kern.gradient.copy() + grad += self.kern.gradient.copy() - self.kern.gradient = grad1+grad2+grad3 + self.kern.gradient = grad def update_gradients_diag(self, dL_dKdiag, X): k1 = self.kern.K(X,self.Xp) From 5f081526513889721614b2e4c4704a58d7f7d3b0 Mon Sep 17 00:00:00 2001 From: Zhenwen Dai Date: Thu, 29 May 2014 14:03:59 +0100 Subject: [PATCH 2/5] [splitkern] some additional implmentation --- GPy/kern/__init__.py | 2 +- GPy/kern/_src/splitKern.py | 69 ++++++++++++++++++++++++++++++++------ 2 files changed, 59 insertions(+), 12 deletions(-) diff --git a/GPy/kern/__init__.py b/GPy/kern/__init__.py index 85c6cf4d..2faf57c6 100644 --- a/GPy/kern/__init__.py +++ b/GPy/kern/__init__.py @@ -14,7 +14,7 @@ from _src.ODE_UYC import ODE_UYC from _src.ODE_st import ODE_st from _src.ODE_t import ODE_t from _src.poly import Poly -from _src.splitKern import SplitKern +from _src.splitKern import SplitKern,DiffGenomeKern # TODO: put this in an init file somewhere #I'm commenting this out because the files were not added. JH. Remember to add the files before commiting diff --git a/GPy/kern/_src/splitKern.py b/GPy/kern/_src/splitKern.py index dc7bef8b..6729b56f 100644 --- a/GPy/kern/_src/splitKern.py +++ b/GPy/kern/_src/splitKern.py @@ -6,21 +6,66 @@ import numpy as np from kern import Kern,CombinationKernel from .independent_outputs import index_to_slices import itertools -from GPy.kern import Linear,RBF +class DiffGenomeKern(Kern): + + def __init__(self, kernel, idx_p, Xp, index_dim=-1, name='DiffGenomeKern'): + self.idx_p = idx_p + self.index_dim=index_dim + self.kern = SplitKern(kernel,Xp, index_dim=index_dim) + super(DiffGenomeKern, self).__init__(input_dim=kernel.input_dim+1, name=name) + self.add_parameter(self.kern) + + def K(self, X, X2=None): + assert X2==None + K = self.kern.K(X,X2) + + slices = index_to_slices(X[:,self.index_dim]) + idx_start = slices[1][0] + idx_end = idx_start+self.idx_p + K[idx_start:idx_end,:] = K[:self.idx_p,:] + K[:,idx_start:idx_end] = K[:,self.idx_p] + + return K + + def Kdiag(self,X): + Kdiag = self.kern.Kdiag(X) + + slices = index_to_slices(X[:,self.index_dim]) + idx_start = slices[1][0] + idx_end = idx_start+self.idx_p + Kdiag[idx_start:idx_end] = Kdiag[:self.idx_p] + + return Kdiag + + def update_gradients_full(self,dL_dK,X,X2=None): + assert X2==None + slices = index_to_slices(X[:,self.index_dim]) + idx_start = slices[1][0] + idx_end = idx_start+self.idx_p + + self.kern.update_gradients_full(dL_dK, X[:self.idx_p],X) + grad_p1 = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dK, X, X[:self.idx_p]) + grad_p2 = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dK, X[:self.idx_p], X[:self.idx_p]) + grad_p3 = self.kern.gradient.copy() + + self.kern.update_gradients_full(dL_dK, X[idx_start:idx_end],X) + grad_n1 = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dK, X, X[idx_start:idx_end]) + grad_n2 = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dK, X[idx_start:idx_end], X[idx_start:idx_end]) + grad_n3 = self.kern.gradient.copy() + + self.kern.update_gradients_full(dL_dK, X) + self.kern.gradient += grad_p1+grad_p2+grad_p3-grad_n1-grad_n2-grad_n3 + + def update_gradients_diag(self, dL_dKdiag, X): + pass class SplitKern(CombinationKernel): - """ - A kernel which can represent several independent functions. this kernel - 'switches off' parts of the matrix where the output indexes are different. - The index of the functions is given by the last column in the input X the - rest of the columns of X are passed to the underlying kernel for - computation (in blocks). - - :param kernels: either a kernel, or list of kernels to work with. If it is - a list of kernels the indices in the index_dim, index the kernels you gave! - """ def __init__(self, kernel, Xp, index_dim=-1, name='SplitKern'): assert isinstance(index_dim, int), "The index dimension must be an integer!" self.kern = kernel @@ -89,6 +134,8 @@ class SplitKern_cross(Kern): def __init__(self, kernel, Xp, name='SplitKern_cross'): assert isinstance(kernel, Kern) self.kern = kernel + if not isinstance(Xp,np.ndarray): + Xp = np.array([[Xp]]) self.Xp = Xp super(SplitKern_cross, self).__init__(input_dim=kernel.input_dim, active_dims=None, name=name) From ea1326c8c8db385ac9249371b696e87253a587f6 Mon Sep 17 00:00:00 2001 From: zhenwen Date: Thu, 29 May 2014 14:54:11 +0100 Subject: [PATCH 3/5] [splitkern] some more changes --- GPy/kern/_src/splitKern.py | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/GPy/kern/_src/splitKern.py b/GPy/kern/_src/splitKern.py index 6729b56f..54a08914 100644 --- a/GPy/kern/_src/splitKern.py +++ b/GPy/kern/_src/splitKern.py @@ -13,7 +13,7 @@ class DiffGenomeKern(Kern): self.idx_p = idx_p self.index_dim=index_dim self.kern = SplitKern(kernel,Xp, index_dim=index_dim) - super(DiffGenomeKern, self).__init__(input_dim=kernel.input_dim+1, name=name) + super(DiffGenomeKern, self).__init__(input_dim=kernel.input_dim+1, active_dims=None, name=name) self.add_parameter(self.kern) def K(self, X, X2=None): @@ -21,10 +21,12 @@ class DiffGenomeKern(Kern): K = self.kern.K(X,X2) slices = index_to_slices(X[:,self.index_dim]) - idx_start = slices[1][0] + idx_start = slices[1][0].start idx_end = idx_start+self.idx_p + K_c = K[idx_start:idx_end,idx_start:idx_end].copy() K[idx_start:idx_end,:] = K[:self.idx_p,:] - K[:,idx_start:idx_end] = K[:,self.idx_p] + K[:,idx_start:idx_end] = K[:,:self.idx_p] + K[idx_start:idx_end,idx_start:idx_end] = K_c return K @@ -32,7 +34,7 @@ class DiffGenomeKern(Kern): Kdiag = self.kern.Kdiag(X) slices = index_to_slices(X[:,self.index_dim]) - idx_start = slices[1][0] + idx_start = slices[1][0].start idx_end = idx_start+self.idx_p Kdiag[idx_start:idx_end] = Kdiag[:self.idx_p] @@ -41,7 +43,7 @@ class DiffGenomeKern(Kern): def update_gradients_full(self,dL_dK,X,X2=None): assert X2==None slices = index_to_slices(X[:,self.index_dim]) - idx_start = slices[1][0] + idx_start = slices[1][0].start idx_end = idx_start+self.idx_p self.kern.update_gradients_full(dL_dK, X[:self.idx_p],X) @@ -59,7 +61,7 @@ class DiffGenomeKern(Kern): grad_n3 = self.kern.gradient.copy() self.kern.update_gradients_full(dL_dK, X) - self.kern.gradient += grad_p1+grad_p2+grad_p3-grad_n1-grad_n2-grad_n3 + self.kern.gradient += grad_p1+grad_p2-2*grad_p3-grad_n1-grad_n2+2*grad_n3 def update_gradients_diag(self, dL_dKdiag, X): pass From 4c538efb64192470449b6c7247d0b445f4404e5a Mon Sep 17 00:00:00 2001 From: Zhenwen Dai Date: Thu, 29 May 2014 17:10:14 +0100 Subject: [PATCH 4/5] DiffGenomeKern bug fix --- GPy/kern/_src/splitKern.py | 28 +++++++++++++++++----------- 1 file changed, 17 insertions(+), 11 deletions(-) diff --git a/GPy/kern/_src/splitKern.py b/GPy/kern/_src/splitKern.py index 54a08914..624bae58 100644 --- a/GPy/kern/_src/splitKern.py +++ b/GPy/kern/_src/splitKern.py @@ -6,6 +6,7 @@ import numpy as np from kern import Kern,CombinationKernel from .independent_outputs import index_to_slices import itertools +from .rbf import RBF class DiffGenomeKern(Kern): @@ -13,6 +14,7 @@ class DiffGenomeKern(Kern): self.idx_p = idx_p self.index_dim=index_dim self.kern = SplitKern(kernel,Xp, index_dim=index_dim) +# self.kern = RBF(1) super(DiffGenomeKern, self).__init__(input_dim=kernel.input_dim+1, active_dims=None, name=name) self.add_parameter(self.kern) @@ -46,22 +48,24 @@ class DiffGenomeKern(Kern): idx_start = slices[1][0].start idx_end = idx_start+self.idx_p - self.kern.update_gradients_full(dL_dK, X[:self.idx_p],X) + self.kern.update_gradients_full(dL_dK[idx_start:idx_end,:], X[:self.idx_p],X) grad_p1 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dK, X, X[:self.idx_p]) + self.kern.update_gradients_full(dL_dK[:,idx_start:idx_end], X, X[:self.idx_p]) grad_p2 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dK, X[:self.idx_p], X[:self.idx_p]) + self.kern.update_gradients_full(dL_dK[idx_start:idx_end,idx_start:idx_end], X[:self.idx_p],X[idx_start:idx_end]) grad_p3 = self.kern.gradient.copy() + self.kern.update_gradients_full(dL_dK[idx_start:idx_end,idx_start:idx_end], X[idx_start:idx_end], X[:self.idx_p]) + grad_p4 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dK, X[idx_start:idx_end],X) + self.kern.update_gradients_full(dL_dK[idx_start:idx_end,:], X[idx_start:idx_end],X) grad_n1 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dK, X, X[idx_start:idx_end]) + self.kern.update_gradients_full(dL_dK[:,idx_start:idx_end], X, X[idx_start:idx_end]) grad_n2 = self.kern.gradient.copy() - self.kern.update_gradients_full(dL_dK, X[idx_start:idx_end], X[idx_start:idx_end]) + self.kern.update_gradients_full(dL_dK[idx_start:idx_end,idx_start:idx_end], X[idx_start:idx_end], X[idx_start:idx_end]) grad_n3 = self.kern.gradient.copy() self.kern.update_gradients_full(dL_dK, X) - self.kern.gradient += grad_p1+grad_p2-2*grad_p3-grad_n1-grad_n2+2*grad_n3 + self.kern.gradient += grad_p1+grad_p2-grad_p3-grad_p4-grad_n1-grad_n2+2*grad_n3 def update_gradients_diag(self, dL_dKdiag, X): pass @@ -92,7 +96,7 @@ class SplitKern(CombinationKernel): assert len(slices2)<=2, 'The Split kernel only support two different indices' target = np.zeros((X.shape[0], X2.shape[0])) # diagonal blocks - [[target.__setitem__((s,s2), self.kern.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[i], slices2[i])] for i in xrange(min(len(slices),len(slices)))] + [[target.__setitem__((s,s2), self.kern.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[i], slices2[i])] for i in xrange(min(len(slices),len(slices2)))] if len(slices)>1: [target.__setitem__((s,s2), self.kern_cross.K(X[s,:],X2[s2,:])) for s,s2 in itertools.product(slices[1], slices2[0])] if len(slices2)>1: @@ -115,17 +119,19 @@ class SplitKern(CombinationKernel): target[:] += self.kern.gradient if X2 is None: + assert dL_dK.shape==(X.shape[0],X.shape[0]) [[collate_grads(dL_dK[s,ss], X[s], X[ss]) for s,ss in itertools.product(slices_i, slices_i)] for slices_i in slices] if len(slices)>1: [collate_grads(dL_dK[s,ss], X[s], X[ss], True) for s,ss in itertools.product(slices[0], slices[1])] [collate_grads(dL_dK[s,ss], X[s], X[ss], True) for s,ss in itertools.product(slices[1], slices[0])] else: + assert dL_dK.shape==(X.shape[0],X2.shape[0]) slices2 = index_to_slices(X2[:,self.index_dim]) - [[collate_grads(dL_dK[s,s2],X[s],X2[s2]) for s,s2 in itertools.product(slices[i], slices2[i])] for i in xrange(min(len(slices),len(slices)))] + [[collate_grads(dL_dK[s,s2],X[s],X2[s2]) for s,s2 in itertools.product(slices[i], slices2[i])] for i in xrange(min(len(slices),len(slices2)))] if len(slices)>1: - [collate_grads(dL_dK[s,ss], X[s], X2[s2], True) for s,s2 in itertools.product(slices[1], slices2[0])] + [collate_grads(dL_dK[s,s2], X[s], X2[s2], True) for s,s2 in itertools.product(slices[1], slices2[0])] if len(slices2)>1: - [collate_grads(dL_dK[s,ss], X[s], X2[s2], True) for s,s2 in itertools.product(slices[0], slices2[1])] + [collate_grads(dL_dK[s,s2], X[s], X2[s2], True) for s,s2 in itertools.product(slices[0], slices2[1])] self.kern.gradient = target def update_gradients_diag(self, dL_dKdiag, X): From 47ba2542c2e9c65d5cc797051dd1ed755422d0c0 Mon Sep 17 00:00:00 2001 From: Zhenwen Dai Date: Thu, 29 May 2014 17:11:26 +0100 Subject: [PATCH 5/5] minor changes --- GPy/kern/_src/splitKern.py | 2 -- 1 file changed, 2 deletions(-) diff --git a/GPy/kern/_src/splitKern.py b/GPy/kern/_src/splitKern.py index 624bae58..dfaf5108 100644 --- a/GPy/kern/_src/splitKern.py +++ b/GPy/kern/_src/splitKern.py @@ -6,7 +6,6 @@ import numpy as np from kern import Kern,CombinationKernel from .independent_outputs import index_to_slices import itertools -from .rbf import RBF class DiffGenomeKern(Kern): @@ -14,7 +13,6 @@ class DiffGenomeKern(Kern): self.idx_p = idx_p self.index_dim=index_dim self.kern = SplitKern(kernel,Xp, index_dim=index_dim) -# self.kern = RBF(1) super(DiffGenomeKern, self).__init__(input_dim=kernel.input_dim+1, active_dims=None, name=name) self.add_parameter(self.kern)