diff --git a/GPy/kern/linear.py b/GPy/kern/linear.py
index 7d817f62..ef6b72bb 100644
--- a/GPy/kern/linear.py
+++ b/GPy/kern/linear.py
@@ -81,6 +81,13 @@ class linear(kernpart):
             self._K_computations(X, X2)
             target += np.sum(self._dot_product*dL_dK)
 
+    def dKdiag_dtheta(self,dL_dKdiag, X, target):
+        tmp = dL_dKdiag[:,None]*X**2
+        if self.ARD:
+            target += tmp.sum(0)
+        else:
+            target += tmp.sum()
+
     def dK_dX(self,dL_dK,X,X2,target):
         target += (((X2[:, None, :] * self.variances)) * dL_dK[:,:, None]).sum(0)
 
@@ -92,13 +99,6 @@ class linear(kernpart):
         self._psi_computations(Z,mu,S)
         target += np.sum(self.variances*self.mu2_S,1)
 
-    def dKdiag_dtheta(self,dL_dKdiag, X, target):
-        tmp = dL_dKdiag[:,None]*X**2
-        if self.ARD:
-            target += tmp.sum(0)
-        else:
-            target += tmp.sum()
-
     def dpsi0_dtheta(self,dL_dpsi0,Z,mu,S,target):
         self._psi_computations(Z,mu,S)
         tmp = dL_dpsi0[:, None] * self.mu2_S
@@ -134,6 +134,7 @@ class linear(kernpart):
         self._psi_computations(Z,mu,S)
         psi2 = self.ZZ*np.square(self.variances)*self.mu2_S[:, None, None, :]
         target += psi2.sum(-1)
+        #TODO: this could be faster using np.tensordot
 
     def dpsi2_dtheta(self,dL_dpsi2,Z,mu,S,target):
         self._psi_computations(Z,mu,S)
diff --git a/GPy/testing/bgplvm_tests.py b/GPy/testing/bgplvm_tests.py
index 80e6fecd..dda92b90 100644
--- a/GPy/testing/bgplvm_tests.py
+++ b/GPy/testing/bgplvm_tests.py
@@ -58,6 +58,7 @@ class BGPLVMTests(unittest.TestCase):
         m.randomize()
         self.assertTrue(m.checkgrad())
 
+    @unittest.skip('psi2 cross terms are NotImplemented for this combination')
     def test_linear_bias_kern(self):
         N, M, Q, D = 10, 3, 2, 4
         X = np.random.rand(N, Q)