R/test_pars.R

In power4mome: Power Analysis for Moderation and Mediation

#' @title Test All Free Parameters
#'
#' @description Test all free parameters,
#' including user-defined parameters,
#' for a `power4test` object.
#'
#' @details
#' This function is to be used in
#' [power4test()] for testing all
#' free and user-defined model
#' parameters, by
#' setting it to the `test_fun`
#' argument.
#'
#' For models fitted by `lavaan`,
#' it uses [lavaan::parameterEstimates()]
#' to do the test. If bootstrapping was
#' requested (by setting `se = "boot"`),
#' then it supports bootstrap
#' confidence intervals returned by
#' [lavaan::parameterEstimates()].
#'
#' It has preliminary, though limited,
#' supported for models fitted by
#' [stats::lm()] (through
#' [lmhelprs::many_lm()]). Tests are
#' conducted by ordinary least squares
#' confidence intervals based on
#' the *t* statistic, reported by
#' [stats::confint()] applied to
#' the output of [stats::lm()].
#'
#' # Find the names of parameters
#'
#' To use the argument `pars`, the
#' names as appeared in the function
#' [coef()] must be used. For the
#' output of `lavaan`, this can
#' usually be inferred from the
#' parameter syntax (e.g., `y~x`,
#' no space). If not sure, call
#' [coef()] on the output of `lavaan`.
#' If a parameter is labelled, then
#' the label should be used in `par`.
#'
#' If not sure, the function
#' [find_par_names()] can be used to
#' find valid names.
#'
#' @return
#' In its normal usage, it returns
#' the output returned by
#' [lavaan::parameterEstimates()]
#' or [lmhelprs::lm_list_to_partable()],
#' with the following modifications:
#'
#' - `est`: The parameter estimates,
#'  even if standardized estimates
#'  are requested (not `est.std`).
#'
#' - `cilo` and `cihi`: The
#'  lower and upper limits of the
#'  confidence interval (95% by
#'  default), respectively (not
#'  `ci.lower` and `ci.upper`).
#'
#' - `sig`: Whether a test by confidence
#'  interval is significant (`1`) or
#'  not significant (`0`).
#'
#' - `test_label`: A column of labels
#'  generated by
#'  [lavaan::lav_partable_labels()],
#'  which are usually the labels used by
#'  `coef()` to label the parameters.
#'
#' @inheritParams test_k_indirect_effects
#'
#' @inheritParams test_indirect_effect
#'
#' @param fit The fit object, to be
#' passed to [lavaan::parameterEstimates()],
#' [lavaan::standardizedSolution()],
#' or [lmhelprs::lm_list_to_partable()].
#'
#' @param standardized Logical. If `TRUE`,
#' [lavaan::standardizedSolution()] will
#' be used. Can be used only with models
#' fitted by `lavaan`.
#'
#' @param pars Optional. If set to
#' a character vector, only parameters
#' with `test_label` equal to values in
#' `pars` will be returned. See the
#' help page
#' on valid names.
#'
#' @param op Optional. If set to a
#' character vector, only parameters with
#' operators (e.g., `"~"`, `"=~"`) will
#' be returned. If both `pars` and `op`
#' are specified, only parameters meeting
#' *both* requirements will be returned.
#'
#' @param remove.nonfree Logical. If
#' `TRUE`, the default, only free
#' parameters will be returned. Ignored
#' if `standardized` is `TRUE` or
#' if the model is not fitted by
#' `lavaan`.
#'
#' @param ... Additional arguments to
#' be passed to [lavaan::parameterEstimates()],
#' [lavaan::standardizedSolution()],
#' or [lmhelprs::lm_list_to_partable()].
#'
#' @seealso [power4test()]
#'
#' @examples
#'
#' # Specify the model
#'
#' mod <-
#' "
#' m ~ x
#' y ~ m + x
#' "
#'
#' # Specify the population values
#'
#' mod_es <-
#' "
#' y ~ m: l
#' m ~ x: m
#' y ~ x: n
#' "
#'
#' # Simulate the data
#'
#' sim_only <- power4test(nrep = 2,
#'                        model = mod,
#'                        pop_es = mod_es,
#'                        n = 100,
#'                        do_the_test = FALSE,
#'                        iseed = 1234)
#'
#' # Do the tests in each replication
#'
#' test_out <- power4test(object = sim_only,
#'                        test_fun = test_parameters)
#'
#' print(test_out,
#'       test_long = TRUE)
#'
#' # Do the tests in each replication: Standardized solution
#' # Delta method SEs will be used to do the tests
#'
#' test_out <- power4test(object = sim_only,
#'                        test_fun = test_parameters,
#'                        test_args = list(standardized = TRUE))
#'
#' print(test_out,
#'       test_long = TRUE)
#'
#' # Do the tests in each replication: Parameters with the selected operator
#'
#' test_out <- power4test(object = sim_only,
#'                        test_fun = test_parameters,
#'                        test_args = list(op = "~"))
#'
#' print(test_out,
#'       test_long = TRUE)
#'
#' @export
test_parameters <- function(fit = fit,
                            standardized = FALSE,
                            pars = NULL,
                            op = NULL,
                            remove.nonfree = TRUE,
                            check_post_check = TRUE,
                            ...,
                            omnibus = c("no", "all_sig", "at_least_one_sig", "at_least_k_sig"),
                            at_least_k = 1,
                            fit_name = "fit",
                            get_map_names = FALSE,
                            get_test_name = FALSE) {
  omnibus <- match.arg(omnibus)
  map_names <- c(fit = fit_name)
  args <- list(...)
  if (get_map_names) {
    return(map_names)
  }
  if (get_test_name) {
    tmp <- character(0)
    if (!is.null(pars)) {
      tmp0 <- paste0("pars: ",
                     paste0(pars,
                            collapse = ","))
      tmp <- c(tmp,
               tmp0)
    }
    if (!is.null(op)) {
      tmp0 <- paste0("op: ",
                     paste0(op,
                            collapse = ","))
      tmp <- c(tmp,
               tmp0)
    }
    if (length(tmp) >= 1) {
      tmp <- paste0("(",
                    paste0(tmp,
                           collapse = "; "),
                    ")")
    } else {
      tmp <- character(0)
    }
    if (standardized) {
      return(paste("test_parameters: CIs (standardized)", tmp))
    } else {
      return(paste("test_parameters: CIs", tmp))
    }
  }

  if (inherits(fit, "lm_list")) {
    fit_type <- "lm_list"
  } else if (inherits(fit, "lavaan")) {
    fit_type <- "lavaan"
  } else {
    stop("fit is not a supported object.")
  }

  if (fit_type == "lm_list") {
    if (is.null(op) && is.null(pars)) {
      op <- "~"
    }
  }

  if (inherits(fit, "lavaan")) {
    fit_ok <- lavaan::lavInspect(fit, "converged") &&
              (suppressWarnings(lavaan::lavInspect(fit, "post.check") ||
               !check_post_check))
  } else {
    fit_ok <- TRUE
  }
  if (standardized) {
    if (fit_type != "lavaan") {
      stop('Standardized solution supported only for `lavaan` output.')
    }
    est <- lavaan::standardizedSolution(object = fit,
                                        pvalue = TRUE,
                                        ci = TRUE)
    if (!fit_ok) {
      est$est.std <- as.numeric(NA)
      est$ci.lower <- as.numeric(NA)
      est$ci.upper <- as.numeric(NA)
      est$pvalue <- as.numeric(NA)
      est$se <- as.numeric(NA)
    }
  } else {
    if (fit_type == "lm_list") {
      # TODO:
      # - Find a better way to handle level
      if (!is.null(args$level)) {
        ci_args <- list(level = args$level)
      } else {
        ci_args <- eval(formals(lmhelprs::lm_list_to_partable)$ci_args)
      }
      est <- lmhelprs::lm_list_to_partable(object = fit,
                                           ci = TRUE,
                                           ci_args = ci_args)
      est <- est[, c("lhs", "op", "rhs", "est", "se", "pvalue",
                     "ci.lower", "ci.upper")]
    } else {
      est <- lavaan::parameterEstimates(object = fit,
                                        pvalue = TRUE,
                                        ci = TRUE,
                                        remove.nonfree = remove.nonfree,
                                        ...)
      if (!fit_ok) {
        est$est.std <- as.numeric(NA)
        est$ci.lower <- as.numeric(NA)
        est$ci.upper <- as.numeric(NA)
        est$pvalue <- as.numeric(NA)
        est$se <- as.numeric(NA)
      }
    }
  }
  enames <- colnames(est)
  enames <- gsub("ci.lower",
                 "cilo",
                 x = enames,
                 fixed = TRUE)
  enames <- gsub("ci.upper",
                 "cihi",
                 x = enames,
                 fixed = TRUE)
  if (standardized) {
    enames <- gsub("est.std",
                   "est",
                   x = enames,
                   fixed = TRUE)
  }
  colnames(est) <- enames

  if (!fit_ok) {
    est$sig <- as.numeric(NA)
  } else {
    est$sig <- ifelse((est$cilo > 0) | (est$cihi < 0),
                      yes = 1,
                      no = 0)
  }
  test_label <- lavaan::lav_partable_labels(est)
  out <- cbind(test_label = test_label,
               est)
  if (!is.null(op)) {
    j <- which(out$op %in% op)
    if (!isTRUE(length(j) > 0)) {
      stop("'op' set but not found in the test results.")
    }
    out <- out[j, ]
  }
  if (!is.null(pars)) {
    j <- out$test_label %in% pars
    if (!is.null(out$label)) {
      j <- j | (out$label %in% pars)
    }
    j <- which(j)
    if (!isTRUE(length(j) > 0)) {
      stop("'pars' set but not found in the test results.")
    }
    out <- out[j, ]
  }
  if (omnibus == "no") {
    attr(out, "test_label") <- "test_label"
    class(out) <- class(est)
    return(out)
  } else {
    out2 <- out[1, ]
    out2[1, ] <- as.numeric(NA)
    tmp <- switch(omnibus,
                  all_sig = "All sig",
                  at_least_one_sig = "1+ sig",
                  at_least_k_sig = paste0(at_least_k,
                                          "+ sig)"))
    out2[1, "test_label"] <- tmp
    tmp <- switch(omnibus,
                  all_sig = as.numeric(isTRUE(all(out$sig == 1))),
                  at_least_one_sig = as.numeric(isTRUE(any(out$sig == 1))),
                  at_least_k_sig = as.numeric(isTRUE(sum(out$sig == 1) >= at_least_k)))
    out2$sig <- tmp
    if (any(is.na(out2$sig))) {
      out2$sig <- as.numeric(NA)
    }
    attr(out2, "test_label") <- "test_label"
    return(out2)
  }
}

#' @param object A `power4test` object.
#'
#' @param fit_name The name of the fit
#' results for which the parameter names
#' will be displayed. Default is `"fit"`.
#'
#' @examples
#'
#' # Finding valid parameter names
#'
#' find_par_names(sim_only)
#'
#' @rdname test_parameters
#' @export
find_par_names <- function(object,
                           fit_name = "fit") {
  if (!inherits(object, "power4test")) {
    stop("Only support a 'power4test' object.")
  }
  out <- tryCatch(methods::getMethod("coef",
                                     signature = "lavaan",
                                     where = asNamespace("lavaan"))(object$sim_all[[1]]$extra[[fit_name]]),
                  error = function(e) e)
  if (inherits(out, "error")) {
    stop("Error in getting the coefficients.")
  }
  names(out)
}