diff --git a/pybasicbayes/distributions/gaussian.py b/pybasicbayes/distributions/gaussian.py
index 66ceed683f92f8e331bf31fdb5d589427f02da77..f283986f623439653a8dd78270a17a3584df8bab 100644
--- a/pybasicbayes/distributions/gaussian.py
+++ b/pybasicbayes/distributions/gaussian.py
@@ -344,10 +344,11 @@ class Gaussian(
             - D*self.kappa_0/self.kappa_mf - self.kappa_0*self.nu_mf*
             np.dot(self.mu_mf -
                 self.mu_0,np.linalg.solve(self.sigma_mf,self.mu_mf - self.mu_0))) \
-            + invwishart_log_partitionfunction(self.sigma_0,self.nu_0) \
+            - invwishart_log_partitionfunction(self.sigma_0,self.nu_0) \
             + (self.nu_0 - D - 1)/2*loglmbdatilde - 1/2*self.nu_mf \
             * np.linalg.solve(self.sigma_mf,self.sigma_0).trace()
 
+
         return p_avgengy + q_entropy
 
     def expected_log_likelihood(self, x=None, stats=None):
diff --git a/pybasicbayes/util/stats.py b/pybasicbayes/util/stats.py
index 7b704570de5102826762250e7d201df33aa38ffa..9ea9b479a449c970aaffc7cab0e6f1a124461c6d 100644
--- a/pybasicbayes/util/stats.py
+++ b/pybasicbayes/util/stats.py
@@ -317,6 +317,7 @@ def invwishart_entropy(sigma,nu,chol=None):
     return invwishart_log_partitionfunction(sigma,nu,chol)-(nu-D-1)/2*Elogdetlmbda + nu*D/2
 
 def invwishart_log_partitionfunction(sigma,nu,chol=None):
+    # In Bishop B.79 notation, this is -log B(W, nu), where W = sigma^{-1}
     D = sigma.shape[0]
     chol = np.linalg.cholesky(sigma) if chol is None else chol
     return -1*(nu*np.log(chol.diagonal()).sum() - (nu*D/2*np.log(2) + D*(D-1)/4*np.log(np.pi) \