From 19407293dcaf7b642ab4335fea985aa7ee398ec7 Mon Sep 17 00:00:00 2001 From: Nicolo Fusi Date: Fri, 17 May 2013 16:29:32 +0100 Subject: [PATCH] cross-terms --- GPy/examples/dimensionality_reduction.py | 4 +- GPy/kern/bias.py | 5 +- GPy/kern/kern.py | 130 ++++++----------------- GPy/kern/kernpart.py | 2 - GPy/kern/white.py | 4 +- 5 files changed, 40 insertions(+), 105 deletions(-) diff --git a/GPy/examples/dimensionality_reduction.py b/GPy/examples/dimensionality_reduction.py index dcda4f42..9b51947b 100644 --- a/GPy/examples/dimensionality_reduction.py +++ b/GPy/examples/dimensionality_reduction.py @@ -17,11 +17,11 @@ def BGPLVM(seed=default_seed): D = 4 # generate GPLVM-like data X = np.random.rand(N, Q) - k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001) + k = GPy.kern.rbf(Q) + GPy.kern.white(Q, 0.00001) K = k.K(X) Y = np.random.multivariate_normal(np.zeros(N), K, D).T - k = GPy.kern.linear(Q, ARD=True) + GPy.kern.white(Q) + k = GPy.kern.rbf(Q, ARD=True) + GPy.kern.linear(Q, ARD=True) + GPy.kern.rbf(Q, ARD=True) + GPy.kern.white(Q) # k = GPy.kern.rbf(Q) + GPy.kern.rbf(Q) + GPy.kern.white(Q) # k = GPy.kern.rbf(Q) + GPy.kern.bias(Q) + GPy.kern.white(Q, 0.00001) # k = GPy.kern.rbf(Q, ARD = False) + GPy.kern.white(Q, 0.00001) diff --git a/GPy/kern/bias.py b/GPy/kern/bias.py index b5883f87..09f0afa9 100644 --- a/GPy/kern/bias.py +++ b/GPy/kern/bias.py @@ -55,8 +55,9 @@ class bias(kernpart): target += self.variance def psi1(self, Z, mu, S, target): - target += self.variance - + self._psi1 = self.variance + target += self._psi1 + def psi2(self, Z, mu, S, target): target += self.variance**2 diff --git a/GPy/kern/kern.py b/GPy/kern/kern.py index 0e425e38..c9582ac8 100644 --- a/GPy/kern/kern.py +++ b/GPy/kern/kern.py @@ -316,32 +316,19 @@ class kern(parameterised): # compute the "cross" terms # TODO: input_slices needed crossterms = 0 + for p1, p2 in itertools.combinations(self.parts, 2): - prod = np.multiply - # # white doesn;t combine with anything - # if p1.name == 'white' or p2.name == 'white': - # pass - # # rbf X bias - # elif p1.name == 'bias' and p2.name == 'rbf': - # target += p1.variance * (p2._psi1[:, :, None] + p2._psi1[:, None, :]) - # elif p2.name == 'bias' and p1.name == 'rbf': - # target += p2.variance * (p1._psi1[:, :, None] + p1._psi1[:, None, :]) - # # linear X bias - # elif p1.name == 'bias' and p2.name == 'linear': - # tmp = np.zeros((mu.shape[0], Z.shape[0])) - # p2.psi1(Z, mu, S, tmp) - # target += p1.variance * (tmp[:, :, None] + tmp[:, None, :]) - # elif p2.name == 'bias' and p1.name == 'linear': - # tmp = np.zeros((mu.shape[0], Z.shape[0])) - # p1.psi1(Z, mu, S, tmp) - # target += p2.variance * (tmp[:, :, None] + tmp[:, None, :]) - # # rbf X linear - # elif p1.name == 'linear' and p2.name == 'rbf': - # raise NotImplementedError # TODO - # elif p2.name == 'linear' and p1.name == 'rbf': - # raise NotImplementedError # TODO - # else: - # raise NotImplementedError, "psi2 cannot be computed for this kernel" + + # TODO psi1 this must be faster/better/precached/more nice + tmp1 = np.zeros((mu.shape[0], Z.shape[0])) + p1.psi1(Z, mu, S, tmp1) + tmp2 = np.zeros((mu.shape[0], Z.shape[0])) + p2.psi1(Z, mu, S, tmp2) + + prod = np.multiply(tmp1, tmp2) + crossterms += prod[:,:,None] + prod[:, None, :] + + target += crossterms return target def dpsi2_dtheta(self, dL_dpsi2, Z, mu, S): @@ -350,71 +337,34 @@ class kern(parameterised): # compute the "cross" terms # TODO: better looping, input_slices - for i1, i2 in itertools.combinations(range(len(self.parts)), 2): + for i1, i2 in itertools.permutations(range(len(self.parts)), 2): p1, p2 = self.parts[i1], self.parts[i2] # ipsl1, ipsl2 = self.input_slices[i1], self.input_slices[i2] ps1, ps2 = self.param_slices[i1], self.param_slices[i2] - # white doesn;t combine with anything - if p1.name == 'white' or p2.name == 'white': - pass - # rbf X bias - elif p1.name == 'bias' and p2.name == 'rbf': - p2.dpsi1_dtheta(dL_dpsi2.sum(1) * p1.variance * 2., Z, mu, S, target[ps2]) - p1.dpsi1_dtheta(dL_dpsi2.sum(1) * p2._psi1 * 2., Z, mu, S, target[ps1]) - elif p2.name == 'bias' and p1.name == 'rbf': - p1.dpsi1_dtheta(dL_dpsi2.sum(1) * p2.variance * 2., Z, mu, S, target[ps1]) - p2.dpsi1_dtheta(dL_dpsi2.sum(1) * p1._psi1 * 2., Z, mu, S, target[ps2]) - # linear X bias - elif p1.name == 'bias' and p2.name == 'linear': - p2.dpsi1_dtheta(dL_dpsi2.sum(1) * p1.variance * 2., Z, mu, S, target[ps2]) # [ps1]) - psi1 = np.zeros((mu.shape[0], Z.shape[0])) - p2.psi1(Z, mu, S, psi1) - p1.dpsi1_dtheta(dL_dpsi2.sum(1) * psi1 * 2., Z, mu, S, target[ps1]) - elif p2.name == 'bias' and p1.name == 'linear': - p1.dpsi1_dtheta(dL_dpsi2.sum(1) * p2.variance * 2., Z, mu, S, target[ps1]) - psi1 = np.zeros((mu.shape[0], Z.shape[0])) - p1.psi1(Z, mu, S, psi1) - p2.dpsi1_dtheta(dL_dpsi2.sum(1) * psi1 * 2., Z, mu, S, target[ps2]) - # rbf X linear - elif p1.name == 'linear' and p2.name == 'rbf': - raise NotImplementedError # TODO - elif p2.name == 'linear' and p1.name == 'rbf': - raise NotImplementedError # TODO - else: - raise NotImplementedError, "psi2 cannot be computed for this kernel" + tmp = np.zeros((mu.shape[0], Z.shape[0])) + p1.psi1(Z, mu, S, tmp) + p2.dpsi1_dtheta((tmp[:,None,:]*dL_dpsi2).sum(1)*2., Z, mu, S, target[ps2]) return self._transform_gradients(target) def dpsi2_dZ(self, dL_dpsi2, Z, mu, S): target = np.zeros_like(Z) [p.dpsi2_dZ(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + #target *= 2 # compute the "cross" terms # TODO: we need input_slices here. - for p1, p2 in itertools.combinations(self.parts, 2): - # white doesn;t combine with anything - if p1.name == 'white' or p2.name == 'white': - pass - # rbf X bias - elif p1.name == 'bias' and p2.name == 'rbf': - p2.dpsi1_dX(dL_dpsi2.sum(1).T * p1.variance, Z, mu, S, target) - elif p2.name == 'bias' and p1.name == 'rbf': - p1.dpsi1_dZ(dL_dpsi2.sum(1).T * p2.variance, Z, mu, S, target) - # linear X bias - elif p1.name == 'bias' and p2.name == 'linear': - p2.dpsi1_dZ(dL_dpsi2.sum(1).T * p1.variance, Z, mu, S, target) - elif p2.name == 'bias' and p1.name == 'linear': - p1.dpsi1_dZ(dL_dpsi2.sum(1).T * p2.variance, Z, mu, S, target) - # rbf X linear - elif p1.name == 'linear' and p2.name == 'rbf': - raise NotImplementedError # TODO - elif p2.name == 'linear' and p1.name == 'rbf': - raise NotImplementedError # TODO - else: - raise NotImplementedError, "psi2 cannot be computed for this kernel" + for p1, p2 in itertools.permutations(self.parts, 2): + if p1.name == 'linear' and p2.name == 'linear': + raise NotImplementedError("We don't handle linear/linear cross-terms") + tmp = np.zeros((mu.shape[0], Z.shape[0])) + p1.psi1(Z, mu, S, tmp) + tmp2 = np.zeros_like(target) + p2.dpsi1_dZ((tmp[:,None,:]*dL_dpsi2).sum(1).T, Z, mu, S, tmp2) + target += tmp2 - return target * 2. + return target * 2 def dpsi2_dmuS(self, dL_dpsi2, Z, mu, S): target_mu, target_S = np.zeros((2, mu.shape[0], mu.shape[1])) @@ -422,27 +372,13 @@ class kern(parameterised): # compute the "cross" terms # TODO: we need input_slices here. - for p1, p2 in itertools.combinations(self.parts, 2): - # white doesn;t combine with anything - if p1.name == 'white' or p2.name == 'white': - pass - # rbf X bias - elif p1.name == 'bias' and p2.name == 'rbf': - p2.dpsi1_dmuS(dL_dpsi2.sum(1).T * p1.variance * 2., Z, mu, S, target_mu, target_S) - elif p2.name == 'bias' and p1.name == 'rbf': - p1.dpsi1_dmuS(dL_dpsi2.sum(1).T * p2.variance * 2., Z, mu, S, target_mu, target_S) - # linear X bias - elif p1.name == 'bias' and p2.name == 'linear': - p2.dpsi1_dmuS(dL_dpsi2.sum(1).T * p1.variance * 2., Z, mu, S, target_mu, target_S) - elif p2.name == 'bias' and p1.name == 'linear': - p1.dpsi1_dmuS(dL_dpsi2.sum(1).T * p2.variance * 2., Z, mu, S, target_mu, target_S) - # rbf X linear - elif p1.name == 'linear' and p2.name == 'rbf': - raise NotImplementedError # TODO - elif p2.name == 'linear' and p1.name == 'rbf': - raise NotImplementedError # TODO - else: - raise NotImplementedError, "psi2 cannot be computed for this kernel" + for p1, p2 in itertools.permutations(self.parts, 2): + if p1.name == 'linear' and p2.name == 'linear': + raise NotImplementedError("We don't handle linear/linear cross-terms") + + tmp = np.zeros((mu.shape[0], Z.shape[0])) + p1.psi1(Z, mu, S, tmp) + p2.dpsi1_dmuS((tmp[:,None,:]*dL_dpsi2).sum(1).T*2., Z, mu, S, target_mu, target_S) return target_mu, target_S diff --git a/GPy/kern/kernpart.py b/GPy/kern/kernpart.py index 30a1cc3d..7de150e9 100644 --- a/GPy/kern/kernpart.py +++ b/GPy/kern/kernpart.py @@ -54,5 +54,3 @@ class kernpart(object): raise NotImplementedError def dK_dX(self,X,X2,target): raise NotImplementedError - - diff --git a/GPy/kern/white.py b/GPy/kern/white.py index be6aad45..d5701cd9 100644 --- a/GPy/kern/white.py +++ b/GPy/kern/white.py @@ -18,7 +18,8 @@ class white(kernpart): self.Nparam = 1 self.name = 'white' self._set_params(np.array([variance]).flatten()) - + self._psi1 = 0 # TODO: more elegance here + def _get_params(self): return self.variance @@ -81,4 +82,3 @@ class white(kernpart): def dpsi2_dmuS(self,dL_dpsi2,Z,mu,S,target_mu,target_S): pass -