R/confint_indirect.R

In manymome: Mediation, Moderation and Moderated-Mediation After Model Fitting

#' @title Confidence Interval of
#' Indirect Effect or Conditional
#' Indirect Effect
#'
#' @description Return the
#' confidence interval of the indirect
#' effect or conditional indirect effect
#' stored in the output of
#' [indirect_effect()] or
#' [cond_indirect()].
#'
#' @details It extracts and returns the
#' stored confidence interval
#' if available.
#'
#' The type of confidence interval
#' depends on the call used to
#' compute the effect. This function
#' merely retrieves the stored estimates,
#' which could be generated by
#' nonparametric bootstrapping,
#' Monte Carlo simulation, or other
#' methods to be supported in
#' the future, and uses them to form the
#' percentile confidence interval.
#'
#' @param object The output of
#' [indirect_effect()] or
#' [cond_indirect()].
#'
#' @param parm Ignored because the
#' stored object always has only one
#' parameter.
#'
#' @param level The level of confidence,
#' default is .95, returning the 95%
#' confidence interval.
#'
#' @param boot_type If bootstrap
#' confidence interval is to be formed,
#' the type of bootstrap confidence
#' interval. The supported types
#' are `"perc"` (percentile bootstrap
#' confidence interval, the recommended
#' method) and `"bc"`
#' (bias-corrected, or BC, bootstrap
#' confidence interval). If not supplied,
#' the stored `boot_type` will be used.
#'
#' @param ...  Additional arguments.
#' Ignored by the function.
#'
#' @return A numeric vector of
#' two elements, the limits of
#' the confidence interval.
#'
#' @seealso [indirect_effect()] and
#' [cond_indirect()]
#'
#' @examples
#'
#' dat <- modmed_x1m3w4y1
#'
#' # Indirect Effect
#'
#' library(lavaan)
#' mod1 <-
#' "
#' m1 ~ x
#' m2 ~ m1
#' y  ~ m2 + x
#' "
#' fit <- sem(mod1, dat,
#'            meanstructure = TRUE, fixed.x = FALSE,
#'            se = "none", baseline = FALSE)
#' # R should be at least 2000 or 5000 in real research.
#' out1 <- indirect_effect(x = "x", y = "y",
#'                         m = c("m1", "m2"),
#'                         fit = fit,
#'                         boot_ci = TRUE, R = 45, seed = 54151,
#'                         parallel = FALSE,
#'                         progress = FALSE)
#' out1
#' confint(out1)
#'
#'
#' @export


confint.indirect <- function(object,
                             parm,
                             level = .95,
                             boot_type,
                             ...) {
    if (missing(boot_type)) {
        ci_boot_type <- object$boot_type
      } else {
        ci_boot_type <- boot_type
      }
    has_ci <- FALSE
    if (isTRUE(!is.null(object$boot_ci))) {
        has_ci <- TRUE
        old_ci <- object$boot_ci
        ci_type <- "boot"
        ind_i <- object$boot_indirect
        if ((level == object$level) &&
            (ci_boot_type == object$boot_type)) {
            new_ci <- FALSE
          } else {
            new_ci <- TRUE
          }
      }
    if (isTRUE(!is.null(object$mc_ci))) {
        has_ci <- TRUE
        old_ci <- object$mc_ci
        ci_type <- "mc"
        ind_i <- object$mc_indirect
        if (level == object$level) {
            new_ci <- FALSE
          } else {
            new_ci <- TRUE
          }
      }
    if (has_ci) {
        if (new_ci) {
            out0 <- boot_ci_internal(t0 = object$indirect,
                            t = ind_i,
                            level = level,
                            boot_type = ifelse(ci_type == "boot",
                                               ci_boot_type,
                                               "perc"),
                            add_names = FALSE)
          } else {
            out0 <- old_ci
          }
      } else {
        warning("Confidence interval not in the object.")
        out0 <- c(NA, NA)
      }
    # Borrowed from stats::confint()
    probs <- c((1 - level) / 2, 1 - (1 - level) / 2)
    cnames <- paste(format(100 * probs,
                           trim = TRUE,
                           scientific = FALSE,
                           digits = 2), "%")
    if (has_ci) {
        if (ci_type == "boot") {
            tmp <- switch(ci_boot_type,
                          perc = "Percentile: ",
                          bc = "Bias-Corrected: ")
            cnames <- paste0(tmp, cnames)
          }
        if (ci_type == "mc") {
            cnames <- paste0("Monte Carlo: ", cnames)
          }
      }
    rnames <- paste0(object$y, "~", object$x)
    out <- array(data = out0,
                 dim = c(1, 2),
                 dimnames = list(rnames, cnames))
    out
  }