From 1cdbf756a15dd6f2a113b1b0ef8b35b9cfaf22b3 Mon Sep 17 00:00:00 2001 From: mawilson1234 Date: Wed, 12 Apr 2023 10:05:22 -0400 Subject: [PATCH 1/4] parallelize sampling from multivariate normal distributions with multiple repetitions --- R/bridge_sampler.R | 18 ++++++++------ R/bridge_sampler_normal.R | 49 +++++++++++++++++++++++++++++++------ R/bridge_sampler_warp3.R | 51 ++++++++++++++++++++++++++++++--------- 3 files changed, 91 insertions(+), 27 deletions(-) diff --git a/R/bridge_sampler.R b/R/bridge_sampler.R index 0c1335f..f7e2c54 100644 --- a/R/bridge_sampler.R +++ b/R/bridge_sampler.R @@ -263,9 +263,10 @@ bridge_sampler.stanfit <- function(samples = NULL, stanfit_model = samples, param_types = rep("real", ncol(samples_4_fit)), transTypes = transTypes, repetitions = repetitions, cores = cores, - packages = "rstan", maxiter = maxiter, silent = silent, - verbose = verbose, - r0 = 0.5, tol1 = 1e-10, tol2 = 1e-4)) + packages = c("rstan", "mvtnorm"), + maxiter = maxiter, silent = silent, + verbose = verbose, r0 = 0.5, tol1 = 1e-10, + tol2 = 1e-4)) } else { bridge_output <- do.call(what = paste0(".bridge.sampler.", method), args = list(samples_4_fit = samples_4_fit, @@ -278,9 +279,10 @@ bridge_sampler.stanfit <- function(samples = NULL, stanfit_model = samples, transTypes = transTypes, repetitions = repetitions, varlist = "stanfit", envir = sys.frame(sys.nframe()), - cores = cores, packages = "rstan", maxiter = maxiter, - silent = silent, verbose = verbose, - r0 = 0.5, tol1 = 1e-10, tol2 = 1e-4)) + cores = cores, packages = c("rstan", "mvtnorm"), + maxiter = maxiter, silent = silent, + verbose = verbose, r0 = 0.5, tol1 = 1e-10, + tol2 = 1e-4)) } return(bridge_output) @@ -490,7 +492,7 @@ bridge_sampler.stanreg <- bridge_output <- bridge_sampler(samples = samples, log_posterior = .stan_log_posterior, data = list(stanfit = sf), lb = lb, ub = ub, repetitions = repetitions, method = method, cores = cores, - use_neff = use_neff, packages = "rstan", + use_neff = use_neff, packages = c("rstan", "mvtnorm"), maxiter = maxiter, silent = silent, verbose = verbose) } else { @@ -500,7 +502,7 @@ bridge_sampler.stanreg <- repetitions = repetitions, varlist = "stanfit", envir = sys.frame(sys.nframe()), method = method, cores = cores, use_neff = use_neff, - packages = "rstan", maxiter = maxiter, + packages = c("rstan", "mvtnorm"), maxiter = maxiter, silent = silent, verbose = verbose) } return(bridge_output) diff --git a/R/bridge_sampler_normal.R b/R/bridge_sampler_normal.R index 6cf550c..f96c9c7 100644 --- a/R/bridge_sampler_normal.R +++ b/R/bridge_sampler_normal.R @@ -38,17 +38,26 @@ # sample from multivariate normal distribution and evaluate for posterior samples and generated samples q12 <- dmvnorm(samples_4_iter, mean = m, sigma = V, log = TRUE) - gen_samples <- vector(mode = "list", length = repetitions) - q22 <- vector(mode = "list", length = repetitions) - for (i in seq_len(repetitions)) { - gen_samples[[i]] <- rmvnorm(n_post, mean = m, sigma = V) - colnames(gen_samples[[i]]) <- colnames(samples_4_iter) - q22[[i]] <- dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE) - } + # gen_samples <- vector(mode = "list", length = repetitions) + # q22 <- vector(mode = "list", length = repetitions) + # for (i in seq_len(repetitions)) { + # gen_samples[[i]] <- rmvnorm(n_post, mean = m, sigma = V) + # colnames(gen_samples[[i]]) <- colnames(samples_4_iter) + # q22[[i]] <- dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE) + # } # evaluate log of likelihood times prior for posterior samples and generated samples q21 <- vector(mode = "list", length = repetitions) if (cores == 1) { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + gen_samples <- vector(mode = "list", length = repetitions) + q22 <- vector(mode = "list", length = repetitions) + for (i in seq_len(repetitions)) { + gen_samples[[i]] <- rmvnorm(n_post, mean = m, sigma = V) + colnames(gen_samples[[i]]) <- colnames(samples_4_iter) + q22[[i]] <- dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE) + } + q11 <- apply(.invTransform2Real(samples_4_iter, lb, ub, param_types), 1, log_posterior, data = data, ...) + .logJacobian(samples_4_iter, transTypes, lb, ub) for (i in seq_len(repetitions)) { @@ -57,6 +66,17 @@ } } else if (cores > 1) { if ( .Platform$OS.type == "unix") { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = + function(x) rmvnorm(n_post, mean = m, sigma = V), + mc.preschedule = FALSE, + mc.cores = cores) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + q22 <- parallel::mclapply(seq_along(gen_samples), FUN = + function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE), + mc.preschedule = FALSE, + mc.cores = cores) + split1 <- .split_matrix(matrix=.invTransform2Real(samples_4_iter, lb, ub, param_types), cores=cores) q11 <- parallel::mclapply(split1, FUN = function(x) apply(x, 1, log_posterior, data = data, ...), @@ -73,8 +93,21 @@ } } else { cl <- parallel::makeCluster(cores, useXDR = FALSE) + on.exit(parallel::stopCluster(cl)) sapply(packages, function(x) parallel::clusterCall(cl = cl, "require", package = x, character.only = TRUE)) + + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + parallel::clusterExport(cl = cl, varlist = c("m", "V", "n_post"), envir = environment()) + gen_samples <- parallel::clusterApplyLB(cl = cl, + x = seq_len(repetitions), + fun = function(i) rmvnorm(n_post, mean = m, sigma = V)) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) + q22 <- parallel::clusterApplyLB(cl = cl, + x = seq_along(gen_samples), + fun = function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE)) + parallel::clusterExport(cl = cl, varlist = varlist, envir = envir) if ( ! is.null(rcppFile)) { @@ -91,7 +124,7 @@ q21[[i]] <- parallel::parRapply(cl = cl, x = .invTransform2Real(gen_samples[[i]], lb, ub, param_types), log_posterior, data = data, ...) + .logJacobian(gen_samples[[i]], transTypes, lb, ub) } - parallel::stopCluster(cl) + # parallel::stopCluster(cl) } } if(verbose) { diff --git a/R/bridge_sampler_warp3.R b/R/bridge_sampler_warp3.R index 559e51a..2dc4394 100644 --- a/R/bridge_sampler_warp3.R +++ b/R/bridge_sampler_warp3.R @@ -39,20 +39,28 @@ # sample from multivariate normal distribution and evaluate for posterior samples and generated samples q12 <- dmvnorm((samples_4_iter - matrix(m, nrow = n_post, ncol = length(m), byrow = TRUE)) %*% - t(solve(L)), sigma = diag(ncol(samples_4_fit)), log = TRUE) - q22 <- vector(mode = "list", length = repetitions) - gen_samples <- vector(mode = "list", length = repetitions) - for (i in seq_len(repetitions)) { - gen_samples[[i]] <- rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))) - colnames(gen_samples[[i]]) <- colnames(samples_4_iter) - q22[[i]] <- dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) - } + t(solve(L)), sigma = diag(ncol(samples_4_fit)), log = TRUE) + # q22 <- vector(mode = "list", length = repetitions) + # gen_samples <- vector(mode = "list", length = repetitions) + # for (i in seq_len(repetitions)) { + # gen_samples[[i]] <- rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))) + # colnames(gen_samples[[i]]) <- colnames(samples_4_iter) + # q22[[i]] <- dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) + # } e <- as.brob( exp(1) ) - + # evaluate log of likelihood times prior for posterior samples and generated samples q21 <- vector(mode = "list", length = repetitions) if (cores == 1) { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + q22 <- vector(mode = "list", length = repetitions) + gen_samples <- vector(mode = "list", length = repetitions) + for (i in seq_len(repetitions)) { + gen_samples[[i]] <- rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))) + colnames(gen_samples[[i]]) <- colnames(samples_4_iter) + q22[[i]] <- dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) + } q11 <- log(e^(apply(.invTransform2Real(samples_4_iter, lb, ub, param_types), 1, log_posterior, data = data,...) + .logJacobian(samples_4_iter, transTypes, lb, ub)) + @@ -72,6 +80,15 @@ } } else if (cores > 1) { if ( .Platform$OS.type == "unix") { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + q22 <- vector(mode = "list", length = repetitions) + gen_samples <- vector(mode = "list", length = repetitions) + gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = + function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + q22 <- parallel::mclapply(seq_along(gen_samples), FUN = + function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) + split1a <- .split_matrix(matrix=.invTransform2Real(samples_4_iter, lb, ub, param_types), cores=cores) split1b <- .split_matrix(matrix=.invTransform2Real( matrix(2*m, nrow = n_post, ncol = length(m), byrow = TRUE) - samples_4_iter, lb, ub, param_types @@ -113,9 +130,21 @@ tmp_mat2, transTypes, lb, ub))) } } else { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples cl <- parallel::makeCluster(cores, useXDR = FALSE) + on.exit(parallel::stopCluster(cl)) sapply(packages, function(x) parallel::clusterCall(cl = cl, "require", package = x, character.only = TRUE)) - + + parallel::clusterExport(cl = cl, varlist = c("n_post", "samples_4_fit"), envir = environment()) + gen_samples <- parallel::clusterApplyLB(cl = cl, + x = seq_len(repetitions), + fun = function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) + q22 <- parallel::clusterApplyLB(cl = cl, + x = seq_along(gen_samples), + fun = function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE)) + parallel::clusterExport(cl = cl, varlist = varlist, envir = envir) if ( ! is.null(rcppFile)) { @@ -152,7 +181,7 @@ .logJacobian(matrix(m, nrow = n_post, ncol = length(m), byrow = TRUE) + gen_samples[[i]] %*% t(L), transTypes, lb, ub))) } - parallel::stopCluster(cl) + # parallel::stopCluster(cl) } } From 6b1b3f45c82fd581ed68812597ca47cb395d64a4 Mon Sep 17 00:00:00 2001 From: mawilson1234 Date: Wed, 12 Apr 2023 10:07:34 -0400 Subject: [PATCH 2/4] remove unnecessary preallocation in warp3 after parallelization --- R/bridge_sampler_warp3.R | 2 -- 1 file changed, 2 deletions(-) diff --git a/R/bridge_sampler_warp3.R b/R/bridge_sampler_warp3.R index 2dc4394..7e6aa5c 100644 --- a/R/bridge_sampler_warp3.R +++ b/R/bridge_sampler_warp3.R @@ -81,8 +81,6 @@ } else if (cores > 1) { if ( .Platform$OS.type == "unix") { # sample from multivariate normal distribution and evaluate for posterior samples and generated samples - q22 <- vector(mode = "list", length = repetitions) - gen_samples <- vector(mode = "list", length = repetitions) gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) From c3b72eefcd704fc02d49ad7045a7e133a30e1608 Mon Sep 17 00:00:00 2001 From: mawilson1234 Date: Wed, 12 Apr 2023 10:09:53 -0400 Subject: [PATCH 3/4] fix missing args to parallel::mclapply --- R/bridge_sampler_warp3.R | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/R/bridge_sampler_warp3.R b/R/bridge_sampler_warp3.R index 7e6aa5c..24aada2 100644 --- a/R/bridge_sampler_warp3.R +++ b/R/bridge_sampler_warp3.R @@ -82,10 +82,14 @@ if ( .Platform$OS.type == "unix") { # sample from multivariate normal distribution and evaluate for posterior samples and generated samples gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = - function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) + function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))), + mc.preschedule = FALSE, + mc.cores = cores) sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) q22 <- parallel::mclapply(seq_along(gen_samples), FUN = - function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) + function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE), + mc.preschedule = FALSE, + mc.cores = cores) split1a <- .split_matrix(matrix=.invTransform2Real(samples_4_iter, lb, ub, param_types), cores=cores) split1b <- .split_matrix(matrix=.invTransform2Real( From ea94a97cf5c73af58317469682733d5695720bae Mon Sep 17 00:00:00 2001 From: mawilson1234 Date: Wed, 12 Apr 2023 10:17:02 -0400 Subject: [PATCH 4/4] only use parallelization if repetitions > 1 and cores > 1 --- R/bridge_sampler_normal.R | 66 +++++++++++++++++++------------------- R/bridge_sampler_warp3.R | 67 +++++++++++++++++++-------------------- 2 files changed, 64 insertions(+), 69 deletions(-) diff --git a/R/bridge_sampler_normal.R b/R/bridge_sampler_normal.R index f96c9c7..72873d8 100644 --- a/R/bridge_sampler_normal.R +++ b/R/bridge_sampler_normal.R @@ -38,18 +38,8 @@ # sample from multivariate normal distribution and evaluate for posterior samples and generated samples q12 <- dmvnorm(samples_4_iter, mean = m, sigma = V, log = TRUE) - # gen_samples <- vector(mode = "list", length = repetitions) - # q22 <- vector(mode = "list", length = repetitions) - # for (i in seq_len(repetitions)) { - # gen_samples[[i]] <- rmvnorm(n_post, mean = m, sigma = V) - # colnames(gen_samples[[i]]) <- colnames(samples_4_iter) - # q22[[i]] <- dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE) - # } - - # evaluate log of likelihood times prior for posterior samples and generated samples - q21 <- vector(mode = "list", length = repetitions) - if (cores == 1) { - # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + if (cores == 1 || repetitions == 1) { + # if cores > 1 and repetitions > 1 we do this below using parallelization gen_samples <- vector(mode = "list", length = repetitions) q22 <- vector(mode = "list", length = repetitions) for (i in seq_len(repetitions)) { @@ -57,7 +47,11 @@ colnames(gen_samples[[i]]) <- colnames(samples_4_iter) q22[[i]] <- dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE) } - + } + + # evaluate log of likelihood times prior for posterior samples and generated samples + q21 <- vector(mode = "list", length = repetitions) + if (cores == 1) { q11 <- apply(.invTransform2Real(samples_4_iter, lb, ub, param_types), 1, log_posterior, data = data, ...) + .logJacobian(samples_4_iter, transTypes, lb, ub) for (i in seq_len(repetitions)) { @@ -67,16 +61,18 @@ } else if (cores > 1) { if ( .Platform$OS.type == "unix") { # sample from multivariate normal distribution and evaluate for posterior samples and generated samples - gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = - function(x) rmvnorm(n_post, mean = m, sigma = V), - mc.preschedule = FALSE, - mc.cores = cores) - sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) - q22 <- parallel::mclapply(seq_along(gen_samples), FUN = - function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE), - mc.preschedule = FALSE, - mc.cores = cores) - + if (repetitions > 1) { + gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = + function(x) rmvnorm(n_post, mean = m, sigma = V), + mc.preschedule = FALSE, + mc.cores = cores) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + q22 <- parallel::mclapply(seq_along(gen_samples), FUN = + function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE), + mc.preschedule = FALSE, + mc.cores = cores) + } + split1 <- .split_matrix(matrix=.invTransform2Real(samples_4_iter, lb, ub, param_types), cores=cores) q11 <- parallel::mclapply(split1, FUN = function(x) apply(x, 1, log_posterior, data = data, ...), @@ -97,17 +93,19 @@ sapply(packages, function(x) parallel::clusterCall(cl = cl, "require", package = x, character.only = TRUE)) - # sample from multivariate normal distribution and evaluate for posterior samples and generated samples - parallel::clusterExport(cl = cl, varlist = c("m", "V", "n_post"), envir = environment()) - gen_samples <- parallel::clusterApplyLB(cl = cl, - x = seq_len(repetitions), - fun = function(i) rmvnorm(n_post, mean = m, sigma = V)) - sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) - parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) - q22 <- parallel::clusterApplyLB(cl = cl, - x = seq_along(gen_samples), - fun = function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE)) - + if (repetitions > 1) { + # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + parallel::clusterExport(cl = cl, varlist = c("m", "V", "n_post"), envir = environment()) + gen_samples <- parallel::clusterApplyLB(cl = cl, + x = seq_len(repetitions), + fun = function(i) rmvnorm(n_post, mean = m, sigma = V)) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) + q22 <- parallel::clusterApplyLB(cl = cl, + x = seq_along(gen_samples), + fun = function(i) dmvnorm(gen_samples[[i]], mean = m, sigma = V, log = TRUE)) + } + parallel::clusterExport(cl = cl, varlist = varlist, envir = envir) if ( ! is.null(rcppFile)) { diff --git a/R/bridge_sampler_warp3.R b/R/bridge_sampler_warp3.R index 24aada2..5b20587 100644 --- a/R/bridge_sampler_warp3.R +++ b/R/bridge_sampler_warp3.R @@ -40,20 +40,8 @@ # sample from multivariate normal distribution and evaluate for posterior samples and generated samples q12 <- dmvnorm((samples_4_iter - matrix(m, nrow = n_post, ncol = length(m), byrow = TRUE)) %*% t(solve(L)), sigma = diag(ncol(samples_4_fit)), log = TRUE) - # q22 <- vector(mode = "list", length = repetitions) - # gen_samples <- vector(mode = "list", length = repetitions) - # for (i in seq_len(repetitions)) { - # gen_samples[[i]] <- rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))) - # colnames(gen_samples[[i]]) <- colnames(samples_4_iter) - # q22[[i]] <- dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) - # } - - e <- as.brob( exp(1) ) - - # evaluate log of likelihood times prior for posterior samples and generated samples - q21 <- vector(mode = "list", length = repetitions) - if (cores == 1) { - # sample from multivariate normal distribution and evaluate for posterior samples and generated samples + if (cores == 1 || repetitions == 1) { + # if cores > 1 or repetitions > 1, we do this below using parallelization q22 <- vector(mode = "list", length = repetitions) gen_samples <- vector(mode = "list", length = repetitions) for (i in seq_len(repetitions)) { @@ -61,7 +49,13 @@ colnames(gen_samples[[i]]) <- colnames(samples_4_iter) q22[[i]] <- dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE) } + } + e <- as.brob( exp(1) ) + + # evaluate log of likelihood times prior for posterior samples and generated samples + q21 <- vector(mode = "list", length = repetitions) + if (cores == 1) { q11 <- log(e^(apply(.invTransform2Real(samples_4_iter, lb, ub, param_types), 1, log_posterior, data = data,...) + .logJacobian(samples_4_iter, transTypes, lb, ub)) + e^(apply(.invTransform2Real(matrix(2*m, nrow = n_post, ncol = length(m), byrow = TRUE) - @@ -80,17 +74,18 @@ } } else if (cores > 1) { if ( .Platform$OS.type == "unix") { - # sample from multivariate normal distribution and evaluate for posterior samples and generated samples - gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = - function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))), - mc.preschedule = FALSE, - mc.cores = cores) - sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) - q22 <- parallel::mclapply(seq_along(gen_samples), FUN = - function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE), - mc.preschedule = FALSE, - mc.cores = cores) - + if (repetitions > 1) { + gen_samples <- parallel::mclapply(seq_len(repetitions), FUN = + function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit))), + mc.preschedule = FALSE, + mc.cores = cores) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + q22 <- parallel::mclapply(seq_along(gen_samples), FUN = + function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE), + mc.preschedule = FALSE, + mc.cores = cores) + } + split1a <- .split_matrix(matrix=.invTransform2Real(samples_4_iter, lb, ub, param_types), cores=cores) split1b <- .split_matrix(matrix=.invTransform2Real( matrix(2*m, nrow = n_post, ncol = length(m), byrow = TRUE) - samples_4_iter, lb, ub, param_types @@ -137,16 +132,18 @@ on.exit(parallel::stopCluster(cl)) sapply(packages, function(x) parallel::clusterCall(cl = cl, "require", package = x, character.only = TRUE)) - parallel::clusterExport(cl = cl, varlist = c("n_post", "samples_4_fit"), envir = environment()) - gen_samples <- parallel::clusterApplyLB(cl = cl, - x = seq_len(repetitions), - fun = function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) - sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) - parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) - q22 <- parallel::clusterApplyLB(cl = cl, - x = seq_along(gen_samples), - fun = function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE)) - + if (repetitions > 1) { + parallel::clusterExport(cl = cl, varlist = c("n_post", "samples_4_fit"), envir = environment()) + gen_samples <- parallel::clusterApplyLB(cl = cl, + x = seq_len(repetitions), + fun = function(x) rmvnorm(n_post, sigma = diag(ncol(samples_4_fit)))) + sapply(seq_along(gen_samples), function(i) colnames(gen_samples[[i]]) <- colnames(samples_4_iter)) + parallel::clusterExport(cl = cl, varlist = "gen_samples", envir = environment()) + q22 <- parallel::clusterApplyLB(cl = cl, + x = seq_along(gen_samples), + fun = function(i) dmvnorm(gen_samples[[i]], sigma = diag(ncol(samples_4_fit)), log = TRUE)) + } + parallel::clusterExport(cl = cl, varlist = varlist, envir = envir) if ( ! is.null(rcppFile)) {