From cb5cbcec98d05f9768297d958a141309dab64ae6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Tomasz=20Wo=C5=BAniak?= Date: Sun, 2 Jun 2024 22:36:47 +1000 Subject: [PATCH] corrected sampler for nu #6 --- R/specify_bvarpanel.R | 2 +- src/sample_mniw.cpp | 10 ++++------ 2 files changed, 5 insertions(+), 7 deletions(-) diff --git a/R/specify_bvarpanel.R b/R/specify_bvarpanel.R index a9d1a02..a6dded0 100644 --- a/R/specify_bvarpanel.R +++ b/R/specify_bvarpanel.R @@ -210,7 +210,7 @@ specify_starting_values_bvarPANEL = R6::R6Class( self$A = matrix(stats::rnorm(K * N, sd = 0.001), K, N) + diag(K)[,1:N] self$V = stats::rWishart(1, K + 1, diag(K))[,,1] self$Sigma = stats::rWishart(1, N + 1, diag(N))[,,1] - self$nu = 30 #stats::rgamma(1, shape = 1, scale = 30) + self$nu = N + 1 + 0.1 self$m = stats::rnorm(1, sd = 0.001) self$w = stats::rgamma(1, 1) self$s = stats::rgamma(1, 1) diff --git a/src/sample_mniw.cpp b/src/sample_mniw.cpp index c4deca3..b6c6663 100644 --- a/src/sample_mniw.cpp +++ b/src/sample_mniw.cpp @@ -136,6 +136,9 @@ double log_kernel_nu ( mat sum_aux_Sigma_c_inv(N, N); for (int c = 0; c < C; c++) { + double ldSc = log_det_sympd(aux_Sigma_c_cpp.slice(c)); + log_kernel_nu -= 0.5 * (aux_nu + N + K + 1) * ldSc; + sum_aux_Sigma_c_inv += aux_Sigma_c_inv.slice(c); } log_kernel_nu -= 0.5 * (aux_nu - N - 1) * trace(aux_Sigma * sum_aux_Sigma_c_inv); @@ -144,11 +147,6 @@ double log_kernel_nu ( log_kernel_nu -= C * R::lgammafn(0.5 * (aux_nu + 1 - n)); } // EDN n loop - for (int c = 0; c < C; c++) { - double ldSc = log_det_sympd(aux_Sigma_c_cpp.slice(c)); - log_kernel_nu -= 0.5 * (aux_nu + N + K + 1) * ldSc; - } // END c loop - return log_kernel_nu; } // END log_kernel_nu @@ -178,7 +176,7 @@ double sample_nu ( Cov_nu += R::psigamma( 0.5 * (aux_nu + 1 - n), 1); } // END n loop Cov_nu *= (C / 4); - Cov_nu -= (C * N * (N + 2)) * (2 * (aux_nu - N -1)); + Cov_nu -= (C * N * aux_nu) * (2 * pow(aux_nu - N - 1, 2)); Cov_nu = sqrt(0.01 / Cov_nu); // Metropolis-Hastings