A vignette about default prior standata and stancode (#176)

* add vignette, additional edits * add bmm version to first line of stancode --------- Co-authored-by: GidonFrischkorn <[email protected]>
venpopov · Mar 25, 2024 · 768a74d · 768a74d
1 parent ad78421
commit 768a74d
Show file tree

Hide file tree

Showing 15 changed files with 338 additions and 30 deletions.
diff --git a/R/bmm_model_IMM.R b/R/bmm_model_IMM.R
@@ -41,9 +41,9 @@
         fixed_parameters = list(mu1 = 0, mu2 = 0, kappa2 = -100),
         default_priors = list(
           mu1 = list(main = "student_t(1, 0, 1)"),
-          kappa = list(main = "normal(2,1)", effects = "normal(0,1)"),
-          a = list(main = "normal(0,1)", effects = "normal(0,1)"),
-          c = list(main = "normal(0,1)", effects = "normal(0,1)")
+          kappa = list(main = "normal(2, 1)", effects = "normal(0, 1)"),
+          a = list(main = "normal(0, 1)", effects = "normal(0, 1)"),
+          c = list(main = "normal(0, 1)", effects = "normal(0, 1)")
         ),
         void_mu = FALSE
       ),
@@ -96,9 +96,9 @@
         fixed_parameters = list(mu1 = 0, mu2 = 0, kappa2 = -100),
         default_priors = list(
           mu1 = list(main = "student_t(1, 0, 1)"),
-          kappa = list(main = "normal(2,1)", effects = "normal(0,1)"),
-          c = list(main = "normal(0,1)", effects = "normal(0,1)"),
-          s = list(main = "normal(0,1)", effects = "normal(0,1)")
+          kappa = list(main = "normal(2, 1)", effects = "normal(0, 1)"),
+          c = list(main = "normal(0, 1)", effects = "normal(0, 1)"),
+          s = list(main = "normal(0, 1)", effects = "normal(0, 1)")
         ),
         void_mu = FALSE
       ),
@@ -152,10 +152,10 @@
         fixed_parameters = list(mu1 = 0, mu2 = 0, kappa2 = -100),
         default_priors = list(
           mu1 = list(main = "student_t(1, 0, 1)"),
-          kappa = list(main = "normal(2,1)", effects = "normal(0,1)"),
-          a = list(main = "normal(0,1)", effects = "normal(0,1)"),
-          c = list(main = "normal(0,1)", effects = "normal(0,1)"),
-          s = list(main = "normal(0,1)", effects = "normal(0,1)")
+          kappa = list(main = "normal(2, 1)", effects = "normal(0, 1)"),
+          a = list(main = "normal(0, 1)", effects = "normal(0, 1)"),
+          c = list(main = "normal(0, 1)", effects = "normal(0, 1)"),
+          s = list(main = "normal(0, 1)", effects = "normal(0, 1)")
         ),
         void_mu = FALSE
       ),
@@ -174,11 +174,11 @@
 #' @name IMM
 #' @details `r model_info(.model_IMMfull(), components =c('domain', 'task', 'name', 'citation'))`
 #' #### Version: `IMMfull`
-#' `r model_info(.model_IMMfull(), components =c('requirements', 'parameters', 'fixed_parameters'))`
+#' `r model_info(.model_IMMfull(), components = c('requirements', 'parameters', 'fixed_parameters', 'links', 'prior'))`
 #' #### Version: `IMMbsc`
-#' `r model_info(.model_IMMbsc(), components =c('requirements', 'parameters', 'fixed_parameters'))`
+#' `r model_info(.model_IMMbsc(), components = c('requirements', 'parameters', 'fixed_parameters', 'links', 'prior'))`
 #' #### Version: `IMMabc`
-#' `r model_info(.model_IMMabc(), components =c('requirements', 'parameters', 'fixed_parameters'))`
+#' `r model_info(.model_IMMabc(), components =c('requirements', 'parameters', 'fixed_parameters', 'links', 'prior'))`
 #'
 #' Additionally, all IMM models have an internal parameter that is fixed to 0 to
 #' allow the model to be identifiable. This parameter is not estimated and is not
@@ -225,8 +225,8 @@
 #'
 #' # specify the full IMM model with explicit column names for non-target features and distances
 #' model1 <- IMMfull(resp_err = "dev_rad",
-#'                   nt_features = paste0('col_nt',1:7),
-#'                   nt_distances = paste0('dist_nt',1:7),
+#'                   nt_features = paste0('col_nt', 1:7),
+#'                   nt_distances = paste0('dist_nt', 1:7),
 #'                   setsize = 'set_size')
 #'
 #' # fit the model

diff --git a/R/bmm_model_mixture3p.R b/R/bmm_model_mixture3p.R
@@ -42,7 +42,7 @@
       fixed_parameters = list(mu1 = 0, mu2 = 0, kappa2 = -100),
       default_priors = list(
         mu1 = list(main = "student_t(1, 0, 1)"),
-        kappa = list(main = "normal(2,1)", effects = "normal(0,1)"),
+        kappa = list(main = "normal(2, 1)", effects = "normal(0, 1)"),
         thetat = list(main = "logistic(0, 1)"),
         thetant = list(main = "logistic(0, 1)")
       ),

diff --git a/R/bmm_model_sdmSimple.R b/R/bmm_model_sdmSimple.R
@@ -32,8 +32,8 @@
       fixed_parameters = list(mu = 0),
       default_priors = list(
         mu = list(main = "student_t(1, 0, 1)"),
-        kappa = list(main = "student_t(5,1.75,0.75)", effects = "normal(0,1)"),
-        c = list(main = "student_t(5,2,0.75)", effects = "normal(0,1)")
+        kappa = list(main = "student_t(5, 1.75, 0.75)", effects = "normal(0, 1)"),
+        c = list(main = "student_t(5, 2, 0.75)", effects = "normal(0, 1)")
       ),
       void_mu = FALSE
     ),

diff --git a/R/helpers-data.R b/R/helpers-data.R
@@ -245,6 +245,7 @@ rad2deg <- function(rad){
 standata.bmmformula <- function(object, data, model, prior = NULL, ...) {
   # check model, formula and data, and transform data if necessary
   formula <- object
+  configure_options(list(...))
   model <- check_model(model, data, formula)
   data <- check_data(model, data, formula)
   formula <- check_formula(model, data, formula)

diff --git a/R/helpers-model.R b/R/helpers-model.R
@@ -279,6 +279,9 @@ model_info.bmmmodel <- function(model, components = 'all') {
   links <- model$links
   links_info <- summarise_links(links)
 
+  priors <- model$default_priors
+  priors_info <- summarise_default_prior(priors)
+
   info_all <-   list(
     domain = paste0("* **Domain:** ", model$domain, "\n\n"),
     task = paste0("* **Task:** ", model$task, "\n\n"),
@@ -288,7 +291,8 @@ model_info.bmmmodel <- function(model, components = 'all') {
     requirements = paste0("* **Requirements:** \n\n  ", model$requirements, "\n\n"),
     parameters = paste0("* **Parameters:** \n\n  ", par_info, "\n"),
     fixed_parameters = paste0("* **Fixed parameters:** \n\n  ", fixed_par_info, "\n"),
-    links = paste0("* **Default parameter links:** \n\n  ", links_info, "\n")
+    links = paste0("* **Default parameter links:** \n\n     - ", links_info, "\n\n"),
+    prior = paste0("* **Default priors:** \n\n", priors_info, "\n")
   )
 
   if (length(components) == 1 && components == 'all') {
@@ -649,5 +653,19 @@ stancode.bmmformula <- function(object, data, model, prior = NULL, ...) {
   fit_args <- combine_args(nlist(config_args, dots, prior))
   fit_args$object <- fit_args$formula
   fit_args$formula <- NULL
-  brms::do_call(brms::stancode, fit_args)
+  code <- brms::do_call(brms::stancode, fit_args)
+  add_bmm_version_to_stancode(code)
+}
+
+
+add_bmm_version_to_stancode <- function(stancode) {
+  version <- packageVersion("bmm")
+  text <- paste0("and bmm ", version)
+  brms_comp <- regexpr("brms.*(?=\\n)", stancode, perl = T)
+  insert_loc <- brms_comp + attr(brms_comp, "match.length") - 1
+  new_stancode <- paste0(substr(stancode, 1, insert_loc),
+                         " ", text,
+                         substr(stancode, insert_loc + 1, nchar(stancode)))
+  class(new_stancode) <- class(stancode)
+  new_stancode
 }
diff --git a/R/helpers-postprocess.R b/R/helpers-postprocess.R
@@ -26,9 +26,15 @@ postprocess_brm.bmmmodel <- function(model, fit, ...) {
   dots <- list(...)
   class(fit) <- c('bmmfit','brmsfit')
   fit$version$bmm <- utils::packageVersion('bmm')
-  fit$bmm <- nlist(model, user_formula = dots$user_formula, configure_opts = dots$configure_opts)
+  fit$bmm <- nlist(model, user_formula = dots$user_formula,
+                   configure_opts = dots$configure_opts)
   attr(fit$data, 'data_name') <- attr(dots$fit_args$data, 'data_name')
+
+  # add bmm version to the stancode
+  fit$model <- add_bmm_version_to_stancode(fit$model)
+
   fit <- NextMethod('postprocess_brm')
+
   # clean up environments stored in the fit object
   reset_env(fit)
 }

diff --git a/R/helpers-prior.R b/R/helpers-prior.R
@@ -258,3 +258,18 @@ combine_prior <- function(prior1, prior2) {
   }
   return(prior)
 }
+
+
+summarise_default_prior <- function(prior_list) {
+  pars <- names(prior_list)
+  prior_info <- ""
+  for (par in pars) {
+    prior_info  <- paste0(prior_info, "   - `", par, "`:\n")
+    types <- names(prior_list[[par]])
+    for (type in types) {
+      prior <- prior_list[[par]][[type]]
+      prior_info  <- paste0(prior_info, "      - `", type, "`: ", prior, "\n")
+    }
+  }
+  prior_info
+}
diff --git a/R/utils.R b/R/utils.R
@@ -357,7 +357,7 @@ order_data_query <- function(model, data, formula) {
         when using brms postprocessing methods that rely on the data order, such as
         generating predictions. Assuming you assigned the result of fit_model to a
         variable called `fit`, you can extract the sorted data from the fitted object
-        with:\n\n   data_sorted <- fit$fit_args$data", width = 80), collapse = "\n")
+        with:\n\n   data_sorted <- fit$data", width = 80), collapse = "\n")
     caution_msg <- crayon::red(caution_msg)
     message(caution_msg)
   }

diff --git a/_pkgdown.yml b/_pkgdown.yml
@@ -79,3 +79,4 @@ articles:
       Guides to using various features of the **bmm** package
    contents:
     - bmm_bmmformula
+    - bmm_extract_info
diff --git a/man/IMM.Rd b/man/IMM.Rd
diff --git a/man/SDM.Rd b/man/SDM.Rd
diff --git a/man/mixture2p.Rd b/man/mixture2p.Rd
diff --git a/man/mixture3p.Rd b/man/mixture3p.Rd
diff --git a/vignettes/articles/bmm_vwm_crt.Rmd b/vignettes/articles/bmm_vwm_crt.Rmd
@@ -59,7 +59,7 @@ The `bmm` package implements several measurement models for analyzing continuous
 
 #### The Interference Measurement Model [@Oberauer_Lin_2017]  {.unnumbered}
 
-  - see `?IMM` and `vignette("imm")`
+  - see `?IMM` and `vignette("IMM")`
 
 #### The Signal Discrimination Model (SDM) by [@Oberauer_2023]  {.unnumbered}