Make vcov and vcov_args arguments explicit

DESCRIPTION

@@ -1,7 +1,7 @@
 Type: Package
 Package: parameters
 Title: Processing of Model Parameters
-Version: 0.22.2.11
+Version: 0.22.2.12
 Authors@R:
     c(person(given = "Daniel",
              family = "Lüdecke",
@@ -194,6 +194,7 @@ Suggests:
     RcppEigen,
     rmarkdown,
     rms,
+    rstan,
     rstanarm,
     sandwich,
     see (>= 0.8.1),

NEWS.md

@@ -9,12 +9,9 @@
   way to test the direction of the effect, which formerly was already (and still
   is) possible with `pd = TRUE` in `model_parameters()`.
 
-* `p_function()` gets a `vcov` and `vcov_args` argument to compute robust
-  standard errors for the confidence curves.
-
-* Functions `p_significance()` and `equivalence_test()` now pass arguments
-  `vcov` and `vcov_args` to `p_value()` and `ci()`, hence, tests can be based
-  on robust standard errors.
+* `p_function()`, `p_significance()` and `equivalence_test()` get a `vcov` and
+  `vcov_args` argument, so that results can be based on robust standard errors
+  and confidence intervals.
 
 * `equivalence_test()` and `p_significance()` work with objects returned by
   `model_parameters()`.

R/equivalence_test.R

@@ -19,10 +19,7 @@
 #' See [`?ggeffects::test_predictions`](https://strengejacke.github.io/ggeffects/reference/test_predictions.html)
 #' for details.
 #' @param verbose Toggle warnings and messages.
-#' @param ... Arguments passed to or from other methods, e.g. `ci()`. Arguments
-#' like `vcov` or `vcov_args` can be used to compute confidence intervals or
-#' p-values using a specific variance-covariance matrix for the standard
-#' errors..
+#' @param ... Arguments passed to or from other methods.
 #' @inheritParams model_parameters.merMod
 #' @inheritParams p_value
 #'
@@ -244,10 +241,21 @@
                                 range = "default",
                                 ci = 0.95,
                                 rule = "classic",
+                                vcov = NULL,
+                                vcov_args = NULL,
                                 verbose = TRUE,
                                 ...) {
   rule <- match.arg(tolower(rule), choices = c("bayes", "classic", "cet"))
-  out <- .equivalence_test_frequentist(x, range, ci, rule, verbose, ...)
+  out <- .equivalence_test_frequentist(
+    x,
+    range = range,
+    ci = ci,
+    rule = rule,
+    vcov = vcov,
+    vcov_args = vcov_args,
+    verbose,
+    ...
+  )
 
   if (is.null(attr(out, "pretty_names", exact = TRUE))) {
     attr(out, "pretty_names") <- format_parameters(x)
@@ -308,6 +316,8 @@
                                     ci = 0.95,
                                     rule = "classic",
                                     effects = c("fixed", "random"),
+                                    vcov = NULL,
+                                    vcov_args = NULL,
                                     verbose = TRUE,
                                     ...) {
   # ==== argument matching ====
@@ -319,7 +329,16 @@
   # ==== equivalent testing for fixed or random effects ====
 
   if (effects == "fixed") {
-    out <- .equivalence_test_frequentist(x, range, ci, rule, verbose, ...)
+    out <- .equivalence_test_frequentist(
+      x,
+      range = range,
+      ci = ci,
+      rule = rule,
+      vcov = vcov,
+      vcov_args = vcov_args,
+      verbose,
+      ...
+    )
   } else {
     out <- .equivalence_test_frequentist_random(x, range, ci, rule, verbose, ...)
   }
@@ -496,6 +515,8 @@
                                           range = "default",
                                           ci = 0.95,
                                           rule = "classic",
+                                          vcov = NULL,
+                                          vcov_args = NULL,
                                           verbose = TRUE,
                                           ...) {
   # ==== define rope range ====
@@ -522,14 +543,14 @@
 
   # ==== requested confidence intervals ====
 
-  params <- conf_int <- .ci_generic(x, ci = ci, ...)
+  params <- conf_int <- .ci_generic(x, ci = ci, vcov = vcov, vcov_args = vcov_args, ...)
   conf_int <- as.data.frame(t(conf_int[, c("CI_low", "CI_high")]))
 
 
   # ==== the "narrower" intervals (1-2*alpha) for CET-rules. ====
 
   alpha <- 1 - ci
-  conf_int2 <- .ci_generic(x, ci = (ci - alpha), ...)
+  conf_int2 <- .ci_generic(x, ci = (ci - alpha), vcov = vcov, vcov_args = vcov_args, ...)
   conf_int2 <- as.data.frame(t(conf_int2[, c("CI_low", "CI_high")]))
 
 
@@ -562,7 +583,7 @@
 
   # ==== (adjusted) p-values for tests ====
 
-  out$p <- .add_p_to_equitest(x, ci, range, ...)
+  out$p <- .add_p_to_equitest(x, ci, range, vcov = vcov, vcov_args = vcov_args, ...)
 
   attr(out, "rope") <- range
   out
@@ -650,12 +671,12 @@
 
 #' @keywords internal
 .equivalence_test_numeric <- function(ci = 0.95,
                                       ci_wide,
                                       ci_narrow,
                                       range_rope,
                                       rule,
                                       dof = Inf,
                                       verbose) {
   final_ci <- NULL
 
   # ==== HDI+ROPE decision rule, by Kruschke ====
@@ -758,7 +779,7 @@
 # same range / limits as the confidence interval, thus indeed representing a
 # normally distributed confidence interval. We then calculate the probability
 # mass of this interval that is inside the ROPE.
 .rope_coverage <- function(ci = 0.95, range_rope, ci_range, dof = Inf) {
   out <- .generate_posterior_from_ci(ci, ci_range, dof = dof)
   # compare: ci_range and range(out)
   # The SGPV refers to the proportion of the confidence interval inside the
@@ -768,7 +789,7 @@
 }
 
 
 .generate_posterior_from_ci <- function(ci = 0.95, ci_range, dof = Inf, precision = 10000) {
   # this function creates an approximate normal distribution that covers the
   # CI-range, i.e. we "simulate" a posterior distribution from a frequentist CI
 
@@ -804,7 +825,7 @@
 }
 
 
-.add_p_to_equitest <- function(model, ci, range, ...) {
+.add_p_to_equitest <- function(model, ci, range, vcov = NULL, vcov_args = NULL, ...) {
   tryCatch(
     {
       params <- insight::get_parameters(model)
@@ -816,7 +837,7 @@
       params$mu <- params$Estimate * -1
 
       # se
-      se <- standard_error(model, ...)
+      se <- standard_error(model, vcov = vcov, vcov_args = vcov_args, ...)
 
       stats::pt((range[1] - params$mu) / se$SE, df = dof, lower.tail = TRUE) +
         stats::pt((range[2] - params$mu) / se$SE, df = dof, lower.tail = FALSE)

R/p_direction.R

@@ -87,9 +87,11 @@ p_direction.lm <- function(x,
                            ci = 0.95,
                            method = "direct",
                            null = 0,
+                           vcov = NULL,
+                           vcov_args = NULL,
                            ...) {
   # generate normal distribution based on CI range
-  result <- .posterior_ci(x, ci, ...)
+  result <- .posterior_ci(x, ci, vcov = vcov, vcov_args = vcov_args, ...)
 
   # copy
   out <- result$out

R/p_significance.R

@@ -28,9 +28,7 @@ bayestestR::p_significance
 #' @inheritParams bayestestR::p_significance
 #' @inheritParams model_parameters.default
 #' @param verbose Toggle warnings and messages.
-#' @param ... Arguments passed to other methods, e.g. `ci()`. Arguments like
-#' `vcov` or `vcov_args` can be used to compute confidence intervals using a
-#' specific variance-covariance matrix for the standard errors.
+#' @param ... Arguments passed to other methods.
 #'
 #' @seealso For more details, see [`bayestestR::p_significance()`]. See also
 #' [`equivalence_test()`], [`p_function()`] and [`bayestestR::p_direction()`]
@@ -143,9 +141,15 @@ bayestestR::p_significance
 #'   plot(result)
 #' }
 #' @export
-p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRUE, ...) {
+p_significance.lm <- function(x,
+                              threshold = "default",
+                              ci = 0.95,
+                              vcov = NULL,
+                              vcov_args = NULL,
+                              verbose = TRUE,
+                              ...) {
   # generate normal distribution based on CI range
-  result <- .posterior_ci(x, ci, ...)
+  result <- .posterior_ci(x, ci, vcov = vcov, vcov_args = vcov_args, ...)
 
   # copy
   out <- result$out
@@ -203,7 +207,7 @@ p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRU
 
 # helper ----------------------------------------------------------------------
 
-.posterior_ci <- function(x, ci, ...) {
+.posterior_ci <- function(x, ci, vcov = NULL, vcov_args = NULL, ...) {
   # first, we need CIs
   if (inherits(x, "parameters_model")) {
     # for model_parameters objects, directly extract CIs
@@ -223,7 +227,7 @@ p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRU
       dof <- Inf
     }
   } else {
-    out <- ci(x, ci = ci, ...)
+    out <- ci(x, ci = ci, vcov = vcov, vcov_args = vcov_args, ...)
     dof <- .safe(insight::get_df(x, type = "wald"), Inf)
   }
   # we now iterate all confidence intervals and create an approximate normal

tests/testthat/test-brms.R

@@ -5,10 +5,11 @@
 skip_if_offline()
 skip_if_not_installed("withr")
 skip_if_not_installed("brms")
+skip_if_not_installed("rstan")
 
 withr::with_options(
   list(parameters_warning_exponentiate = TRUE),
   {
     test_that("mp, footer exp", {
       m <- suppressWarnings(insight::download_model("brms_bernoulli_1"))
       out <- parameters::model_parameters(m, exponentiate = FALSE)