diff --git a/GPy/likelihoods/EP.py b/GPy/likelihoods/EP.py
index 8307b6b4..685195ba 100644
--- a/GPy/likelihoods/EP.py
+++ b/GPy/likelihoods/EP.py
@@ -196,9 +196,8 @@ class EP(likelihood):
                 self.tau_tilde[i] = self.tau_tilde[i] + Delta_tau
                 self.v_tilde[i] = self.v_tilde[i] + Delta_v
                 #Posterior distribution parameters update
-                #LLT = LLT + np.outer(Kmn[:,i],Kmn[:,i])*Delta_tau
-                #L = jitchol(LLT)
-                cholupdate(L,Kmn[:,i]*np.sqrt(Delta_tau))
+                LLT = LLT + np.outer(Kmn[:,i],Kmn[:,i])*Delta_tau
+                L = jitchol(LLT)
                 V,info = linalg.lapack.flapack.dtrtrs(L,Kmn,lower=1)
                 Sigma_diag = np.sum(V*V,-2)
                 si = np.sum(V.T*V[:,i],-1)
@@ -251,6 +250,7 @@ class EP(likelihood):
         R = R0.copy()
         Diag = Diag0.copy()
         Sigma_diag = Knn_diag
+        RPT0 = np.dot(R0,P0.T)
 
         """
         Initial values - Cavity distribution parameters:
@@ -306,13 +306,7 @@ class EP(likelihood):
             Iplus_Dprod_i = 1./(1.+ Diag0 * self.tau_tilde)
             Diag = Diag0 * Iplus_Dprod_i
             P = Iplus_Dprod_i[:,None] * P0
-
-            #Diag = Diag0/(1.+ Diag0 * self.tau_tilde)
-            #P = (Diag / Diag0)[:,None] * P0
-            RPT0 = np.dot(R0,P0.T)
             L = jitchol(np.eye(M) + np.dot(RPT0,((1. - Iplus_Dprod_i)/Diag0)[:,None]*RPT0.T))
-            #L = jitchol(np.eye(M) + np.dot(RPT0,(1./Diag0 - Iplus_Dprod_i/Diag0)[:,None]*RPT0.T))
-            #L = jitchol(np.eye(M) + np.dot(RPT0,(1./Diag0 - Diag/(Diag0**2))[:,None]*RPT0.T))
             R,info = linalg.lapack.flapack.dtrtrs(L,R0,lower=1)
             RPT = np.dot(R,P.T)
             Sigma_diag = Diag + np.sum(RPT.T*RPT.T,-1)