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",

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 @@ bayestestR::equivalence_test
 #' 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 @@ equivalence_test.lm <- function(x,
                                 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 @@ equivalence_test.merMod <- function(x,
                                     ci = 0.95,
                                     rule = "classic",
                                     effects = c("fixed", "random"),
+                                    vcov = NULL,
+                                    vcov_args = NULL,
                                     verbose = TRUE,
                                     ...) {
   # ==== argument matching ====
@@ -319,7 +329,16 @@ equivalence_test.merMod <- function(x,
   # ==== 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 @@ equivalence_test.ggeffects <- function(x,
                                           range = "default",
                                           ci = 0.95,
                                           rule = "classic",
+                                          vcov = NULL,
+                                          vcov_args = NULL,
                                           verbose = TRUE,
                                           ...) {
   # ==== define rope range ====
@@ -522,14 +543,14 @@ equivalence_test.ggeffects <- function(x,
 
   # ==== 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 @@ equivalence_test.ggeffects <- function(x,
 
   # ==== (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
@@ -804,7 +825,7 @@ equivalence_test.ggeffects <- function(x,
 }
 
 
-.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 @@ equivalence_test.ggeffects <- function(x,
       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-p_direction.R

@@ -20,6 +20,20 @@ test_that("p_direction", {
   expect_equal(x$pd, c(1, 0.5567, 0.9997, 0.9309, 1), tolerance = 1e-3)
 })
 
+test_that("p_direction", {
+  skip_if_not_installed("sandwich")
+  data(mtcars)
+  m <- lm(mpg ~ gear + wt + cyl + hp, data = mtcars)
+
+  set.seed(123)
+  x <- p_direction(m, ci = 0.8, vcov = "HC3")
+  expect_equal(x$pd, c(1, 0.6162, 0.9984, 0.8323, 0.8962), tolerance = 1e-3)
+
+  set.seed(123)
+  x <- p_direction(m, null = 0.2, vcov = "HC3")
+  expect_equal(x$pd, c(1, 0.5464, 0.9989, 0.88, 1), tolerance = 1e-3)
+})
+
 test_that("p_direction, glmmTMB", {
   skip_if_not_installed("glmmTMB")
   data(Salamanders, package = "glmmTMB")