diff --git a/examples/ishigami/data/sparse_solver_comparison.py b/examples/ishigami/data/sparse_solver_comparison.py
index 4c2be19b5da1aacc726062dfe77daf5a7f7031e9..8bf8e00ca4b32c46ba7df8651397f969669c735b 100644
--- a/examples/ishigami/data/sparse_solver_comparison.py
+++ b/examples/ishigami/data/sparse_solver_comparison.py
@@ -9,6 +9,7 @@ Created on Sat Sep 10 09:44:05 2022
 import numpy as np
 import joblib
 import os
+import scipy.stats as st
 
 import sys
 sys.path.append("../../../src/bayesvalidrox/")
@@ -70,7 +71,7 @@ def plot_seq_design_diagnostics(meta_model, util_funcs,
     # Plot the evolution of the diagnostic criteria of the
     # Sequential Experimental Design.
     for plotidx, plot in enumerate(plotList):
-        fig, ax = plt.subplots()
+        fig, ax = plt.subplots(figsize=(27, 15))
         seq_dict = seqList[plotidx]
         name_util = list(seq_dict.keys())
 
@@ -105,16 +106,19 @@ def plot_seq_design_diagnostics(meta_model, util_funcs,
 
             # BoxPlot
             def draw_plot(data, labels, edge_color, fill_color, idx):
-                pos = labels - (2*idx-1)
+                pos = labels - (4*idx-6)
                 bp = plt.boxplot(data, positions=pos, labels=labels,
-                                 patch_artist=True, sym='', widths=1.75)
+                                 patch_artist=True, sym='', widths=3)
+
+                ax.plot(pos, np.median(data, axis=0), lw=4, color=fill_color[idx])
+
                 elements = ['boxes', 'whiskers', 'fliers', 'means',
                             'medians', 'caps']
                 for element in elements:
-                    plt.setp(bp[element], color=edge_color[idx])
+                    plt.setp(bp[element], color=edge_color[idx], alpha=0.6)
 
                 for patch in bp['boxes']:
-                    patch.set(facecolor=fill_color[idx])
+                    patch.set(facecolor=fill_color[idx], alpha=0.6)
 
             if meta_model.ExpDesign.n_new_samples != 1:
                 step1 = meta_model.ExpDesign.n_new_samples
@@ -156,20 +160,24 @@ def plot_seq_design_diagnostics(meta_model, util_funcs,
                                     ls='--', lw=2)
 
                 # Plot each UtilFuncs
-                labels = np.arange(n_init_samples, n_total_samples+1, step1)
-                draw_plot(all_errors[:, ::step2], labels, edge_color,
+                # labels = np.arange(n_init_samples, n_total_samples+1, step1)
+                labels = np.array([10, 30, 50, 70, 90, 120, 150, 200])
+                indices = [0, 20, 40, 60, 80, 110, 140, 190]
+                draw_plot(all_errors[:, indices], labels, edge_color,
                           fill_color, idx)
+                # draw_plot(all_errors[:, ::step2], labels, edge_color,
+                #           fill_color, idx)
 
             plt.xticks(labels)
             # Set the major and minor locators
             # ax.xaxis.set_major_locator(ticker.AutoLocator())
             # ax.xaxis.set_minor_locator(ticker.AutoMinorLocator())
-            ax.xaxis.grid(True, which='major', linestyle='-')
-            ax.xaxis.grid(True, which='minor', linestyle='--')
+            # ax.xaxis.grid(True, which='major', linestyle='-')
+            # ax.xaxis.grid(True, which='minor', linestyle='--')
 
             # Shade
-            # for center in np.linspace(25, 200, 25):
-            #     ax.axvspan(center-3, center+2, alpha=0.5, color='grey')
+            for center in [30, 70, 120, 200]:#labels[::2]:
+                ax.axvspan(center-8, center+8, alpha=0.1, color='grey')
 
             # Legend
             legend_elements = []
@@ -182,9 +190,13 @@ def plot_seq_design_diagnostics(meta_model, util_funcs,
             if plot != 'BME' and plot != 'KLD':
                 plt.yscale('log')
             plt.autoscale(True)
-            plt.xlabel('\\# of training samples')
-            plt.ylabel(plot_label)
-            plt.title(plot)
+            # ax.yaxis.set_minor_locator(ticker.LogLocator(numticks=999, subs="auto"))
+            ax.yaxis.grid(True, which='minor', linestyle='--')
+            plt.xlabel('\\# of training samples', fontsize=f_size)
+            plt.ylabel(plot_label, fontsize=f_size)
+            # plt.title(plot)
+            plt.xticks(fontsize=f_size)
+            plt.yticks(fontsize=f_size)
 
             if save_fig:
                 # save the current figure
@@ -201,8 +213,9 @@ def plot_seq_design_diagnostics(meta_model, util_funcs,
 if __name__ == "__main__":
     # Set variables
     model_name = 'Ishigami'
-    solvers = ['BCS', 'FastARD', 'OMP', 'OLS']
-    path = f'/home/farid/bwSyncShare/Scientific_LH2/Promotion/dissertation/surrogate/data/'
+    solvers = ['FastLaplace', 'FastARD', 'OMP', 'OLS']
+    path = f'/home/farid/bwSyncShare/Scientific_LH2/Promotion/dissertation/surrogate/data-ishigami/'
+    f_size = 45
 
     all_loo_errors = {}
     all_valid_errors = {}