Add latest changes from C++ library

Author: apulsipher

Diff for: inst/include/epiworld/math/lfmcmc/lfmcmc-meat-print.hpp

@@ -5,15 +5,18 @@ template<typename TData>
 inline void LFMCMC<TData>::print(size_t burnin) const
 {
 
+    // For each statistic or parameter in the model, we print three values: 
+    // - mean, the 2.5% quantile, and the 97.5% quantile
     std::vector< epiworld_double > summ_params(n_parameters * 3, 0.0);
     std::vector< epiworld_double > summ_stats(n_statistics * 3, 0.0);
 
+    // Compute the number of samples to use based on burnin rate
     size_t n_samples_print = n_samples;
     if (burnin > 0)
     {
         if (burnin >= n_samples)
             throw std::length_error(
-                "The burnin is greater than the number of samples."
+                "The burnin is greater than or equal to the number of samples."
                 );
 
         n_samples_print = n_samples - burnin;
@@ -24,15 +27,16 @@ inline void LFMCMC<TData>::print(size_t burnin) const
         n_samples_print
         );
 
+    // Compute parameter summary values
     for (size_t k = 0u; k < n_parameters; ++k)
     {
 
         // Retrieving the relevant parameter
         std::vector< epiworld_double > par_i(n_samples_print);
         for (size_t i = burnin; i < n_samples; ++i)
         {
-            par_i[i] = accepted_params[i * n_parameters + k];
-            summ_params[k * 3] += par_i[i]/n_samples_dbl;
+            par_i[i-burnin] = accepted_params[i * n_parameters + k];
+            summ_params[k * 3] += par_i[i-burnin]/n_samples_dbl;
         }
 
         // Computing the 95% Credible interval
@@ -45,20 +49,20 @@ inline void LFMCMC<TData>::print(size_t burnin) const
 
     }
 
+    // Compute statistics summary values
     for (size_t k = 0u; k < n_statistics; ++k)
     {
 
         // Retrieving the relevant parameter
         std::vector< epiworld_double > stat_k(n_samples_print);
         for (size_t i = burnin; i < n_samples; ++i)
         {
-            stat_k[i] = accepted_stats[i * n_statistics + k];
-            summ_stats[k * 3] += stat_k[i]/n_samples_dbl;
+            stat_k[i-burnin] = accepted_stats[i * n_statistics + k];
+            summ_stats[k * 3] += stat_k[i-burnin]/n_samples_dbl;
         }
 
         // Computing the 95% Credible interval
         std::sort(stat_k.begin(), stat_k.end());
-
         summ_stats[k * 3 + 1u] = 
             stat_k[std::floor(.025 * n_samples_dbl)];
         summ_stats[k * 3 + 2u] =