From eb5eefdcf0e5e1d9242aebfb09f0db820a30fafb Mon Sep 17 00:00:00 2001
From: William Gatt <williamjgatt@gmail.com>
Date: Wed, 5 Jul 2023 10:07:38 +0000
Subject: [PATCH] Variable number of parameter MCMC convergence plots

This change allows the user to change the number of parameters for which the MCMC convergence diagnostics are shown. The user needs to set a value for options_.convergence.brooksgelman.plotrows for an arbitrary number of parameters to be shown, but defaulting to 3 if no input is supplied.
---
 matlab/convergence_diagnostics/McMCDiagnostics.m | 12 +++++++-----
 matlab/default_option_values.m                   |  2 ++
 2 files changed, 9 insertions(+), 5 deletions(-)

diff --git a/matlab/convergence_diagnostics/McMCDiagnostics.m b/matlab/convergence_diagnostics/McMCDiagnostics.m
index 18887db3b..16f58ed85 100644
--- a/matlab/convergence_diagnostics/McMCDiagnostics.m
+++ b/matlab/convergence_diagnostics/McMCDiagnostics.m
@@ -48,6 +48,7 @@ MAX_nruns = ceil(options_.MaxNumberOfBytes/(npar+2)/8);
 
 load_last_mh_history_file(MetropolisFolder, ModelName);
 NumberOfMcFilesPerBlock = record.LastFileNumber;
+npardisp = options_.convergence.brooksgelman.plotrows;
 
 % Check that the declared number of blocks is consistent with informations saved in mh-history files.
 if ~isequal(nblck,record.Nblck)
@@ -284,12 +285,13 @@ 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/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: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 +308,7 @@ for i = 1:pages
                 plt2 = UDIAG(:,6,k);
                 namnam  = [nam , ' (m3)'];
             end
-            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)
@@ -339,8 +341,8 @@ if reste
     if reste == 1
         nr = 3;
         nc = 1;
-    elseif reste == 2
-        nr = 2;
+    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)');
diff --git a/matlab/default_option_values.m b/matlab/default_option_values.m
index 0c3616b91..e8442aa00 100644
--- a/matlab/default_option_values.m
+++ b/matlab/default_option_values.m
@@ -766,6 +766,8 @@ options_.convergence.geweke.geweke_interval=[0.2 0.5];
 %Raftery/Lewis convergence diagnostics;
 options_.convergence.rafterylewis.indicator=false;
 options_.convergence.rafterylewis.qrs=[0.025 0.005 0.95];
+%Brooks Gelman convergence plots
+options_.convergence.brooksgelman.plotrows=3;
 
 %tolerance for Modified Harmonic Mean estimator
 options_.marginal_data_density.harmonic_mean.tolerance = 0.01;