diff --git a/GPy/kern/Matern32.py b/GPy/kern/Matern32.py
index 9831ae40..c175009d 100644
--- a/GPy/kern/Matern32.py
+++ b/GPy/kern/Matern32.py
@@ -106,7 +106,7 @@ class Matern32(kernpart):
         dK_dX = - np.transpose(3*self.variance*dist*np.exp(-np.sqrt(3)*dist)*ddist_dX,(1,0,2))
         target += np.sum(dK_dX*partial.T[:,:,None],0)
 
-    def dKdiag_dX(self,X,target):
+    def dKdiag_dX(self,partial,X,target):
         pass
 
     def Gram_matrix(self,F,F1,F2,lower,upper):
diff --git a/GPy/kern/Matern52.py b/GPy/kern/Matern52.py
index 2994fc45..26caad1c 100644
--- a/GPy/kern/Matern52.py
+++ b/GPy/kern/Matern52.py
@@ -104,7 +104,7 @@ class Matern52(kernpart):
         dK_dX = -  np.transpose(self.variance*5./3*dist*(1+np.sqrt(5)*dist)*np.exp(-np.sqrt(5)*dist)*ddist_dX,(1,0,2))
         target += np.sum(dK_dX*partial.T[:,:,None],0)
 
-    def dKdiag_dX(self,X,target):
+    def dKdiag_dX(self,partial,X,target):
         pass
 
     def Gram_matrix(self,F,F1,F2,F3,lower,upper):
diff --git a/GPy/kern/exponential.py b/GPy/kern/exponential.py
index 3c9cb192..366ddf3b 100644
--- a/GPy/kern/exponential.py
+++ b/GPy/kern/exponential.py
@@ -102,7 +102,7 @@ class exponential(kernpart):
         dK_dX = - np.transpose(self.variance*np.exp(-dist)*ddist_dX,(1,0,2))
         target += np.sum(dK_dX*partial.T[:,:,None],0)
 
-    def dKdiag_dX(self,X,target):
+    def dKdiag_dX(self,partial,X,target):
         pass
 
     def Gram_matrix(self,F,F1,lower,upper):
diff --git a/GPy/kern/product.py b/GPy/kern/product.py
index cefd7458..92522418 100644
--- a/GPy/kern/product.py
+++ b/GPy/kern/product.py
@@ -82,5 +82,27 @@ class product(kernpart):
         self.k1.dK_dX(partial*K2, X, X2, target)
         self.k2.dK_dX(partial*K1, X, X2, target)
 
-    def dKdiag_dX(self,X,target):
-        pass
+    def dKdiag_dX(self,partial,X,target):
+        target1 = np.zeros((X.shape[0],))
+        target2 = np.zeros((X.shape[0],))
+        self.k1.Kdiag(X,target1)
+        self.k2.Kdiag(X,target2)
+
+        self.k1.dKdiag_dX(partial*target2, X, target)
+        self.k2.dKdiag_dX(partial*target1, X, target)
+
+    def dKdiag_dtheta(self,partial,X,target):
+        """Compute the diagonal of the covariance matrix associated to X."""
+        target1 = np.zeros((X.shape[0],))
+        target2 = np.zeros((X.shape[0],))
+        self.k1.Kdiag(X,target1)
+        self.k2.Kdiag(X,target2)
+
+        k1_target = np.zeros(self.k1.Nparam)
+        k2_target = np.zeros(self.k2.Nparam)
+        self.k1.dKdiag_dtheta(partial*target2, X, k1_target)
+        self.k2.dKdiag_dtheta(partial*target1, X, k2_target)
+
+        target[:self.k1.Nparam] += k1_target
+        target[self.k1.Nparam:] += k2_target
+