Init dqrng

Author: MaartenMarsman

DESCRIPTION

@@ -20,9 +20,9 @@ Description: Bayesian variable selection methods for analyzing the structure of
 License: GPL (>= 2)
 URL: https://maartenmarsman.github.io/bgms/
 BugReports: https://github.com/MaartenMarsman/bgms/issues
-Imports: Rcpp (>= 1.0.7), RcppParallel, Rdpack, methods, coda
+Imports: Rcpp (>= 1.0.7), RcppParallel, Rdpack, methods, coda, dqrng
 RdMacros: Rdpack
-LinkingTo: Rcpp, RcppProgress, RcppArmadillo, RcppParallel
+LinkingTo: Rcpp, RcppProgress, RcppArmadillo, RcppParallel, dqrng
 RoxygenNote: 7.3.2
 Depends: 
     R (>= 2.10)

R/RcppExports.R

@@ -5,16 +5,8 @@ run_bgmCompare_parallel <- function(observations, num_groups, num_obs_categories
     .Call(`_bgms_run_bgmCompare_parallel`, observations, num_groups, num_obs_categories, sufficient_blume_capel, sufficient_pairwise, num_categories, main_alpha, main_beta, pairwise_scale, difference_scale, difference_selection_alpha, difference_selection_beta, difference_prior, iter, burnin, na_impute, missing_data_indices, is_ordinal_variable, baseline_category, difference_selection, main_effect_indices, pairwise_effect_indices, target_accept, nuts_max_depth, learn_mass_matrix, projection, group_membership, group_indices, interaction_index_matrix, inclusion_probability, num_chains, nThreads)
 }
 
-run_gibbs_sampler_for_bgmCompare <- function(chain_id, observations, num_groups, num_obs_categories, sufficient_blume_capel, sufficient_pairwise, num_categories, main_alpha, main_beta, pairwise_scale, difference_scale, difference_selection_alpha, difference_selection_beta, difference_prior, iter, burnin, na_impute, missing_data_indices, is_ordinal_variable, baseline_category, difference_selection, main_effect_indices, pairwise_effect_indices, target_accept, nuts_max_depth, learn_mass_matrix, projection, group_membership, group_indices, interaction_index_matrix, inclusion_probability) {
-    .Call(`_bgms_run_gibbs_sampler_for_bgmCompare`, chain_id, observations, num_groups, num_obs_categories, sufficient_blume_capel, sufficient_pairwise, num_categories, main_alpha, main_beta, pairwise_scale, difference_scale, difference_selection_alpha, difference_selection_beta, difference_prior, iter, burnin, na_impute, missing_data_indices, is_ordinal_variable, baseline_category, difference_selection, main_effect_indices, pairwise_effect_indices, target_accept, nuts_max_depth, learn_mass_matrix, projection, group_membership, group_indices, interaction_index_matrix, inclusion_probability)
-}
-
-run_bgm_parallel <- function(observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix, num_chains, nThreads) {
-    .Call(`_bgms_run_bgm_parallel`, observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix, num_chains, nThreads)
-}
-
-run_gibbs_sampler_for_bgm <- function(chain_id, observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix) {
-    .Call(`_bgms_run_gibbs_sampler_for_bgm`, chain_id, observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix)
+run_bgm_parallel <- function(observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix, num_chains, nThreads, seed) {
+    .Call(`_bgms_run_bgm_parallel`, observations, num_categories, interaction_scale, edge_prior, inclusion_probability, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, lambda, interaction_index_matrix, iter, burnin, num_obs_categories, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, is_ordinal_variable, reference_category, edge_selection, update_method, pairwise_effect_indices, target_accept, sufficient_pairwise, hmc_num_leapfrogs, nuts_max_depth, learn_mass_matrix, num_chains, nThreads, seed)
 }
 
 sample_omrf_gibbs <- function(no_states, no_variables, no_categories, interactions, thresholds, iter) {

R/bgm.R

@@ -336,7 +336,8 @@ bgm = function(
     nuts_max_depth = 10,
     learn_mass_matrix = FALSE,
     chains = 4,
-    cores = parallel::detectCores()
+    cores = parallel::detectCores(),
+    seed = NULL
 ) {
   # Check update method
   update_method_input = update_method
@@ -488,6 +489,18 @@ bgm = function(
     }
   }
 
+  #Setting the seed
+  if (missing(seed) || is.null(seed)) {
+    # Draw a random seed if none provided
+    seed <- sample.int(.Machine$integer.max, 1)
+  }
+
+  if (!is.numeric(seed) || length(seed) != 1 || is.na(seed) || seed < 0) {
+    stop("Argument 'seed' must be a single non-negative integer.")
+  }
+
+  seed <- as.integer(seed)
+
   out = run_bgm_parallel(
     observations = x, num_categories = num_categories,
     interaction_scale = interaction_scale, edge_prior = edge_prior,
@@ -507,7 +520,7 @@ bgm = function(
     target_accept = target_accept, sufficient_pairwise = sufficient_pairwise,
     hmc_num_leapfrogs = hmc_num_leapfrogs, nuts_max_depth = nuts_max_depth,
     learn_mass_matrix = learn_mass_matrix, num_chains = chains,
-    nThreads = cores
+    nThreads = cores, seed = seed
   )
 
   # Main output handler in the wrapper function

R/bgmCompare2.R

@@ -24,7 +24,8 @@ bgmCompare2 = function(
     nuts_max_depth = 10,
     learn_mass_matrix = FALSE,
     chains = 4,
-    cores = parallel::detectCores()
+    cores = parallel::detectCores(),
+    seed = NULL
 ) {
   # Check update method
   update_method_input = update_method
@@ -190,6 +191,15 @@ bgmCompare2 = function(
     projection = matrix(projection, ncol = 1) / sqrt(2)
   }
 
+  if (!is.null(seed)) {
+    if (!is.numeric(seed) || any(is.na(seed)) || any(seed < 0)) {
+      stop("Argument 'seed' must be a non-negative integer or vector of non-negative integers.")
+    }
+    # Force to integer type
+    seed <- as.integer(seed)
+    dqrng::dqset.seed(seed)
+  }
+
   # Call the Rcpp function
   out = run_bgmCompare_parallel(
     observations = observations,

src/Makevars

@@ -1,4 +1,4 @@
-CXX_STD = CXX11
+CXX_STD = CXX14
 
 ## Pull in the include-paths for RcppParallel
 PKG_CPPFLAGS = \

src/Makevars.win

@@ -1,4 +1,4 @@
-CXX_STD = CXX11
+CXX_STD = CXX14
 
 PKG_CPPFLAGS = \
   $(shell "$(R_HOME)\bin\Rscript.exe" -e "cat(RcppParallel::CxxFlags())")