Anova

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: bgms
 Type: Package
 Title: Bayesian Analysis of Networks of Binary and/or Ordinal Variables
-Version: 0.1.4.2
+Version: 0.1.4.3
 Date: 2024-12-03
 Authors@R: c(
     person("Maarten", "Marsman", , "m.marsman@uva.nl", role = c("aut", "cre"),
@@ -20,7 +20,7 @@ URL: https://maartenmarsman.github.io/bgms/
 BugReports: https://github.com/MaartenMarsman/bgms/issues
 Imports: Rcpp (>= 1.0.7), Rdpack, methods
 RdMacros: Rdpack
-LinkingTo: Rcpp, RcppProgress
+LinkingTo: Rcpp, RcppProgress, RcppArmadillo
 RoxygenNote: 7.3.2
 Depends: 
     R (>= 2.10)

NEWS.md

@@ -1,8 +1,16 @@
+# bgms 0.1.4.3
+
+bgmCompare now allows for network comparison for two or more groups
+fixed a bug with selecting group differences of blume-capel parameters
+fixed a bug with handling the samples of blume-capel parameters
+fixed a bug with handling threshold estimation when missing categories and main_model = "Free".
+
 # bgms 0.1.4.2
 
 Fixed a bug with adjusting the variance of the proposal distributions.
 Fixed a bug with recoding data under the "collapse" condition.
 When selection = true, we run 2 * burnin iterations instead of 1 * burnin in the burnin phase. This helps ensure that the Markov chain used for estimating the pseudoposterior starts with good parameter values and that proposals are properly calibrated. In rare cases, the Markov chain could get stuck before. The default setting for the burnin is also changed from 1000 to 500.
+Changed the maximum standard deviation of the adaptive proposal from 20 back to 2.
 
 # bgms 0.1.4.1
 

R/RcppExports.R

@@ -13,7 +13,7 @@ gibbs_sampler <- function(observations, indicator, interactions, thresholds, no_
     .Call(`_bgms_gibbs_sampler`, observations, indicator, interactions, thresholds, no_categories, interaction_scale, proposal_sd, proposal_sd_blumecapel, edge_prior, theta, beta_bernoulli_alpha, beta_bernoulli_beta, dirichlet_alpha, Index, iter, burnin, n_cat_obs, sufficient_blume_capel, threshold_alpha, threshold_beta, na_impute, missing_index, variable_bool, reference_category, save, display_progress, edge_selection)
 }
 
-compare_gibbs_sampler <- function(observations_gr1, observations_gr2, no_categories_gr1, no_categories_gr2, interaction_scale, pairwise_difference_scale, main_difference_scale, pairwise_difference_prior, main_difference_prior, inclusion_probability_difference, pairwise_beta_bernoulli_alpha, pairwise_beta_bernoulli_beta, main_beta_bernoulli_alpha, main_beta_bernoulli_beta, Index, iter, burnin, n_cat_obs_gr1, n_cat_obs_gr2, sufficient_blume_capel_gr1, sufficient_blume_capel_gr2, threshold_alpha, threshold_beta, na_impute, missing_index_gr1, missing_index_gr2, ordinal_variable, reference_category, independent_thresholds, save = FALSE, display_progress = FALSE, difference_selection = TRUE) {
-    .Call(`_bgms_compare_gibbs_sampler`, observations_gr1, observations_gr2, no_categories_gr1, no_categories_gr2, interaction_scale, pairwise_difference_scale, main_difference_scale, pairwise_difference_prior, main_difference_prior, inclusion_probability_difference, pairwise_beta_bernoulli_alpha, pairwise_beta_bernoulli_beta, main_beta_bernoulli_alpha, main_beta_bernoulli_beta, Index, iter, burnin, n_cat_obs_gr1, n_cat_obs_gr2, sufficient_blume_capel_gr1, sufficient_blume_capel_gr2, threshold_alpha, threshold_beta, na_impute, missing_index_gr1, missing_index_gr2, ordinal_variable, reference_category, independent_thresholds, save, display_progress, difference_selection)
+compare_anova_gibbs_sampler <- function(observations, main_effect_indices, pairwise_effect_indices, projection, num_categories, num_groups, group_indices, interaction_scale, pairwise_difference_scale, main_difference_scale, pairwise_difference_prior, main_difference_prior, inclusion_probability_difference, pairwise_beta_bernoulli_alpha, pairwise_beta_bernoulli_beta, main_beta_bernoulli_alpha, main_beta_bernoulli_beta, Index, iter, burnin, num_obs_categories, sufficient_blume_capel, prior_threshold_alpha, prior_threshold_beta, na_impute, missing_data_indices, is_ordinal_variable, baseline_category, independent_thresholds, save_main = FALSE, save_pairwise = FALSE, save_indicator = FALSE, display_progress = FALSE, difference_selection = TRUE) {
+    .Call(`_bgms_compare_anova_gibbs_sampler`, observations, main_effect_indices, pairwise_effect_indices, projection, num_categories, num_groups, group_indices, interaction_scale, pairwise_difference_scale, main_difference_scale, pairwise_difference_prior, main_difference_prior, inclusion_probability_difference, pairwise_beta_bernoulli_alpha, pairwise_beta_bernoulli_beta, main_beta_bernoulli_alpha, main_beta_bernoulli_beta, Index, iter, burnin, num_obs_categories, sufficient_blume_capel, prior_threshold_alpha, prior_threshold_beta, na_impute, missing_data_indices, is_ordinal_variable, baseline_category, independent_thresholds, save_main, save_pairwise, save_indicator, display_progress, difference_selection)
 }
 

R/bgm.R

@@ -130,9 +130,9 @@
 #' \code{beta_bernoulli_beta = 1}.
 #' @param dirichlet_alpha The shape of the Dirichlet prior on the node-to-block
 #' allocation probabilities for the Stochastic Block model.
-#' @param na.action How do you want the function to handle missing data? If
-#' \code{na.action = "listwise"}, listwise deletion is used. If
-#' \code{na.action = "impute"}, missing data are imputed iteratively during the
+#' @param na_action How do you want the function to handle missing data? If
+#' \code{na_action = "listwise"}, listwise deletion is used. If
+#' \code{na_action = "impute"}, missing data are imputed iteratively during the
 #' MCMC procedure. Since imputation of missing data can have a negative impact
 #' on the convergence speed of the MCMC procedure, it is recommended to run the
 #' MCMC for more iterations. Also, since the numerical routines that search for
@@ -276,7 +276,11 @@
 #'  #Restore user par() settings
 #'  par(op)
 #' }
+#'
 #' @importFrom utils packageVersion
+#' @importFrom Rcpp evalCpp
+#' @importFrom Rdpack reprompt
+#'
 #' @export
 bgm = function(x,
                variable_type = "ordinal",
@@ -292,7 +296,7 @@ bgm = function(x,
                beta_bernoulli_alpha = 1,
                beta_bernoulli_beta = 1,
                dirichlet_alpha = 1,
-               na.action = c("listwise", "impute"),
+               na_action = c("listwise", "impute"),
                save = FALSE,
                display_progress = TRUE) {
 
@@ -342,12 +346,12 @@ bgm = function(x,
   if(burnin <= 0)
     stop("Parameter ``burnin'' needs to be a positive integer.")
 
-  #Check na.action -------------------------------------------------------------
-  na.action_input = na.action
-  na.action = try(match.arg(na.action), silent = TRUE)
-  if(inherits(na.action, what = "try-error"))
-    stop(paste0("The na.action argument should equal listwise or impute, not ",
-                na.action_input,
+  #Check na_action -------------------------------------------------------------
+  na_action_input = na_action
+  na_action = try(match.arg(na_action), silent = TRUE)
+  if(inherits(na_action, what = "try-error"))
+    stop(paste0("The na_action argument should equal listwise or impute, not ",
+                na_action_input,
                 "."))
   #Check save ------------------------------------------------------------------
   save_input = save
@@ -364,7 +368,7 @@ bgm = function(x,
 
   #Format the data input -------------------------------------------------------
   data = reformat_data(x = x,
-                       na.action = na.action,
+                       na_action = na_action,
                        variable_bool = variable_bool,
                        reference_category = reference_category)
   x = data$x
@@ -425,11 +429,33 @@ bgm = function(x,
     for(variable2 in (variable1 + 1):no_variables) {
       cntr =  cntr + 1
       Index[cntr, 1] = cntr
-      Index[cntr, 2] = variable1
-      Index[cntr, 3] = variable2
+      Index[cntr, 2] = variable1 - 1
+      Index[cntr, 3] = variable2 - 1
     }
   }
 
+  #Preparing the output --------------------------------------------------------
+  arguments = list(
+    no_variables = no_variables,
+    no_cases = nrow(x),
+    na_impute = na_impute,
+    variable_type = variable_type,
+    iter = iter,
+    burnin = burnin,
+    interaction_scale = interaction_scale,
+    threshold_alpha = threshold_alpha,
+    threshold_beta = threshold_beta,
+    edge_selection = edge_selection,
+    edge_prior = edge_prior,
+    inclusion_probability = theta,
+    beta_bernoulli_alpha = beta_bernoulli_alpha ,
+    beta_bernoulli_beta =  beta_bernoulli_beta,
+    dirichlet_alpha = dirichlet_alpha,
+    na_action = na_action,
+    save = save,
+    version = packageVersion("bgms")
+  )
+
   #The Metropolis within Gibbs sampler -----------------------------------------
   out = gibbs_sampler(observations = x,
                       indicator = indicator,
@@ -459,29 +485,6 @@ bgm = function(x,
                       display_progress = display_progress,
                       edge_selection = edge_selection)
 
-
-  #Preparing the output --------------------------------------------------------
-  arguments = list(
-    no_variables = no_variables,
-    no_cases = nrow(x),
-    na_impute = na_impute,
-    variable_type = variable_type,
-    iter = iter,
-    burnin = burnin,
-    interaction_scale = interaction_scale,
-    threshold_alpha = threshold_alpha,
-    threshold_beta = threshold_beta,
-    edge_selection = edge_selection,
-    edge_prior = edge_prior,
-    inclusion_probability = theta,
-    beta_bernoulli_alpha = beta_bernoulli_alpha ,
-    beta_bernoulli_beta =  beta_bernoulli_beta,
-    dirichlet_alpha = dirichlet_alpha,
-    na.action = na.action,
-    save = save,
-    version = packageVersion("bgms")
-  )
-
   if(save == FALSE) {
     if(edge_selection == TRUE) {
       indicator = out$indicator