MIT License
ParaMonte: plain powerful parallel Monte Carlo library.
Copyright (C) 2012-present, The Computational Data Science Lab
https://github.com/cdslaborg/paramonte
References

Sampling a 4-dimensional MultiVariate Normal distribution (MVN) via the ParaMonte library's ParaDRAM routine

NOTE
If you are viewing an HTML version of this MATLAB live script on the web, you can download the corresponding MATLAB live script *.mlx file to this HTML page at,
https://github.com/cdslaborg/paramontex/blob/fbeca6745684c798ff28c1bf57cfae0c190db478/MATLAB/mlx/sampling_multivariate_normal_distribution_via_paradram/sampling_multivariate_normal_distribution_via_paradram.mlx
This code, along with other MATLAB live scripts, is located at,
https://github.com/cdslaborg/paramontex/blob/fbeca6745684c798ff28c1bf57cfae0c190db478/MATLAB/mlx
Once you download the file, open it in MATLAB to view and interact with its contents, which is the same as what you see on this page.
NOTE
This ParaDRAM sampler example is also available in Python language as a Jupyter Notebook, on this page,
https://nbviewer.jupyter.org/github/cdslaborg/paramontex/blob/fbeca6745684c798ff28c1bf57cfae0c190db478/Python/Jupyter/sampling_multivariate_normal_distribution_via_paradram/sampling_multivariate_normal_distribution_via_paradram.ipynb
This Jupyter Notebook, along with other Python Jupyter Notebooks, is located at,
https://github.com/cdslaborg/paramontex/blob/fbeca6745684c798ff28c1bf57cfae0c190db478/Python/Jupyter
IMPORTANT - macOS users
If you are a macOS (Darwin) user and you have downloaded the prebuilt ParaMonte libraries from the GitHub release page to use on your macOS system, read this troubleshooting page before continuing.

The logic of Monte Carlo sampling is very simple

==========================================================
Suppose we have a mathematical objective function defined on an ndim-dimensional domain. Typically, this domain is defined on the space of real numbers. In general, understanding and visualizing the structure of such objective functions is an extremely difficult task, if not impossible. As a better alternative, you employ Monte Carlo techniques that can randomly explore the structure of the objective function and find the function's extrema.
In the following example, one of the simplest such objective functions, the Multivariate Normal Distribution (MVN), is constructed and sampled using the ParaMonte library samplers, here, the ParaDRAM sampler (Delayed-Rejection Adaptive Metropolis-Hastings Markov Chain Monte Carlo sampler).
Suppose we want to sample random points from a MultiVariate Normal (MVN) Probability Density Function (PDF). The PDF equation of MVN is the following,

Setting up the path to the ParaMonte::MATLAB library

================================================================
First, we will clean up the MATLAB environment and make sure the path to the ParaMonte library is in MATLAB's path list.
clc;
clear all;
close all;
format compact; format long;
%%%%%%%%%%%% IMPORTANT %%%%%%%%%%%%%
% Set the path to the ParaMonte library:
% Change the following path to the ParaMonte library root directory,
% otherwise, make sure the path to the ParaMonte library is already added
% to MATLAB's path list.
pmlibRootDir = './';
addpath(genpath(pmlibRootDir),"-begin");
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% change MATLAB's working directory to the folder containing this script
% if MATLAB Live Scripts did not create a temporary folder, we would not
% have all of these problems!
try
setwdFilePath = websave("setwd.m","https://github.com/cdslaborg/paramontex/blob/fbeca6745684c798ff28c1bf57cfae0c190db478/MATLAB/mlx/setwd.m");
run(setwdFilePath); % This is a MATLAB script that you can download from the same GitHub location given in the above.
catch % alas, we will have to run the simulations in MATLAB Live Script's temporary folder
filePath = mfilename('fullpath');
[currentDir,fileName,fileExt] = fileparts(filePath); cd(currentDir);
cd(fileparts(mfilename('fullpath'))); % Change working directory to source code directory.
end

Constructing the objective function

==========================================
The following MATLAB function getLogFunc() returns the natural logarithm of the Probability Density Function (PDF) of the MultiVariate Normal (MVN) distribution,
NDIM = 4; % the number of dimensions of the domain of the objective function
MEAN = [-10; 15.; 20.; 0.0]; % This is the mean of the multivariate normal (MVN) distribution. THIS IS A COLUMN VECTOR.
COVMAT = [ 1.0,.45,-.3,0.0 ...
; .45,1.0,0.3,-.2 ...
; -.3,0.3,1.0,0.6 ...
; 0.0,-.2,0.6,1.0 ...
]; % This is the covariance matrix of the MVN distribution.
INVCOV = inv(COVMAT); % This is the inverse of the covariance matrix of the MVN distribution.
MVN_COEF = NDIM * log( 1. / sqrt(2.*pi) ) + log( sqrt(det(INVCOV)) ); % This is the log of the coefficient used in the definition of the MVN.
getLogFunc = @(point) MVN_COEF - 0.5 * dot( (MEAN-point)' , INVCOV*(MEAN-point) ); % return the logarithm of the objective function: log(MVN)
NOTE
Since the mathematical objective functions (e.g., probability density functions) can take extremely small or large values, we often work with their natural logarithms instead. This is the reason behind the naming convention used in the ParaMonte library for the user's objective functions: getLogFunc, implying that the user must provide a function that returns the natural logarithm of the target objective function.
See MATLAB Live Script for an in-depth discussion of why we need to work with the logarithm of mathematical objective functions in optimization and sampling problems.
IMPORTANT
Note that in MATLAB all vectors and arrays are, by default, column-major, and so is the input value, point, to getLogFunc(). This means that the size of point is (ndim,1). Thus, all other vectors that interact with point, for example, the mean vector, must be also of the same size as point: (ndim,1).

Running a ParaDRAM simulation on a single processor

==================================================================
We will sample random points from this objective function by calling the ParaDRAM sampler (Parallel Delayed-Rejection Adaptive Metropolis-Hastings Markov Chain Monte Carlo sampler) of the ParaMonte library.
Note: To run the ParaMonte samplers in parallel, visit the ParaMonte library's documentation website.
We will first create an instance of the paramonte class,
pm = paramonte();

Checking the ParaMonte library's version

==========================================
Before runnnig any simulations, it is good to check which version of the ParaMonte library you have and are working with,
pm.version.interface.get() % get the ParaMonte::MATLAB (interface) version
ans = "ParaMonte MATLAB Interface Version 2.3.0"
pm.version.interface.dump()
ans = "2.3.0"
pm.version.kernel.get() % the ParaMonte::Kernel version (the version of the ParaDRAM sampler)
ans = "ParaMonte MATLAB Kernel Version 1.5.0"
pm.version.kernel.dump()
ans = "1.5.0"
We can also verify the existence of the essential components of the current installation of the ParaMonte library on the system,
pm.verify();
::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: :::: :::: _/_/_/_/ _/_/ _/_/ _/ _/ _/_/_/_/_/ _/ _/ _/ _/_/_/_/ _/ /_/_/ _/_/_/_/ _/ _/ _/ _/_/ _/_/_/ _/_/_/ _/_/_/ _/_/_/ _/ _/ _/_/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/_/_/_/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/ _/_/_/ _/_/_/_/ _/ _/_/_/_/ _/_/_/ _/_/_/ _/_/ _/ _/ _/_/ _/_/_/ ParaMonte plain powerful parallel Monte Carlo library Version 2.3.0 :::: :::: ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: ParaMonte - NOTE: Intel MPI library for 64-bit architecture detected at: ParaMonte - NOTE: ParaMonte - NOTE: "C:\Program Files (x86)\IntelSWTools\compilers_and_libraries_2020.1.216\windows\mpi\intel64\bin" ParaMonte - NOTE: ParaMonte - NOTE: To perform ParaMonte simulations in parallel on a single node, run the ParaMonte - NOTE: following two commands, in the form and order specified, on a MATLAB-aware ParaMonte - NOTE: mpiexec-aware command-line interface, such as the Intel Parallel Studio's command-prompt: ParaMonte - NOTE: ParaMonte - NOTE: "C:\Program Files (x86)\IntelSWTools\compilers_and_libraries_2020.1.216\windows\mpi\intel64\bin\mpivars.bat" ParaMonte - NOTE: ParaMonte - NOTE: mpiexec -localonly -n 3 matlab -batch "main_mpi" ParaMonte - NOTE: ParaMonte - NOTE: where, ParaMonte - NOTE: ParaMonte - NOTE: 0. the first command defines the essential environment variables and, ParaMonte - NOTE: the second command runs in the simulation in parallel, in which, ParaMonte - NOTE: 1. you should replace the default number 3 with the number of ParaMonte - NOTE: processors you wish to assign to your simulation task, ParaMonte - NOTE: 2. main_mpi.m is the MATLAB file which serves as the entry point to ParaMonte - NOTE: your simulation, where you call the ParaMonte sampler routines. ParaMonte - NOTE: ATTN: Notice the MATLAB filename must appear without the file extension. ParaMonte - NOTE: ATTN: Notice the MATLAB filename must be enclosed with quotation marks. ParaMonte - NOTE: 3. -localonly indicates a parallel simulation on only a single node (this ParaMonte - NOTE: flag will obviate the need for MPI library credentials registration). ParaMonte - NOTE: For more information, visit: ParaMonte - NOTE: ParaMonte - NOTE: https://software.intel.com/en-us/get-started-with-mpi-for-windows ParaMonte - NOTE: ParaMonte - NOTE: Note that the above two commands must be executed on a command-line that recognizes ParaMonte - NOTE: both MATLAB and mpiexec applications, such as the Intel Parallel Studio's command-prompt. ParaMonte - NOTE: For more information, in particular, on how to register to run Hydra services ParaMonte - NOTE: for multi-node simulations on Windows servers, visit: ParaMonte - NOTE: ParaMonte - NOTE: https://www.cdslab.org/paramonte ParaMonte - NOTE: To check for the ParaMonte or the MPI library installations status or, ParaMonte - NOTE: to display the above messages in the future, use the following ParaMonte - NOTE: commands on the MATLAB command-line: ParaMonte - NOTE: ParaMonte - NOTE: pm = paramonte(); ParaMonte - NOTE: pm.verify();

Setting up the ParaDRAM simulation specifications

====================================================
Now, we create an instance of the ParaDRAM sampler class,
pmpd = pm.ParaDRAM();
The simplest scenario is to run the simulation with the default specifications that are appropriately determined by the ParaDRAM sampler. However, for the purpose of clarity reproducibility, we will specify a few simulation specs,
For a complete list of all ParaDRAM simulation specifications, visit the ParaDRAM simulation specifications webpage.

The simulation output file names and storage path

We will store all output files in a directory with the same name as this MATLAB Live Script's name.
By default, all ParaDRAM simulation output files are named properly (see this page for a review of ParaDRAM simulation output files). In general, it is a good habit to specify an output file prefix for any ParaDRAM simulation. This way, the restart of an incomplete simulation, if it happens for some reason, would be seamless and doable by just rerunning the simulation, without any change of any sort to any parts of the codes or the simulation.
We will prefix all simulation output files the same as the filename of this script,
pmpd.spec.outputFileName = "./out/mvn_serial";
pmpd.spec.helpme("outputFileName");
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#outputfilename
The above outputFileName will result in the generation of a folder named ./out/, which will keep the simulation output files, all of which are prefixed with mvn_serial.
NOTE
If the specified prefix ends with a forward slash / on Linux/macOS or with a backslash \ on Windows, then the prefix provided will serve as the directory in which the output files will be stored.

Ensuring the reproducibility of the results

We will also initialize the random seed to a fixed number to ensure the reproducibility of the simulation results in the future,
pmpd.spec.randomSeed = 3751;
pmpd.spec.helpme("randomSeed");
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#randomseed
For the sake of illustration, let's also ask for a modest refinement (decorrelation) of the final sample,
pmpd.spec.sampleRefinementCount = 1;
pmpd.spec.helpme("sampleRefinementCount");
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#samplerefinementcount

Setting the output sample size to a non-default value

Not so important, but we will also specify a chainSize here, 50000, the number of unique points to be sampled from the objective function. This number is smaller than the default 100,000 unique sample points of the sampler.
pmpd.spec.chainSize = 50000;
pmpd.spec.helpme("chainSize");
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#chainsize

Requesting an ASCII-formatted output restart file (for illustration purposes)

For the sake of illustration, we will also request an ASCII-format restart output file,
pmpd.spec.restartFileFormat = "ascii";
pmpd.spec.helpme("restartFileFormat");
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#restartfileformat
Note, however, that you do not want to request an ACSII-formatted output restart file in your real production simulation as it will slow down your simulation. Just do not set this property and it will be automatically taken care of.

The default domain of the objective function

The default starting point of the sampler for exploring the objective function is the center of the domain of the objective function. The default domain of the objective function chosen by the sampler is along each dimension of the objective function. This is a fine default assumption for this particular simulation and in most other scenarios. If needed, you can change them via the following attributes,
pmpd.spec.domainLowerLimitVec = -1.e300 * ones(NDIM,1); % set the lower limits of the domain of the objective function to negative infinity (This is the default value)
pmpd.spec.domainUpperLimitVec = +1.e300 * ones(NDIM,1); % set the upper limits of the domain of the objective function to negative infinity (This is the default value)
The ParaDRAM sampler guarantees to not sample any points outside the hypercube defined by the above limits.

Calling the ParaDRAM sampler

===============================
Note that none of the above specification settings were necessary, but is in general a good habit to define them prior to the simulation.
We now call the ParaDRAM sampler routine to run the sampling simulation. For the same of this simulation, we will also make the output files overwritable (before doing, The ParaDRAM sampler method takes only two arguments. We can call it like the following,
pmpd.spec.overwriteRequested = true; % overwrite the simulation output files if they already exist
pmpd.runSampler ( NDIM ... This is the number of dimensions of the objective function
, getLogFunc ... This is the objective function
);
ParaDRAM - NOTE: Running the ParaDRAM sampler in serial mode... ParaDRAM - NOTE: To run the ParaDRAM sampler in parallel mode visit: cdslab.org/pm ************************************************************************************************************************************ ************************************************************************************************************************************ **** **** **** **** **** ParaMonte **** **** Plain Powerful Parallel **** **** Monte Carlo Library **** **** **** **** Version 1.5.0 **** **** **** **** Build: Thu Dec 17 06:20:11 2020 **** **** **** **** Department of Physics **** **** Computational & Data Science Lab **** **** Data Science Program, College of Science **** **** The University of Texas at Arlington **** **** **** **** originally developed at **** **** **** **** Multiscale Modeling Group **** **** Center for Computational Oncology (CCO) **** **** Oden Institute for Computational Engineering and Sciences **** **** Department of Aerospace Engineering and Engineering Mechanics **** **** Department of Neurology, Dell-Seton Medical School **** **** Department of Biomedical Engineering **** **** The University of Texas at Austin **** **** **** **** For questions and further information, please contact: **** **** **** **** Amir Shahmoradi **** **** **** **** shahmoradi@utexas.edu **** **** amir.shahmoradi@uta.edu **** **** ashahmoradi@gmail.com **** **** **** **** cdslab.org/pm **** **** **** **** https://www.cdslab.org/paramonte/ **** **** **** **** **** ************************************************************************************************************************************ ************************************************************************************************************************************ ************************************************************************************************************************************ **** **** **** Setting up the ParaDRAM simulation environment **** **** **** ************************************************************************************************************************************ ParaDRAM - NOTE: Interfacing MATLAB with ParaDRAM... ParaDRAM - NOTE: Variable outputFileName detected among the input variables to ParaDRAM: ParaDRAM - NOTE: ./out/mvn_serial ParaDRAM - NOTE: ParaDRAM - NOTE: Absolute path to the current working directory: ParaDRAM - NOTE: D:\temp\libparamonte_MATLAB ParaDRAM - NOTE: ParaDRAM - NOTE: Generating the requested directory for the ParaDRAM output files: ParaDRAM - NOTE: .\out\ ParaDRAM - NOTE: ParaDRAM - NOTE: ParaDRAM output files will be prefixed with: ParaDRAM - NOTE: .\out\mvn_serial ParaDRAM - NOTE: Searching for previous runs of ParaDRAM... ParaDRAM - NOTE: No pre-existing ParaDRAM run detected. ParaDRAM - NOTE: Starting a fresh ParaDRAM run... ParaDRAM - NOTE: Generating the output report file: ParaDRAM - NOTE: .\out\mvn_serial_process_1_report.txt ParaDRAM - NOTE: Running the simulation in serial on 1 process. ParaDRAM - NOTE: Please see the output report and progress files for further realtime simulation details. Accepted/Total Func. Call Dynamic/Overall Acc. Rate Elapsed/Remained Time [s] ========================= ========================= ========================= 89 / 999 0.0884 / 0.0885 0.0860 / 48.2286 228 / 1999 0.1338 / 0.1112 0.1270 / 27.7239 394 / 2999 0.1618 / 0.1281 0.1740 / 21.9072 574 / 3999 0.1800 / 0.1410 0.2190 / 18.8577 784 / 4999 0.2134 / 0.1555 0.2610 / 16.3844 979 / 5999 0.2019 / 0.1633 0.3040 / 15.2221 1226 / 6999 0.2412 / 0.1744 0.3500 / 13.9241 1442 / 7999 0.2191 / 0.1800 0.3880 / 13.0655 1671 / 8999 0.2320 / 0.1858 0.4320 / 12.4944 1909 / 9999 0.2509 / 0.1923 0.4750 / 11.9661 2170 / 10999 0.2605 / 0.1985 0.5180 / 11.4175 2423 / 11999 0.2614 / 0.2037 0.5600 / 10.9959 2690 / 12999 0.2685 / 0.2087 0.6100 / 10.7283 2989 / 13999 0.2834 / 0.2140 0.6570 / 10.3333 3284 / 14999 0.2977 / 0.2196 0.7000 / 9.9577 3571 / 15999 0.2643 / 0.2224 0.7470 / 9.7123 3850 / 16999 0.2796 / 0.2258 0.7850 / 9.4098 4120 / 17999 0.2704 / 0.2283 0.8220 / 9.1537 4429 / 18999 0.3000 / 0.2320 0.8660 / 8.9105 4719 / 19999 0.2847 / 0.2347 0.9080 / 8.7127 5012 / 20999 0.2906 / 0.2373 0.9470 / 8.5003 5276 / 21999 0.2793 / 0.2392 1.0000 / 8.4769 5570 / 22999 0.2888 / 0.2414 1.0370 / 8.2718 5850 / 23999 0.2890 / 0.2434 1.0770 / 8.1281 6119 / 24999 0.2757 / 0.2447 1.1250 / 8.0677 6384 / 25999 0.2788 / 0.2460 1.1700 / 7.9935 6693 / 26999 0.3042 / 0.2481 1.2170 / 7.8746 6993 / 27999 0.2936 / 0.2498 1.2560 / 7.7244 7296 / 28999 0.2939 / 0.2513 1.2940 / 7.5739 7629 / 29999 0.3288 / 0.2539 1.3340 / 7.4090 7940 / 30999 0.3063 / 0.2556 1.3820 / 7.3208 8276 / 31999 0.3233 / 0.2577 1.4200 / 7.1590 8549 / 32999 0.2812 / 0.2584 1.4610 / 7.0839 8845 / 33999 0.2956 / 0.2595 1.5050 / 7.0026 9166 / 34999 0.3070 / 0.2608 1.5520 / 6.9141 9477 / 35999 0.3025 / 0.2620 1.5940 / 6.8158 9783 / 36999 0.3006 / 0.2630 1.6360 / 6.7254 10068 / 37999 0.3002 / 0.2640 1.6800 / 6.6633 10341 / 38999 0.2832 / 0.2645 1.7210 / 6.6002 10635 / 39999 0.3040 / 0.2655 1.7690 / 6.5479 10908 / 40999 0.2756 / 0.2657 1.8110 / 6.4902 11218 / 41999 0.3037 / 0.2666 1.8550 / 6.4130 11542 / 42999 0.3219 / 0.2679 1.8920 / 6.3042 11836 / 43999 0.2979 / 0.2686 1.9340 / 6.2360 12136 / 44999 0.3044 / 0.2694 1.9800 / 6.1775 12459 / 45999 0.3108 / 0.2703 2.0240 / 6.0986 12753 / 46999 0.3004 / 0.2709 2.0670 / 6.0370 13042 / 47999 0.3020 / 0.2716 2.1120 / 5.9849 13312 / 48999 0.2794 / 0.2718 2.1520 / 5.9309 13634 / 49999 0.3156 / 0.2726 2.1910 / 5.8441 13959 / 50999 0.3196 / 0.2735 2.2370 / 5.7758 14273 / 51999 0.3088 / 0.2742 2.2940 / 5.7422 14569 / 52999 0.2964 / 0.2746 2.3450 / 5.7029 14850 / 53999 0.2922 / 0.2750 2.3940 / 5.6666 15148 / 54999 0.3050 / 0.2755 2.4430 / 5.6208 15432 / 55999 0.2925 / 0.2758 2.4860 / 5.5687 15719 / 56999 0.2989 / 0.2762 2.5320 / 5.5219 16010 / 57999 0.3099 / 0.2768 2.5780 / 5.4732 16335 / 58999 0.3139 / 0.2774 2.6280 / 5.4161 16641 / 59999 0.3047 / 0.2779 2.6740 / 5.3604 16968 / 60999 0.3319 / 0.2788 2.7210 / 5.2970 17313 / 61999 0.3274 / 0.2796 2.7650 / 5.2203 17627 / 62999 0.3119 / 0.2801 2.8090 / 5.1589 17938 / 63999 0.3251 / 0.2808 2.8490 / 5.0922 18246 / 64999 0.3112 / 0.2812 2.8930 / 5.0348 18556 / 65999 0.3141 / 0.2817 2.9360 / 4.9752 18900 / 66999 0.3404 / 0.2826 2.9750 / 4.8954 19189 / 67999 0.2948 / 0.2828 3.0180 / 4.8459 19520 / 68999 0.3317 / 0.2835 3.0620 / 4.7812 19856 / 69999 0.3245 / 0.2841 3.1130 / 4.7259 20174 / 70999 0.3123 / 0.2845 3.1570 / 4.6674 20494 / 71999 0.3184 / 0.2850 3.2070 / 4.6172 20828 / 72999 0.3271 / 0.2855 3.2490 / 4.5506 21124 / 73999 0.3090 / 0.2859 3.2910 / 4.4987 21445 / 74999 0.3237 / 0.2864 3.3400 / 4.4474 21777 / 75999 0.3241 / 0.2869 3.3880 / 4.3908 22082 / 76999 0.3144 / 0.2872 3.4300 / 4.3365 22387 / 77999 0.3166 / 0.2876 3.4760 / 4.2874 22704 / 78999 0.3140 / 0.2879 3.5310 / 4.2452 23036 / 79999 0.3218 / 0.2883 3.5800 / 4.1904 23370 / 80999 0.3169 / 0.2887 3.6350 / 4.1421 23693 / 81999 0.3343 / 0.2893 3.6940 / 4.1016 24007 / 82999 0.3237 / 0.2897 3.7450 / 4.0548 24335 / 83999 0.3297 / 0.2901 3.7980 / 4.0056 24645 / 84999 0.3102 / 0.2904 3.8420 / 3.9527 24953 / 85999 0.3093 / 0.2906 3.8930 / 3.9077 25259 / 86999 0.3136 / 0.2909 3.9430 / 3.8621 25600 / 87999 0.3275 / 0.2913 4.0020 / 3.8144 25910 / 88999 0.3196 / 0.2916 4.0540 / 3.7692 26232 / 89999 0.3266 / 0.2920 4.1160 / 3.7294 26546 / 90999 0.3154 / 0.2922 4.1770 / 3.6905 26869 / 91999 0.3183 / 0.2925 4.2270 / 3.6389 27207 / 92999 0.3277 / 0.2929 4.2860 / 3.5906 27521 / 93999 0.3204 / 0.2932 4.3370 / 3.5424 27851 / 94999 0.3233 / 0.2935 4.3800 / 3.4833 28139 / 95999 0.2973 / 0.2936 4.4260 / 3.4385 28451 / 96999 0.3091 / 0.2937 4.4790 / 3.3924 28756 / 97999 0.3100 / 0.2939 4.5280 / 3.3451 29087 / 98999 0.3239 / 0.2942 4.5760 / 3.2901 29403 / 99999 0.3102 / 0.2943 4.6170 / 3.2342 29714 / 100999 0.3184 / 0.2946 4.6660 / 3.1855 30025 / 101999 0.3164 / 0.2948 4.7160 / 3.1375 30353 / 102999 0.3249 / 0.2951 4.7650 / 3.0843 30642 / 103999 0.3024 / 0.2952 4.8100 / 3.0387 30934 / 104999 0.3111 / 0.2953 4.8580 / 2.9942 31266 / 105999 0.3280 / 0.2956 4.9040 / 2.9384 31580 / 106999 0.3158 / 0.2958 4.9490 / 2.8867 31902 / 107999 0.3237 / 0.2961 4.9940 / 2.8331 32222 / 108999 0.3166 / 0.2963 5.0440 / 2.7830 32529 / 109999 0.3150 / 0.2964 5.0880 / 2.7327 32840 / 110999 0.3173 / 0.2966 5.1440 / 2.6879 33158 / 111999 0.3247 / 0.2969 5.1950 / 2.6387 33474 / 112999 0.3135 / 0.2970 5.2440 / 2.5889 33804 / 113999 0.3257 / 0.2973 5.2880 / 2.5336 34140 / 114999 0.3229 / 0.2975 5.3440 / 2.4826 34446 / 115999 0.3026 / 0.2975 5.3900 / 2.4338 34763 / 116999 0.3155 / 0.2977 5.4460 / 2.3870 35078 / 117999 0.3191 / 0.2979 5.4940 / 2.3371 35373 / 118999 0.3063 / 0.2979 5.5400 / 2.2908 35671 / 119999 0.3178 / 0.2981 5.5840 / 2.2431 35972 / 120999 0.3108 / 0.2982 5.6300 / 2.1955 36296 / 121999 0.3306 / 0.2985 5.6740 / 2.1423 36634 / 122999 0.3364 / 0.2988 5.7230 / 2.0880 36964 / 123999 0.3287 / 0.2990 5.7700 / 2.0349 37290 / 124999 0.3171 / 0.2992 5.8220 / 1.9844 37614 / 125999 0.3224 / 0.2994 5.8680 / 1.9323 37958 / 126999 0.3330 / 0.2996 5.9130 / 1.8759 38281 / 127999 0.3181 / 0.2998 5.9530 / 1.8224 38610 / 128999 0.3199 / 0.2999 5.9990 / 1.7697 38919 / 129999 0.3154 / 0.3000 6.0560 / 1.7243 39241 / 130999 0.3194 / 0.3002 6.1040 / 1.6736 39567 / 131999 0.3313 / 0.3004 6.1590 / 1.6240 39885 / 132999 0.3142 / 0.3005 6.2060 / 1.5739 40219 / 133999 0.3438 / 0.3008 6.2510 / 1.5202 40540 / 134999 0.3181 / 0.3010 6.3090 / 1.4722 40870 / 135999 0.3321 / 0.3012 6.3570 / 1.4201 41168 / 136999 0.3027 / 0.3012 6.4060 / 1.3743 41538 / 137999 0.3561 / 0.3016 6.4530 / 1.3146 41845 / 138999 0.3079 / 0.3017 6.4990 / 1.2666 42168 / 139999 0.3206 / 0.3018 6.5480 / 1.2162 42520 / 140999 0.3426 / 0.3021 6.6010 / 1.1612 42846 / 141999 0.3167 / 0.3022 6.6490 / 1.1102 43157 / 142999 0.3188 / 0.3023 6.7010 / 1.0625 43480 / 143999 0.3307 / 0.3025 6.7510 / 1.0123 43841 / 144999 0.3586 / 0.3029 6.8060 / 0.9561 44161 / 145999 0.3184 / 0.3030 6.8590 / 0.9069 44474 / 146999 0.3216 / 0.3031 6.9000 / 0.8573 44810 / 147999 0.3200 / 0.3032 6.9360 / 0.8033 45155 / 148999 0.3375 / 0.3035 6.9870 / 0.7497 45474 / 149999 0.3181 / 0.3036 7.0300 / 0.6997 45787 / 150999 0.3207 / 0.3037 7.0770 / 0.6512 46110 / 151999 0.3200 / 0.3038 7.1260 / 0.6012 46436 / 152999 0.3216 / 0.3039 7.1750 / 0.5507 46757 / 153999 0.3122 / 0.3040 7.2230 / 0.5010 47098 / 154999 0.3348 / 0.3042 7.2670 / 0.4478 47425 / 155999 0.3192 / 0.3042 7.3150 / 0.3972 47726 / 156999 0.3108 / 0.3043 7.3620 / 0.3508 48072 / 157999 0.3364 / 0.3045 7.4060 / 0.2970 48388 / 158999 0.3253 / 0.3046 7.4580 / 0.2485 48724 / 159999 0.3258 / 0.3048 7.5020 / 0.1965 49045 / 160999 0.3076 / 0.3048 7.5530 / 0.1471 49383 / 161999 0.3339 / 0.3050 7.5990 / 0.0949 49693 / 162999 0.3125 / 0.3050 7.6490 / 0.0473 50000 / 163943 0.3236 / 0.3051 7.6920 / 0.0000 ParaDRAM - NOTE: Computing the statistical properties of the Markov chain... ParaDRAM - NOTE: Computing the final decorrelated sample size... ParaDRAM - NOTE: Generating the output sample file: ParaDRAM - NOTE: .\out\mvn_serial_process_1_sample.txt ParaDRAM - NOTE: Computing the statistical properties of the final refined sample... ParaDRAM - NOTE: Mission Accomplished. ParaDRAM - NOTE: To read the generated output files, try the following: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.readReport() % to read the summary report from the output report file. ParaDRAM - NOTE: pmpd.readSample() % to read the final i.i.d. sample from the output sample file. ParaDRAM - NOTE: pmpd.readChain() % to read the uniquely-accepted points from the output chain file. ParaDRAM - NOTE: pmpd.readMarkovChain() % to read the Markov Chain. NOT recommended for extremely-large chains. ParaDRAM - NOTE: pmpd.readRestart() % to read the contents of an ASCII-format output restart file. ParaDRAM - NOTE: pmpd.readProgress() % to read the contents of an output progress file. ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte

The simulation output files

================================
Once the simulation is finished, the ParaDRAM routine generates 5 output files, each of which contains information about certain aspects of the simulation,
outPath = fullfile(currentDir,"out");
outFileList = ["mvn_serial_process_1_report.txt", "mvn_serial_process_1_progress.txt", "mvn_serial_process_1_sample.txt", "mvn_serial_process_1_chain.txt", "mvn_serial_process_1_restart.txt"];
for file = outFileList(end:-1:1)
open(fullfile(outPath,file))
end

General Guidelines on the ParaMonte visualization tools

===================================================================
In the following sections, we will discuss some of the postprocessing and visualization tools that automatically ship with the ParaMonte library. However, this is just the tip of the iceberg. For more detailed instructions on how to use the individual visualization tools of the ParaMonte library, visist the following MATLAB live scripts,
line / scatter plots
density plots
other plots

Reading the simulation output sample file

==================================================
You may have noticed that, upon finishing the simulation, the sampler provides hints on how to postprocess the results. We can read the generated output sample, contained in the file suffixed with *_sample.txt via,
pmpd.readSample();
ParaDRAM - WARNING: The ParaDRAM input simulation specification `pmpd.spec.outputDelimiter` is not set. ParaDRAM - WARNING: This information is essential for successful reading of the requested sample file(s). ParaDRAM - WARNING: Proceeding with the default assumption of comma-delimited sample file contents... ParaDRAM - NOTE: 1 files detected matching the pattern: "D:\temp\libparamonte_MATLAB\out\mvn_serial*" ParaDRAM - NOTE: processing file: "D:\temp\libparamonte_MATLAB\out\mvn_serial_process_1_sample.txt" ParaDRAM - NOTE: reading the file contents... ParaDRAM - NOTE: done in 1.645700 seconds. ParaDRAM - NOTE: ndim = 4, count = 11869 ParaDRAM - NOTE: computing the sample correlation matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 1.126500 seconds. ParaDRAM - NOTE: computing the sample covariance matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.219380 seconds. ParaDRAM - NOTE: computing the sample autocorrelation... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: done in 1.283500 seconds. ParaDRAM - NOTE: adding the graphics tools... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: creating the line3 plot object from scratch... ParaDRAM - NOTE: creating the scatter3 plot object from scratch... ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch... ParaDRAM - NOTE: creating the histogram plot object from scratch... ParaDRAM - NOTE: creating the histogram2 plot object from scratch... ParaDRAM - NOTE: creating the histfit plot object from scratch... ParaDRAM - NOTE: creating the contour plot object from scratch... ParaDRAM - NOTE: creating the contourf plot object from scratch... ParaDRAM - NOTE: creating the contour3 plot object from scratch... ParaDRAM - NOTE: creating the grid plot object from scratch... ParaDRAM - NOTE: The processed sample files are now stored in the newly-created component "pmpd.sampleList" of the ParaDRAM - NOTE: ParaDRAM object as a cell array. For example, to access the contents of the first (or the only) sample ParaDRAM - NOTE: file, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.sampleList{1}.df ParaDRAM - NOTE: ParaDRAM - NOTE: To access the plotting tools, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.sampleList{1}.plot.<PRESS TAB TO SEE THE LIST OF PLOTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For example, ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.sampleList{1}.plot.line.make() % to make 2D line plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.scatter.make() % to make 2D scatter plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.lineScatter.make() % to make 2D line-scatter plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.line3.make() % to make 3D line plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.scatter3.make() % to make 3D scatter plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.lineScatter3.make() % to make 3D line-scatter plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.contour3.make() % to make 3D kernel-density contour plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.contourf.make() % to make 2D kernel-density filled-contour plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.contour.make() % to make 2D kernel-density plots. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.histogram2.make() % to make 2D histograms. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.histogram.make() % to make 1D histograms. ParaDRAM - NOTE: pmpd.sampleList{1}.plot.grid.make() % to make GridPlot ParaDRAM - NOTE: ParaDRAM - NOTE: To plot or inspect the variable autocorrelations or the correlation/covariance matrices, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.sampleList{1}.stats.<PRESS TAB TO SEE THE LIST OF COMPONENTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte
sample = pmpd.sampleList{1};
By default, the ParaDRAM sampler generates a highly refined decorrelated final sample, essentially refining the Markov chain for as long as there is any residual auto-correlation in any of the individual variables of the Markov chain that has been sampled. This is an extremely strong and conservative refinement strategy (and it is so for a good reason: to ensure independent and identically distributed (i.i.d.) samples from the objective function). However, this extreme refinement can be mitigated and controlled by the setting the following ParaDRAM simulation specification variable prior to the simulation run,
pmpd.spec.sampleRefinementCount = 1; % refine each of the compact and verbose (Markov) chains only once.
pmpd.spec.helpme("sampleRefinementCount")
See this page: https://www.cdslab.org/paramonte/notes/usage/paradram/specifications/#samplerefinementcount
Rerun the simulation with this specification set and you will notice an increase in the final sample size by a few factors. Note that the above specification has to be set before running the simulation in the above, if you really want to use it.

Visualizing the refined sample

To quickly visualize the generated sample as a histogram, try,
sample.plot.histogram.make()
To plot all of the sampled variables, try,
sample.plot.histogram.reset(); % reset the plot properties to the default
ParaDRAM - NOTE: resetting the properties of the histogram plot...
sample.df.Properties.VariableNames
ans = 1×5 cell array {'SampleLogFunc'} {'SampleVariable1'} {'SampleVariable2'} {'SampleVariable3'} {'SampleVariable4'}
sample.plot.histogram
ans =
DensityPlot with properties: xcolumns: "SampleVariable1" dfref: [11869×5 table] figure: [1×1 struct] currentFig: [1×1 struct] rows: {} axes: [1×1 struct] histogram: [1×1 struct] target: [1×1 Target_class] legend: [1×1 struct]
sample.plot.histogram.xcolumns = sample.df.Properties.VariableNames(2:end);
sample.plot.histogram.legend.kws.location = "northeast";
sample.plot.histogram.legend.enabled = true;
sample.plot.histogram.make();
Or, we could also used the column indices of the variabels, instead of the column names,
sample.plot.histogram.reset() % reset the plot properties to the default
ParaDRAM - NOTE: resetting the properties of the histogram plot...
sample.plot.histogram.legend.enabled = false; % let's also remove the legends
sample.plot.histogram.make("xcolumns",[2,3,4,5]); % the column indices of all sampled variables in the output sample file.
We can also rebuild the histogram object from scratch, which includes reading again the input data used in the plotting from the dataFrame, by specifying the input parameter "hard",
sample.plot.histogram.reset("hard")
ParaDRAM - NOTE: creating the histogram plot object from scratch...
Make a trace plot of the refined sample,
sample.plot.line.reset();
ParaDRAM - NOTE: resetting the properties of the line plot...
sample.plot.line.colormap.values = "autumn"; % let's request autumn color-mapping
sample.plot.line.make();
We can simply take the axis logarithmic scale,
sample.plot.line.axes.kws.xscale = "log";
sample.plot.line.make();
or, we could do it via,
set(gca,'xscale','log','yscale','linear');
We can also include more variables in a single trace plot, for example, the the first, second, and the fourth dimensions of the objective function,
sample.plot.line.reset(); % reset all settings including colormap to the default values.
ParaDRAM - NOTE: resetting the properties of the line plot...
sample.plot.line.make("ycolumns",[2,3,5]); % these are the column numbers in the dataFrame sample.df
set(gca,'xscale','log','yscale','linear');
By default, the color of the line in the trace-plot will represent the value returned by getLogFunc() at the given sampled point.
sample.plot.line.ccolumns
ans = "SampleLogFunc"
To turn the color mapping off, you can instead try,
sample.plot.line.reset(); % reset to the default settings
ParaDRAM - NOTE: resetting the properties of the line plot...
sample.plot.line.colormap.enabled = false; % disable varying-color lines
sample.plot.line.legend.enabled = true; % display legend
sample.plot.line.legend.kws.location = "best";
sample.plot.line.make("ycolumns",3:5); % these are the column numbers in the dataFrame sample.df
set(gca,'xscale','log','yscale','linear');
We can also change the properties of the lines as we wish,
sample.plot.line.reset();
ParaDRAM - NOTE: resetting the properties of the line plot...
sample.plot.line.colormap.enabled = false;
sample.plot.line.legend.enabled = true; % display legend
sample.plot.line.plot.kws.linewidth = 1.0;
sample.plot.line.plot.kws.linestyle = ":"; % we can add any new proprty to plot_kws, that is an input argument to MATLAB's plot().
sample.plot.line.make("ycolumns",["SampleLogFunc",3:5]); % we can mix column names with column numbers
set(gca,'xscale','log','yscale','linear');
To make a scatter trace-plot of the sampled points, try,
sample.plot.scatter.make();
How about plotting other variables?
for variableName = sample.df.Properties.VariableNames(2:end)
sample.plot.scatter.make("xcolumns",variableName,"ycolumns","SampleLogFunc");
end
Let's turn off the scatter plot's varying colors and give it single color,
sample.plot.scatter.reset("hard");
ParaDRAM - NOTE: creating the scatter plot object from scratch...
sample.plot.scatter.colormap.enabled = false;
sample.plot.scatter.legend.enabled = true;
sample.plot.scatter.make("xcolumns",[3,4],"ycolumns","SampleLogFunc");
How about lineScatter plots?
sample.plot.lineScatter.colormap.values = "autumn";
sample.plot.lineScatter.make();
plot it on logarithmic axis,
sample.plot.lineScatter.make();
set(gca,'xscale','log','yscale','linear');
Let's turn off the varying-color,
sample.plot.reset("lineScatter")
ParaDRAM - NOTE: resetting the properties of the line plot... ParaDRAM - NOTE: resetting the properties of the scatter plot... ParaDRAM - NOTE: resetting the properties of the lineScatter plot... ParaDRAM - NOTE: resetting the properties of the line3 plot... ParaDRAM - NOTE: resetting the properties of the scatter3 plot... ParaDRAM - NOTE: resetting the properties of the lineScatter3 plot... ParaDRAM - NOTE: resetting the properties of the histogram plot... ParaDRAM - NOTE: resetting the properties of the histogram2 plot... ParaDRAM - NOTE: resetting the properties of the histfit plot... ParaDRAM - NOTE: resetting the properties of the contour plot... ParaDRAM - NOTE: resetting the properties of the contourf plot... ParaDRAM - NOTE: resetting the properties of the contour3 plot... ParaDRAM - NOTE: resetting the properties of the grid plot...
sample.plot.lineScatter.ccolumns
ans = "SampleLogFunc"
sample.plot.lineScatter.colormap
ans = struct with fields:
enabled: 1 values: "winter"
sample.plot.lineScatter.colormap.enabled = false; % no color mapping
sample.plot.lineScatter.make();
set(gca,'xscale','log','yscale','linear');
Make 3D line plots instead of 2D line plots,
sample.plot.line3.make();
Let's change the variables to plot,
disp(sample.plot.line3.xcolumns);
SampleVariable1
sample.plot.line3.xcolumns = {}; % reset the x-axis variable to the MCMC counts
sample.plot.line3.ycolumns = "SampleVariable1"; % ,"SampleVariable2"; % reset the y-axis variable to the first variable. We could include more variables.
sample.plot.line3.zcolumns % we will keep the z-axis as is: SampleLogFunc
ans = "SampleLogFunc"
sample.plot.line3.make();
Let's make a bivariate grid plot of the sample,
sample.plot.grid.reset("hard");
ParaDRAM - NOTE: creating the grid plot object from scratch...
sample.plot.grid.make();
generating subplot #1: (1,1) out of 25 generating subplot #2: (1,2) out of 25 generating subplot #3: (1,3) out of 25 generating subplot #4: (1,4) out of 25 generating subplot #5: (1,5) out of 25 generating subplot #6: (2,1) out of 25 generating subplot #7: (2,2) out of 25 generating subplot #8: (2,3) out of 25 generating subplot #9: (2,4) out of 25 generating subplot #10: (2,5) out of 25 generating subplot #11: (3,1) out of 25 generating subplot #12: (3,2) out of 25 generating subplot #13: (3,3) out of 25 generating subplot #14: (3,4) out of 25 generating subplot #15: (3,5) out of 25 generating subplot #16: (4,1) out of 25 generating subplot #17: (4,2) out of 25 generating subplot #18: (4,3) out of 25 generating subplot #19: (4,4) out of 25 generating subplot #20: (4,5) out of 25 generating subplot #21: (5,1) out of 25 generating subplot #22: (5,2) out of 25 generating subplot #23: (5,3) out of 25 generating subplot #24: (5,4) out of 25 generating subplot #25: (5,5) out of 25
In the above plot, the lower-triangle contour subplots represent the contours of the density of sampled points.
Let's change the diagonal plots of the grid plot to histfit,
sample.plot.grid.reset() % reset the grid plot to the original settings
ParaDRAM - NOTE: resetting the properties of the grid plot...
sample.plot.grid.plotType.diag.value = "histfit"; % by default, MATLAB's histfit(), fits a Gaussian to the histogram data.
sample.plot.grid.template.histfit.histfit % show the properties of the histfit plot
ans = struct with fields:
enabled: 0 kws: [1×1 struct] dist: "Normal" nbins: []
sample.plot.grid.colormap.values = "autumn";
sample.plot.grid.make("columns",[2,3,4,5]); % ignore the first column: SampleLogFunc as it is not Gaussian
generating subplot #1: (1,1) out of 16 generating subplot #2: (1,2) out of 16 generating subplot #3: (1,3) out of 16 generating subplot #4: (1,4) out of 16 generating subplot #5: (2,1) out of 16 generating subplot #6: (2,2) out of 16 generating subplot #7: (2,3) out of 16 generating subplot #8: (2,4) out of 16 generating subplot #9: (3,1) out of 16 generating subplot #10: (3,2) out of 16 generating subplot #11: (3,3) out of 16 generating subplot #12: (3,4) out of 16 generating subplot #13: (4,1) out of 16 generating subplot #14: (4,2) out of 16 generating subplot #15: (4,3) out of 16 generating subplot #16: (4,4) out of 16
for i = 1:4
sample.plot.grid.currentFig.subplotList{i,i}.currentFig.histfit{1}(1).FaceColor = "red"; % make histogram color red
end
We can also add targets to the subplots of the gridplot via addTarget() method of the grid object,
sample.plot.grid.reset("hard") % reset the plot settings to the original conditions
ParaDRAM - NOTE: creating the grid plot object from scratch...
sample.plot.grid.make("columns",1:5);
generating subplot #1: (1,1) out of 25 generating subplot #2: (1,2) out of 25 generating subplot #3: (1,3) out of 25 generating subplot #4: (1,4) out of 25 generating subplot #5: (1,5) out of 25 generating subplot #6: (2,1) out of 25 generating subplot #7: (2,2) out of 25 generating subplot #8: (2,3) out of 25 generating subplot #9: (2,4) out of 25 generating subplot #10: (2,5) out of 25 generating subplot #11: (3,1) out of 25 generating subplot #12: (3,2) out of 25 generating subplot #13: (3,3) out of 25 generating subplot #14: (3,4) out of 25 generating subplot #15: (3,5) out of 25 generating subplot #16: (4,1) out of 25 generating subplot #17: (4,2) out of 25 generating subplot #18: (4,3) out of 25 generating subplot #19: (4,4) out of 25 generating subplot #20: (4,5) out of 25 generating subplot #21: (5,1) out of 25 generating subplot #22: (5,2) out of 25 generating subplot #23: (5,3) out of 25 generating subplot #24: (5,4) out of 25 generating subplot #25: (5,5) out of 25
sample.plot.grid.addTarget() % no arguments will lead to adding the mode of SampleLogFunc as the target to the subplots
Instead of plotting the mode of SampleLogFunc, we can also plot other quantities of the variables, such as mean or median of the variables,
sample.plot.reset("grid") % reset the plot settings to the original conditions
ParaDRAM - NOTE: resetting the properties of the line plot... ParaDRAM - NOTE: resetting the properties of the scatter plot... ParaDRAM - NOTE: resetting the properties of the lineScatter plot... ParaDRAM - NOTE: resetting the properties of the line3 plot... ParaDRAM - NOTE: resetting the properties of the scatter3 plot... ParaDRAM - NOTE: resetting the properties of the lineScatter3 plot... ParaDRAM - NOTE: resetting the properties of the histogram plot... ParaDRAM - NOTE: resetting the properties of the histogram2 plot... ParaDRAM - NOTE: resetting the properties of the histfit plot... ParaDRAM - NOTE: resetting the properties of the contour plot... ParaDRAM - NOTE: resetting the properties of the contourf plot... ParaDRAM - NOTE: resetting the properties of the contour3 plot... ParaDRAM - NOTE: resetting the properties of the grid plot...
sample.plot.grid.make("columns",1:5);
generating subplot #1: (1,1) out of 25 generating subplot #2: (1,2) out of 25 generating subplot #3: (1,3) out of 25 generating subplot #4: (1,4) out of 25 generating subplot #5: (1,5) out of 25 generating subplot #6: (2,1) out of 25 generating subplot #7: (2,2) out of 25 generating subplot #8: (2,3) out of 25 generating subplot #9: (2,4) out of 25 generating subplot #10: (2,5) out of 25 generating subplot #11: (3,1) out of 25 generating subplot #12: (3,2) out of 25 generating subplot #13: (3,3) out of 25 generating subplot #14: (3,4) out of 25 generating subplot #15: (3,5) out of 25 generating subplot #16: (4,1) out of 25 generating subplot #17: (4,2) out of 25 generating subplot #18: (4,3) out of 25 generating subplot #19: (4,4) out of 25 generating subplot #20: (4,5) out of 25 generating subplot #21: (5,1) out of 25 generating subplot #22: (5,2) out of 25 generating subplot #23: (5,3) out of 25 generating subplot #24: (5,4) out of 25 generating subplot #25: (5,5) out of 25
sample.plot.grid.addTarget("mean") % target will represent the means of the variables on each subplot
Now, let's say that we want another set of target values to appear on these plots, for example, the last point that we have in our sample.
sample.plot.grid.reset() % reset the plot settings to the original conditions
ParaDRAM - NOTE: resetting the properties of the grid plot...
sample.plot.grid.make("columns",1:3);
generating subplot #1: (1,1) out of 9 generating subplot #2: (1,2) out of 9 generating subplot #3: (1,3) out of 9 generating subplot #4: (2,1) out of 9 generating subplot #5: (2,2) out of 9 generating subplot #6: (2,3) out of 9 generating subplot #7: (3,1) out of 9 generating subplot #8: (3,2) out of 9 generating subplot #9: (3,3) out of 9
sample.plot.grid.addTarget( [] ... skip the default target type
, "values", sample.df{end,:} ... specify the target values directly
) % add the last sample as the extra target, with all colors set to magenta
sample.plot.grid.hide("upper","colorbar") % hides the upper traingle and colorbar. Undo it by sample.plot.grid.show("upper","colorbar").
There are many more capabilities of the grid plot that were not discussed here. Take a look at the visualization links given in the above. You can also always get help via the helpme() (uncomment the following lines for help),
% sample.plot.grid.helpme("all")
% sample.plot.grid.helpme("plot")
% sample.plot.grid.helpme("addTarget")

Ensuring diminishing adaptation in the adaptive MCMC sampling

=============================================================================
By default, the ParaDRAM sampler adapts the proposal distribution of the sampler throughout the entire simulation. This breaks the ergodicity condition of the Markov chain, however, if the adaptation tends to continuously and progressively diminish throughout the entire simulation, then we potentially trust the results. The alternative is to simply turn off the adaptation, but that is rarely a good strategy. Instead, a better solution is to keep the adaptivity on, however, simulate for very long times to ensure convergence of the Markov chain to the target density.
A unique feature of the ParaDRAM sampler is that it provides a measure of the amount of adaptation of the proposal distribution of in the Markov chain by measuring how much the proposal distribution progressively changes throughout the simulation. In essence, the computed quantity, represented by the column AdaptationMeasure in the output chain file, represents an upper limit on the total variation distance of each updated proposal distribution with the last used proposal distribution before the update. We can plot the AdaptationMeasure for this sampling problem to ensure that the adaptation of the proposal distribution happens progressively less and less throughout the simulation. If AdaptationMeasure does not diminish monotonically throughout the simulation, it is a strong indication of the lack of convergence of the MCMC sampler to the target objective function.
The adaptation measure is reported in the simulation output *_chain.txt file, which is by default, in compact format, meaning that it only contains the set of uniquely-visited states of the full verbose Markov chain. This compact chain is however good enough to check for the convergence of the adaptive MCMC sampler.
pmpd.readChain(); % read the Markov chain in compact format.
ParaDRAM - WARNING: The ParaDRAM input simulation specification `pmpd.spec.outputDelimiter` is not set. ParaDRAM - WARNING: This information is essential for successful reading of the requested chain file(s). ParaDRAM - WARNING: Proceeding with the default assumption of comma-delimited chain file contents... ParaDRAM - NOTE: 1 files detected matching the pattern: ParaDRAM - NOTE: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial*" ParaDRAM - NOTE: processing file: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial_process_1_chain.txt" ParaDRAM - NOTE: reading the file contents... ParaDRAM - NOTE: done in 0.184050 seconds. ParaDRAM - NOTE: ndim = 4, count = 50000 ParaDRAM - NOTE: computing the sample correlation matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.072364 seconds. ParaDRAM - NOTE: computing the sample covariance matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.020549 seconds. ParaDRAM - NOTE: computing the sample autocorrelation... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: done in 0.233050 seconds. ParaDRAM - NOTE: adding the graphics tools... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: creating the line3 plot object from scratch... ParaDRAM - NOTE: creating the scatter3 plot object from scratch... ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch... ParaDRAM - NOTE: creating the histogram plot object from scratch... ParaDRAM - NOTE: creating the histogram2 plot object from scratch... ParaDRAM - NOTE: creating the histfit plot object from scratch... ParaDRAM - NOTE: creating the contour plot object from scratch... ParaDRAM - NOTE: creating the contourf plot object from scratch... ParaDRAM - NOTE: creating the contour3 plot object from scratch... ParaDRAM - NOTE: creating the grid plot object from scratch... ParaDRAM - NOTE: The processed chain files are now stored in the newly-created component "pmpd.chainList" of the ParaDRAM - NOTE: ParaDRAM object as a cell array. For example, to access the contents of the first (or the only) chain ParaDRAM - NOTE: file, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.chainList{1}.df ParaDRAM - NOTE: ParaDRAM - NOTE: To access the plotting tools, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.chainList{1}.plot.<PRESS TAB TO SEE THE LIST OF PLOTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For example, ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.chainList{1}.plot.line.make() % to make 2D line plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.scatter.make() % to make 2D scatter plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.lineScatter.make() % to make 2D line-scatter plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.line3.make() % to make 3D line plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.scatter3.make() % to make 3D scatter plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.lineScatter3.make() % to make 3D line-scatter plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.contour3.make() % to make 3D kernel-density contour plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.contourf.make() % to make 2D kernel-density filled-contour plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.contour.make() % to make 2D kernel-density plots. ParaDRAM - NOTE: pmpd.chainList{1}.plot.histogram2.make() % to make 2D histograms. ParaDRAM - NOTE: pmpd.chainList{1}.plot.histogram.make() % to make 1D histograms. ParaDRAM - NOTE: pmpd.chainList{1}.plot.grid.make() % to make GridPlot ParaDRAM - NOTE: ParaDRAM - NOTE: To plot or inspect the variable autocorrelations or the correlation/covariance matrices, try: ParaDRAM - NOTE: ParaDRAM - NOTE: pmpd.chainList{1}.stats.<PRESS TAB TO SEE THE LIST OF COMPONENTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte
chain = pmpd.chainList{1};
% first remove the zero values in Adaptation Measure for better
% illustration on log-log plot. set them all to nan, then convert nans to the last non-zero value in sequence.
chain.plot.line.dfref.AdaptationMeasure( chain.plot.line.dfref.AdaptationMeasure == 0 ) = nan;
chain.plot.line.dfref.AdaptationMeasure = fillmissing( chain.plot.line.dfref.AdaptationMeasure, "previous" );
chain.plot.line.surface.kws.lineWidth = 1;
chain.plot.line.make("ycolumns","AdaptationMeasure"); % plot adaptation measure
set(gca,"xscale","log","yscale","log");
A similar plot can be made with the full verbose Markov chain of visited states,
markov = pmpd.readMarkovChain(); % read the chain in full verbose (Markovian) format
ParaDRAM - WARNING: The ParaDRAM input simulation specification `pmpd.spec.outputDelimiter` is not set. ParaDRAM - WARNING: This information is essential for successful reading of the requested chain file(s). ParaDRAM - WARNING: Proceeding with the default assumption of comma-delimited chain file contents... ParaDRAM - NOTE: 1 files detected matching the pattern: ParaDRAM - NOTE: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial*" ParaDRAM - NOTE: processing file: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial_process_1_chain.txt" ParaDRAM - NOTE: reading the file contents... ParaDRAM - NOTE: done in 0.432100 seconds. ParaDRAM - NOTE: ndim = 4, count = 167575 ParaDRAM - NOTE: computing the sample correlation matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.065212 seconds. ParaDRAM - NOTE: computing the sample covariance matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.024562 seconds. ParaDRAM - NOTE: computing the sample autocorrelation... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: done in 0.535960 seconds. ParaDRAM - NOTE: adding the graphics tools... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: creating the line3 plot object from scratch... ParaDRAM - NOTE: creating the scatter3 plot object from scratch... ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch... ParaDRAM - NOTE: creating the histogram plot object from scratch... ParaDRAM - NOTE: creating the histogram2 plot object from scratch... ParaDRAM - NOTE: creating the histfit plot object from scratch... ParaDRAM - NOTE: creating the contour plot object from scratch... ParaDRAM - NOTE: creating the contourf plot object from scratch... ParaDRAM - NOTE: creating the contour3 plot object from scratch... ParaDRAM - NOTE: creating the grid plot object from scratch... ParaDRAM - NOTE: The processed chain files are now stored in the output variable as a cell array. For example, to ParaDRAM - NOTE: access the contents of the first (or the only) chain file stored in an output variable named ParaDRAM - NOTE: OUTPUT_CELL_ARRAY, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.df ParaDRAM - NOTE: ParaDRAM - NOTE: To access the plotting tools, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.<PRESS TAB TO SEE THE LIST OF PLOTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For example, ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line.make() % to make 2D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter.make() % to make 2D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter.make() % to make 2D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line3.make() % to make 3D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter3.make() % to make 3D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter3.make() % to make 3D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour3.make() % to make 3D kernel-density contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contourf.make() % to make 2D kernel-density filled-contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour.make() % to make 2D kernel-density plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram2.make() % to make 2D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram.make() % to make 1D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.grid.make() % to make GridPlot ParaDRAM - NOTE: ParaDRAM - NOTE: To plot or inspect the variable autocorrelations or the correlation/covariance matrices, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.stats.<PRESS TAB TO SEE THE LIST OF COMPONENTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte
markov = markov{1};
% first remove the zero values in Adaptation Measure for better
% illustration on log-log plot. set them all to nan, then convert nans to the last non-zero value in sequence.
markov.plot.line.dfref.AdaptationMeasure(markov.plot.line.dfref.AdaptationMeasure==0) = nan;
markov.plot.line.dfref.AdaptationMeasure = fillmissing(markov.plot.line.dfref.AdaptationMeasure,'previous');
markov.plot.line.ccolumns = []; % set the colormap to count of states;
markov.plot.line.surface.kws.lineWidth = 0.75;
markov.plot.line.make("ycolumns","AdaptationMeasure");
set(gca,'yscale','log');
markov.plot.line.reset(); % reset the line plot to the original settings
ParaDRAM - NOTE: resetting the properties of the line plot...
The above curves are exactly the kind of diminishing adaptation that we would want to see in a ParaDRAM simulation: At the beginning, the sampler appears to struggle for a while to find the shape of the target objective function, but around step 1000, it appears to have found the peak of the target and starts to quickly converge to and adapt the shape of the proposal distribution of the Markov Chain sampler to the shape of the objective function.

Reading the simulation output Markov chain file

=========================================================
There are many other types of plots that can be made. Take a look at the visualization links given in the above. Let's make those plots with the verbose (Markov) Chain instead.
markov = pmpd.readMarkovChain();
ParaDRAM - WARNING: The ParaDRAM input simulation specification `pmpd.spec.outputDelimiter` is not set. ParaDRAM - WARNING: This information is essential for successful reading of the requested chain file(s). ParaDRAM - WARNING: Proceeding with the default assumption of comma-delimited chain file contents... ParaDRAM - NOTE: 1 files detected matching the pattern: ParaDRAM - NOTE: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial*" ParaDRAM - NOTE: processing file: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial_process_1_chain.txt" ParaDRAM - NOTE: reading the file contents... ParaDRAM - NOTE: done in 0.578390 seconds. ParaDRAM - NOTE: ndim = 4, count = 167575 ParaDRAM - NOTE: computing the sample correlation matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.092776 seconds. ParaDRAM - NOTE: computing the sample covariance matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.031010 seconds. ParaDRAM - NOTE: computing the sample autocorrelation... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: done in 0.685230 seconds. ParaDRAM - NOTE: adding the graphics tools... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: creating the line3 plot object from scratch... ParaDRAM - NOTE: creating the scatter3 plot object from scratch... ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch... ParaDRAM - NOTE: creating the histogram plot object from scratch... ParaDRAM - NOTE: creating the histogram2 plot object from scratch... ParaDRAM - NOTE: creating the histfit plot object from scratch... ParaDRAM - NOTE: creating the contour plot object from scratch... ParaDRAM - NOTE: creating the contourf plot object from scratch... ParaDRAM - NOTE: creating the contour3 plot object from scratch... ParaDRAM - NOTE: creating the grid plot object from scratch... ParaDRAM - NOTE: The processed chain files are now stored in the output variable as a cell array. For example, to ParaDRAM - NOTE: access the contents of the first (or the only) chain file stored in an output variable named ParaDRAM - NOTE: OUTPUT_CELL_ARRAY, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.df ParaDRAM - NOTE: ParaDRAM - NOTE: To access the plotting tools, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.<PRESS TAB TO SEE THE LIST OF PLOTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For example, ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line.make() % to make 2D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter.make() % to make 2D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter.make() % to make 2D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line3.make() % to make 3D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter3.make() % to make 3D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter3.make() % to make 3D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour3.make() % to make 3D kernel-density contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contourf.make() % to make 2D kernel-density filled-contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour.make() % to make 2D kernel-density plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram2.make() % to make 2D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram.make() % to make 1D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.grid.make() % to make GridPlot ParaDRAM - NOTE: ParaDRAM - NOTE: To plot or inspect the variable autocorrelations or the correlation/covariance matrices, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.stats.<PRESS TAB TO SEE THE LIST OF COMPONENTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte
markov = markov{1};
Let's make a quick line trace plot of all variables in a single plot. But first, we have to keep in mind that the chain file contains more information than only the sampled points,
markov.df.Properties.VariableNames
ans = 1×11 cell array {'ProcessID'} {'DelayedRejectionStage'} {'MeanAcceptanceRate'} {'AdaptationMeasure'} {'BurninLocation'} {'SampleWeight'} {'SampleLogFunc'} {'SampleVariable1'} {'SampleVariable2'} {'SampleVariable3'} {'SampleVariable4'}
We will only plot the sample columns,
line = markov.plot.line;
line.legend.enabled = true;
line.axes.kws.xscale = "log";
line.colormap.enabled = false;
line.plot.kws.lineWidth = 0.75;
line.legend.kws.location = "northwest";
line.make("ycolumns", markov.df.Properties.VariableNames(8:end) ); % only plot these columns
Let's make a 3D lineScatter plot,
markov.plot.lineScatter3.reset("hard");
ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch...
markov.plot.lineScatter3.make("xcolumns", "SampleVariable1", "ycolumns", "SampleVariable3");
Because of the bad initial starting point of the sampling (which was intentional in this example), the colormap does not show the colors well close to the peak of the objective function. One solution is to cut the initial part of the sample that is causing the color-problem,
markov.count % the total number of sampled points
ans = 167575
markov.plot.lineScatter3.make("rows", 1000:10:markov.count ) % we will skip the first 10 points, and will include every other point (for the purpose of illustration).
Note that in the above and in all other cases in this example, you can always set the properties in two ways, either passing them as keyword arguments to the make() method, or directly setting them as object properties. For example,
markov.plot.lineScatter3.rows = 1000:10:markov.count;

Reading the simulation output Markov chain file in compact format

===============================================================================
The readMarkovChain() method of the ParaDRAM sampler often requires significant amounts of RAM and cen become slow if the simulated Markov chain has a very low acceptance rate or the number of dimensions of the objective funciton is large. In general, the most important output is the *_sample.txt file, and if anything else needs to be directly checked in the Markov chain (, like the adaptation meadure of the proposal to ensure convergence), it suffices to read and use the Markov chain in compact format. This can be done va the readChain() method of the ParaDRAM sampler which reads the contents of the output chain file in the (compact or verbose (Markov)) format that it has been written, without attempting to unroll the compact weighted chain. This compact chain is the sequence of uniquely-accepted points during the simulation.
chain = pmpd.readChain(); % read the Markov chain in the (compact) format as it is in the output file.
ParaDRAM - WARNING: The ParaDRAM input simulation specification `pmpd.spec.outputDelimiter` is not set. ParaDRAM - WARNING: This information is essential for successful reading of the requested chain file(s). ParaDRAM - WARNING: Proceeding with the default assumption of comma-delimited chain file contents... ParaDRAM - NOTE: 1 files detected matching the pattern: ParaDRAM - NOTE: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial*" ParaDRAM - NOTE: processing file: "D:\Dropbox\Projects\20180101_ParaMonte\paramontex\MATLAB\mlx\sampling_multivariate_normal_distribution_via_paradram\out\mvn_serial_process_1_chain.txt" ParaDRAM - NOTE: reading the file contents... ParaDRAM - NOTE: done in 0.300340 seconds. ParaDRAM - NOTE: ndim = 4, count = 50000 ParaDRAM - NOTE: computing the sample correlation matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.134110 seconds. ParaDRAM - NOTE: computing the sample covariance matrix... ParaDRAM - NOTE: creating the heatmap plot object from scratch... ParaDRAM - NOTE: done in 0.025217 seconds. ParaDRAM - NOTE: computing the sample autocorrelation... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: done in 0.278530 seconds. ParaDRAM - NOTE: adding the graphics tools... ParaDRAM - NOTE: creating the line plot object from scratch... ParaDRAM - NOTE: creating the scatter plot object from scratch... ParaDRAM - NOTE: creating the lineScatter plot object from scratch... ParaDRAM - NOTE: creating the line3 plot object from scratch... ParaDRAM - NOTE: creating the scatter3 plot object from scratch... ParaDRAM - NOTE: creating the lineScatter3 plot object from scratch... ParaDRAM - NOTE: creating the histogram plot object from scratch... ParaDRAM - NOTE: creating the histogram2 plot object from scratch... ParaDRAM - NOTE: creating the histfit plot object from scratch... ParaDRAM - NOTE: creating the contour plot object from scratch... ParaDRAM - NOTE: creating the contourf plot object from scratch... ParaDRAM - NOTE: creating the contour3 plot object from scratch... ParaDRAM - NOTE: creating the grid plot object from scratch... ParaDRAM - NOTE: The processed chain files are now stored in the output variable as a cell array. For example, to ParaDRAM - NOTE: access the contents of the first (or the only) chain file stored in an output variable named ParaDRAM - NOTE: OUTPUT_CELL_ARRAY, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.df ParaDRAM - NOTE: ParaDRAM - NOTE: To access the plotting tools, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.<PRESS TAB TO SEE THE LIST OF PLOTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For example, ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line.make() % to make 2D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter.make() % to make 2D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter.make() % to make 2D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.line3.make() % to make 3D line plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.scatter3.make() % to make 3D scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.lineScatter3.make() % to make 3D line-scatter plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour3.make() % to make 3D kernel-density contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contourf.make() % to make 2D kernel-density filled-contour plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.contour.make() % to make 2D kernel-density plots. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram2.make() % to make 2D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.histogram.make() % to make 1D histograms. ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.plot.grid.make() % to make GridPlot ParaDRAM - NOTE: ParaDRAM - NOTE: To plot or inspect the variable autocorrelations or the correlation/covariance matrices, try: ParaDRAM - NOTE: ParaDRAM - NOTE: OUTPUT_CELL_ARRAY{1}.stats.<PRESS TAB TO SEE THE LIST OF COMPONENTS> ParaDRAM - NOTE: ParaDRAM - NOTE: For more information and examples on the usage, visit: ParaDRAM - NOTE: ParaDRAM - NOTE: https://www.cdslab.org/paramonte
chain = chain{1};
The weight of each state in the compact chain is the number of times that proposed states from that state have been rejected before moving to the next accepted state in the chain,
chain.df
ans = 50000×11 table
 ProcessIDDelayedRejectionStageMeanAcceptanceRateAdaptationMeasureBurninLocationSampleWeightSampleLogFuncSampleVariable1SampleVariable2SampleVariable3SampleVariable4
1101.000000000000000011-3.292231700000000e+020000
2100.666666670000000012-3.193877600000000e+02-1.131467100000000-1.0114699000000002.440996800000000-1.638003900000000
3100.600000000000000022-2.928708600000000e+02-2.708689500000000-1.7464301000000002.8874169000000000.428330540000000
4100.571428570000000032-2.907333800000000e+02-4.913327800000000-2.9019973000000002.351047100000000-0.017945322000000
5100.625000000000000041-2.464273300000000e+02-6.311879700000000-1.0692085000000003.3947139000000000.746133930000000
6100.443390130000000056-2.210873300000000e+02-6.030484400000000-0.0797748080000005.5681717000000002.826638500000000
7100.480497460000000061-2.024019700000000e+02-4.3554141000000000.7212299900000007.2478948000000000.972193160000000
8100.5129663700000000.59390401000000071-1.938168700000000e+02-4.5113634000000001.1759642000000006.729581700000000-0.097407638000000
9100.541615400000000081-1.844784800000000e+02-5.2852795000000001.2274968000000006.1575672000000000.904767930000000
10100.567081220000000091-1.324347100000000e+02-8.6517902000000011.35976080000000010.3890470000000002.323448900000000
â‹®
The compact chain format is a very efficient representation of the Markov chain in terms of memory footprint and processing speed. In partictular, the kernel density estimate plots work very well with the compact format of the Markov chain.
We can also make contour plots of bivariate density-maps of data. But let's make these plots for the set of uniquely-sampled points in the chain file (that is, the Markov chain in compact format),
burnin = chain.df.BurninLocation(end);
chain.plot.contour.make( "rows", burnin:chain.count ); % exclude the burnin episoide from the plot
To make filled contour plots of the bivariate density-maps of data,
chain.plot.contourf.make( "rows", burnin:chain.count );
or, make 3D contour plots of the density of the bivariate density-maps of data,
chain.plot.contour3.make( "rows" , burnin:chain.count )

Inspecting the autocorrelation of the MCMC chain and the final sample

====================================================================================
By default, the ParaDRAM sampler will decorrelate and refine the MCMC chain to the highest levels possible until it cannot detect any residual autocorrelation in the refined sample. Now, let's first take a look at the the autocorrelation plots of the raw MCMC chain,
line = markov.stats.autocorr.plot.line;
line.legend.enabled = true;
line.make();
For comparison, here is the autocorrelation plot of the variables of the final refined decorrelated sample,
line = sample.stats.autocorr.plot.line;
line.legend.enabled = true;
line.make();
The above AutoCorrelation plot for the refined sample is reassuring since the sampled states do not appear to be correlated with each other at all. This is because the ParaDRAM routine, by default, applies as many rounds of refinements and decorrelation of the Markov chain as necessary to remove any residual correlations from the final output final sample.
By contrast, unlike the refined sample, the Markov chain is significantly correlated with itself along each dimension as illustrated in the above. The large amount of autocorrelation seen for SampleLogFunc is because of the fact that we started the MCMC sampling from a very bad low-probability location, which is also visible the grid plots in the above.

Inspecting the correlation / covariance matrices of the final sample or the MCMC chain

=======================================================================================================
To construct and visualize the correlation matrix of the sampled points in the chain, try,
sample.stats.cormat.plot.heatmap.make();
For comparison, here is the correlation matrix of the raw (verbose) Markov chain,
markov.stats.cormat.plot.heatmap.make()
The two plots are mostly consistent with each other, an indication that good convergence has occurred, and that the requested chainSize of uniquely-visited states was more than enough to obtain a reasonable number i.i.d. random sample from the objective function.
We can also visualize the covariance matrix of the sample,
sample.stats.covmat.plot.heatmap.make();
For comparison, here is the covariance matrix of the raw Markov chain,
markov.stats.covmat.plot.heatmap.make();

Running a single-chain ParaDRAM simulation in parallel on multiple processors

===============================================================================================
For large scale time-intensive problems, you can also run the ParaDRAM sampler in parallel on multiple processors, and with as many processors as needed. One unique feature of the parallelization capabilities of the ParaMonte library is that, you do not need to code anything in your objective function in parallel. You do not even need to have the MATLAB parallel Toolbox installed on your system. All parallelizations and communications are automatically handled by the ParaMonte library routines. For instructions on how to sample your expensive objective functions in parallel, visit the ParaMonte library's documentation website.

Final Remarks

=================
There are many more functionalities and features of the ParaMonte library that were neither explored nor mentioned in this example Jupyter notebook. You can explore them by checking the existing components of each attribute of the ParaDRAM sampler class and by visiting the ParaMonte library's documentation website.

{% comment %}
{% endcomment %}
web counter