diff --git a/matlab/convergence_diagnostics/McMCDiagnostics.m b/matlab/convergence_diagnostics/McMCDiagnostics.m
index 18887db3b4229fed99e4ed015b938abdfc66c4f2..20e906992fec66ffb0efaab672dba08ac9cd3275 100644
--- a/matlab/convergence_diagnostics/McMCDiagnostics.m
+++ b/matlab/convergence_diagnostics/McMCDiagnostics.m
@@ -49,6 +49,19 @@ MAX_nruns = ceil(options_.MaxNumberOfBytes/(npar+2)/8);
 load_last_mh_history_file(MetropolisFolder, ModelName);
 NumberOfMcFilesPerBlock = record.LastFileNumber;
 
+% Check whether umber of parameters displayed in MCMC convergence
+% diagnostics figures has been declared. If not, default to 3.
+if isfield(options_,'nparMCMCdisp')
+    if options_.nparMCMCdisp == 0
+        warning('Cannot plot 0 variables, reseting options_.nparMCMCdisp to 3 [default]');
+        options_.nparMCMCdisp = 3;
+    end
+    npardisp = options_.nparMCMCdisp;
+else
+    npardisp = 3;
+end
+%------------------
+
 % Check that the declared number of blocks is consistent with informations saved in mh-history files.
 if ~isequal(nblck,record.Nblck)
     disp(['Estimation::mcmc::diagnostics: The number of MCMC chains you declared (' num2str(nblck) ') is inconsistent with the information available in the mh-history files (' num2str(record.Nblck) ' chains)!'])
@@ -284,12 +297,14 @@ end
 UDIAG(:,[2 4 6],:) = UDIAG(:,[2 4 6],:)/nblck;
 skipline()
 clear pmet temp moyenne CSUP CINF csup cinf n linea iter tmp;
-pages = floor(npar/3);
+% pages = floor(npar/3);
+pages = floor(npar/npardisp); % changed from 3 to npardisp
 k = 0;
 for i = 1:pages
     h=dyn_figure(options_.nodisplay,'Name','MCMC univariate convergence diagnostic (Brooks and Gelman,1998)');
     boxplot = 1;
-    for j = 1:3 % Loop over parameters
+%     for j = 1:3 % Loop over parameters
+    for j = 1:npardisp % Loop over parameters  %npardisp instead of 3
         k = k+1;
         [nam,namtex] = get_the_name(k,TeX,M_,estim_params_,options_);
         for crit = 1:3% Loop over criteria
@@ -306,7 +321,8 @@ for i = 1:pages
                 plt2 = UDIAG(:,6,k);
                 namnam  = [nam , ' (m3)'];
             end
-            subplot(3,3,boxplot);
+%           subplot(3,3,boxplot);
+            subplot(npardisp,3,boxplot);  %Added more rows to display more variables
             plot(xx,plt1,'-b');     % Pooled
             hold on;
             plot(xx,plt2,'-r');     % Within (mean)
@@ -342,6 +358,9 @@ if reste
     elseif reste == 2
         nr = 2;
         nc = 3;
+    else % Conditional for additional rows (variables) when not towards the end of the loop
+        nr = npardisp;
+        nc = 3;
     end
     h = dyn_figure(options_.nodisplay,'Name','MCMC univariate convergence diagnostic (Brooks and Gelman, 1998)');
     boxplot = 1;