diff --git a/pyfstat/mcmc_based_searches.py b/pyfstat/mcmc_based_searches.py
index 13d08f64aea48a1ac27745ea0e1faaa0b558fc52..7f5fe106e3481e312f99db0bdac125990a684536 100644
--- a/pyfstat/mcmc_based_searches.py
+++ b/pyfstat/mcmc_based_searches.py
@@ -993,18 +993,20 @@ class MCMCSearch(core.BaseSearchClass):
                 if burnin_idx and add_det_stat_burnin:
                     burn_in_vals = lnl[:, :burnin_idx].flatten()
                     try:
-                        axes[-1].hist(burn_in_vals[~np.isnan(burn_in_vals)],
-                                      bins=50, histtype='step', color='C3')
+                        twoF_burnin = (burn_in_vals[~np.isnan(burn_in_vals)]
+                                       - self.likelihoodcoef)
+                        axes[-1].hist(twoF_burnin, bins=50, histtype='step',
+                                      color='C3')
                     except ValueError:
                         logging.info('Det. Stat. hist failed, most likely all '
                                      'values where the same')
                         pass
                 else:
-                    burn_in_vals = []
+                    twoF_burnin = []
                 prod_vals = lnl[:, burnin_idx:].flatten()
                 try:
-                    axes[-1].hist(prod_vals[~np.isnan(prod_vals)], bins=50,
-                                  histtype='step', color='k')
+                    twoF = prod_vals[~np.isnan(prod_vals)]-self.likelihoodcoef
+                    axes[-1].hist(twoF, bins=50, histtype='step', color='k')
                 except ValueError:
                     logging.info('Det. Stat. hist failed, most likely all '
                                  'values where the same')
@@ -1014,7 +1016,7 @@ class MCMCSearch(core.BaseSearchClass):
                 else:
                     axes[-1].set_xlabel(r'$\widetilde{2\mathcal{F}}$')
                 axes[-1].set_ylabel(r'$\textrm{Counts}$')
-                combined_vals = np.append(burn_in_vals, prod_vals)
+                combined_vals = np.append(twoF_burnin, twoF)
                 if len(combined_vals) > 0:
                     minv = np.min(combined_vals)
                     maxv = np.max(combined_vals)