R/print_sbt_std_boot.R

In semboottools: Bootstrapping Helpers for Structural Equation Modelling

#' @title Print a 'sbt_std_boot' Object
#'
#' @description Print method for a
#' 'sbt_std_boot' object, which
#' is the output of
#' [standardizedSolution_boot()].
#'
#' @details
#' The default format of the printout,
#' `"lavaan.printer"`,
#' is a compact version of the lavaan-style
#' printout, generated by `lavaan.printer`.
#' Alternatively, users can request a format
#' similar to that of the printout
#' of the summary of a `lavaan` output
#' by setting `output` to `"text"`. This
#' format can be used if `"lavaan.printer"`
#' failed.
#'
#' Users can also print the content just
#' as a data frame by setting `output`
#' to `"table"`. Not easy to read much
#' more compact.
#'
#' For the `"text"` or `"lavaan.printer"` format, users can
#' also select whether
#' only the standardized solution is
#' printed (the default) or whether
#' the standardized solution is appended
#' to the right of the printout.
#'
#' @param x Object of the class
#' `sbt_std_boot`, the output of
#' [standardizedSolution_boot()].
#'
#' @param ... Optional arguments to be
#' passed to [print()] methods.
#'
#' @param nd The number of digits
#' after the decimal place. Default
#' is 3.
#'
#' @param output String. How the results
#' are printed. If set to `"table"`,
#' the results are printed in a table
#' format similar to that of
#' [lavaan::standardizedSolution()].
#' If set to `"text"`, the results will be
#' printed in a text format similar to
#' the printout of the output of
#' [summary()] of
#' a 'lavaan'-class object. If set
#' to `"lavaan.printer"`, the default,
#' `lavaan.printer`
#' will be used to print a more compact
#' version of the `"text"` output.
#'
#' @param standardized_only Logical.
#' If `TRUE`, the default, only the
#' results for the standardized solution
#' will be printed. If `FALSE`, then
#' the standardized solution is printed
#' alongside the unstandardized solution,
#' as in the printout of the output
#' of [summary()] of a 'lavaan'-class
#' object.
#'
#' @param boot_ci_only Logical. Whether
#' only bootstrap confidence intervals
#' are printed. If `FALSE`, the default,
#' the delta method confidence intervals
#' by [lavaan::standardizedSolution()]
#' are also printed.
#'
#' @param drop_cols The name(s) of the
#' column(s) to drop
#' if output format is `"lavaan.printer"`.
#' Default is `"Z"`, to fit the print
#' out to the usual screen width of 80.
#'
#' @seealso [standardizedSolution_boot()]
#'
#' @examples
#' library(lavaan)
#' set.seed(5478374)
#' n <- 50
#' x <- runif(n) - .5
#' m <- .40 * x + rnorm(n, 0, sqrt(1 - .40))
#' y <- .30 * m + rnorm(n, 0, sqrt(1 - .30))
#' dat <- data.frame(x = x, y = y, m = m)
#' model <-
#' '
#' m ~ a*x
#' y ~ b*m
#' ab := a*b
#' '
#'
#' # Should set bootstrap to at least 2000 in real studies
#' fit <- sem(model, data = dat, fixed.x = FALSE,
#'            se = "boot",
#'            bootstrap = 50)
#' std_out <- standardizedSolution_boot(fit)
#' std_out
#' print(std_out, standardized_only = FALSE)
#'
#' @return
#'  `x` is returned invisibly. Called for its side effect.
#'
#' @author Shu Fai Cheung
#' <https://orcid.org/0000-0002-9871-9448>
#'
#' @export

print.sbt_std_boot <- function(x,
                               ...,
                               nd = 3,
                               output = c("lavaan.printer", "text", "table"),
                               standardized_only = TRUE,
                               boot_ci_only = FALSE,
                               drop_cols = "Z") {
    output <- match.arg(output)
    x_call <- attr(x, "call")
    if (output == "table") {
        NextMethod()
        return(invisible(x))
      }

    ptable <- attr(x, "partable")
    est0 <- attr(x, "est")
    est1 <- est0
    est1$id <- seq_len(nrow(est1))
    if (standardized_only) {
      i0 <- colnames(est1) %in% c("est", "se", "z", "pvalue",
                                  "ci.lower", "ci.upper")
      est1 <- est1[, !i0]
      est1 <- merge(est1,
                    x)
      i0 <- colnames(ptable) %in% c("est", "se",
                                    "user", "free",
                                    "ustart", "plabel",
                                    "start",
                                    "id")
      est1 <- merge(est1, ptable[, !i0])
      est1 <- est1[order(est1$id), ]
      est1$id <- NULL
      # if (boot_ci_only) {
      #   est1$ci.lower <- NULL
      #   est1$ci.upper <- NULL
      #   est1$se <- NULL
      #   est1$z <- NULL
      #   est1$pvalue <- NULL
      # }
    } else {
      # Unstandardized and Standardized
      # Always use boot CI for standardized
      i0 <- colnames(x) %in% c("se", "z", "pvalue",
                              "ci.lower", "ci.upper")
      est1 <- merge(est1,
                    x[, !i0])
      i0 <- colnames(ptable) %in% c("est", "se",
                                    "user", "free",
                                    "ustart", "plabel",
                                    "start",
                                    "id")
      est1 <- merge(est1, ptable[, !i0])
      est1 <- est1[order(est1$id), ]
      est1$id <- NULL
    }
    class(est1) <- class(est0)
    pe_attrib <- attr(x, "pe_attrib")
    tmp <- !(names(pe_attrib) %in% names(attributes(est1)))
    attributes(est1) <- c(attributes(est1),
                          pe_attrib[tmp])
    class(est1) <- c("lavaan.parameterEstimates", class(est1))

    # If both est and est.std are present:
    # - !standardized_only
    # If only est.std is present:
    # - if any of these are present: se, pvalue, ci.lower, ci.upper
    #   - both_ci
    # - else
    #   - boot_ci only
    if (output == "text") {
      if (!standardized_only) {
        # This is for diagnostic purpose
        # Therefore, not as informative as the
        # standardized_only printout
        tmp <- colnames(est1)
        tmp[tmp == "est.std"] <- "Standardized"
        tmp[tmp == "boot.ci.lower"] <- "ci.std.lower"
        tmp[tmp == "boot.ci.upper"] <- "ci.std.upper"
        tmp[tmp == "boot.se"] <- "Std.Err.std"
        tmp[tmp == "boot.p"] <- "pvalue.std"
        colnames(est1) <- tmp
        print(est1, ..., nd = nd)
        return(invisible(x))
      } else {
        level <- attr(x, "level")
        est2 <- est1
        tmp <- colnames(est2)
        tmp[tmp == "est.std"] <- "est"
        colnames(est2) <- tmp

        if (boot_ci_only) {
          est2$ci.lower <- NULL
          est2$ci.upper <- NULL
          est2$se <- est2$boot.se
          est2$z <- NULL
          est2$pvalue <- est2$boot.p
          est2$boot.se <- NULL
          est2$boot.p <- NULL
        }

        out <- utils::capture.output(print(est2, nd = nd))
        i <- grepl("Parameter Estimates:", out, fixed = TRUE)
        out[i] <- "Standardized Estimates Only"
        i <- grepl("  Standard errors  ", out, fixed = TRUE)
        # j <- unlist(gregexpr("Standard", out[i]))[1]
        j <- nchar(out[i])
        tmp <- "  Standard errors (boot.se)"
        tmp2 <- "Bootstrap"
        st0 <- paste0(tmp,
                      paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                            collapse = ""),
                      tmp2)
        tmp <- "  Confidence interval (boot.ci.)"
        tmp2 <- "Bootstrap"
        st1 <- paste0(tmp,
                      paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                            collapse = ""),
                      tmp2)
        j <- nchar(out[i])
        tmp <- "  Confidence Level (boot.ci.)"
        tmp2 <- paste0(formatC(level * 100, digits = 1, format = "f"),
                      "%")
        st2 <- paste0(tmp,
                      paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                            collapse = ""),
                      tmp2)
        tmp <- "  Bootstrap CI Type (boot.ci.)"
        tmp2 <- switch(attr(x, "boot_ci_type"),
                      perc = "Percentile",
                      bc = "Bias-Corrected",
                      bca.simple = "Bias-Corrected")
        st2b <- paste0(tmp,
                      paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                              collapse = ""),
                      tmp2)
        if (!is.null(est1$boot.p)) {
          tmp <- "  Bootstrap P-Value (boot.p)"
          tmp2 <- "Asymmetric P-Value"
          st2c <- paste0(tmp,
                        paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                                collapse = ""),
                        tmp2)
        } else {
          st2c <- NULL
        }
        tmp <- "  Standardization Type"
        tmp2 <- attr(x, "type")
        st3 <- paste0(tmp,
                      paste0(rep(" ", j - nchar(tmp) - nchar(tmp2)),
                            collapse = ""),
                      tmp2)
        out <- c(out[seq_len(which(i))],
                st0,
                st1,
                st2,
                st2b,
                st2c,
                st3,
                out[-seq_len(which(i))])
        out <- gsub("    Estimate  Std.Err",
                    "Standardized  Std.Err",
                    out)
        if (boot_ci_only) {
          i <- grepl("  Standard errors  ", out, fixed = TRUE)
          out <- out[!i]
          out <- gsub(" Std.Err ",
                      " boot.se ",
                      out)
          out <- gsub(" P(>|z|) ",
                      "  boot.p ",
                      out,
                      fixed = TRUE)
        }
        cat(out, sep = "\n")
        return(invisible(x))
      }
    } else {
      # lavaan.printer
      level <- attr(x, "level")
      boot_ci_type <- attr(x, "boot_ci_type")
      has_boot_p <- !is.null(est1$boot.p)
      boot_est_std_i <- attr(x, "boot_est_std")
      std_type <- attr(x, "type")
      # TODO:
      # - Revise store_boot() to store more information
      if (!is.null(boot_est_std_i)) {
        R <- sum(stats::complete.cases(boot_est_std_i))
      } else {
        R <- NA
      }
      if (!standardized_only) {
        # This is for diagnostic purpose
        # Therefore, not as informative as the
        # standardized_only printout
        est2 <- lavaan.printer::parameterEstimates_table_list(est1,
                          rename_cols = c("est" = "Est",
                                          "est.std" = "Std",
                                          "P(>|z|)" = "p",
                                          "S.E." = "SE",
                                          "boot.ci.lower" = "bCI.Lo",
                                          "boot.ci.upper" = "bCI.Up",
                                          "boot.se" = "bSE",
                                          "boot.p" = "bp"),
                          header_funs = list(hdr_std_bootstrap),
                          header_funs_args = list(list(level = level,
                                                       boot_ci_type = boot_ci_type,
                                                       has_boot_p = has_boot_p,
                                                       R = R,
                                                       std_type = std_type)),
                          footer_funs = list(ftr_ustd_std_bootstrap),
                          footer_funs_args = list(list(cnames = colnames(est1))))
        lavaan.printer::print_parameterEstimates_table_list(est2,
                                                            nd = nd,
                                                            drop = drop_cols)
        return(invisible(x))
      } else {

        est2 <- est1
        tmp <- colnames(est2)
        tmp[tmp == "est.std"] <- "est"
        colnames(est2) <- tmp

        if (boot_ci_only) {
          est2$ci.lower <- NULL
          est2$ci.upper <- NULL
          est2$se <- NULL
          est2$z <- NULL
          est2$pvalue <- NULL
        }
        est3 <- lavaan.printer::parameterEstimates_table_list(est2,
                          rename_cols = c("est" = "Std",
                                          "P(>|z|)" = "p",
                                          "S.E." = "SE",
                                          "boot.ci.lower" = "bCI.Lo",
                                          "boot.ci.upper" = "bCI.Up",
                                          "boot.se" = "bSE",
                                          "boot.p" = "bp"),
                          header_funs = list(hdr_std_bootstrap),
                          header_funs_args = list(list(level = level,
                                                       boot_ci_type = boot_ci_type,
                                                       has_boot_p = has_boot_p,
                                                       R = R,
                                                       std_type = std_type)),
                          footer_funs = list(ftr_std_bootstrap),
                          footer_funs_args = list(list(cnames = colnames(est2))))
        lavaan.printer::print_parameterEstimates_table_list(est3,
                                                            nd = nd,
                                                            drop = drop_cols)
        return(invisible(x))
      }
    }
  }


#' @noRd

hdr_std_bootstrap <- function(x,
                              level,
                              boot_ci_type,
                              has_boot_p,
                              R,
                              std_type) {
  out0 <- data.frame(Field = "Valid Bootstrap Samples:",
                      Value = as.character(R))
  out1 <- data.frame(Field = "Level of Confidence:",
                      Value = sprintf("%3.1f%%", level * 100))
  out2 <- data.frame(Field = "CI Type:",
                      Value = switch(boot_ci_type,
                                    perc = "Percentile",
                                    bc = "Bias-Corrected",
                                    bca.simple = "Bias-Corrected"))
  if (has_boot_p) {
    out3 <- data.frame(Field = "P-Value:",
                        Value = "Asymmetric")
  } else {
    out3 <- NULL
  }
  out4 <- data.frame(Field = "Standardization Type:",
                     Value = std_type)
  out <- rbind(out0,
                out1,
                out2,
                out3,
                out4)
  colnames(out) <- NULL
  attr(out, "section_title") <- "Bootstrapping:"
  out
}

ftr_ustd_std_bootstrap <- function(x,
                                   cnames) {
  out0 <- character(0)
  out0 <- c(out0,
            "- Estimate: Unstandardized estimates.")
  if ("se" %in% cnames) {
    out0 <- c(out0,
              "- SE: Standard errors of unstandardized estimates.")
  }
  if ("pvalue" %in% cnames) {
    out0 <- c(out0,
              "- p: P-values of unstandardized estimates.")
  }
  if ("ci.lower" %in% cnames) {
    out0 <- c(out0,
              "- CI.Lo, CI.Up: Confidence intervals of unstandardized estimates.")
  }
  out0 <- c(out0,
            "- Std: Standardized estimates.")
  if ("boot.se" %in% cnames) {
    out0 <- c(out0,
              "- bSE: Bootstrap standard errors of standardized estimates.")
  }
  if ("boot.ci.lower" %in% cnames) {
    out0 <- c(out0,
              "- bCI.Lo, bCI.Up: Bootstrap confidence intervals of standardized estimates.")
  }
  if ("boot.p" %in% cnames) {
    out0 <- c(out0,
              "- bp: Bootstrap p-values of standardized estimates.")
  }
  attr(out0, "section_title") <- "Footnote:"
  attr(out0, "print_fun") <- "cat"
  attr(out0, "strwrap_args") <- list(exdent = 2)
  out0
}

ftr_std_bootstrap <- function(x,
                              cnames) {
  out0 <- character(0)
  out0 <- c(out0,
            "- Std: Standardized estimates.")
  if ("se" %in% cnames) {
    out0 <- c(out0,
              "- SE: Delta method standard errors.")
  }
  if ("pvalue" %in% cnames) {
    out0 <- c(out0,
              "- p: Delta method p-values.")
  }
  if ("ci.lower" %in% cnames) {
    out0 <- c(out0,
              "- CI.Lo, CI.Up: Delta method confidence intervals.")
  }
  if ("boot.se" %in% cnames) {
    out0 <- c(out0,
               "- bSE: Bootstrap standard errors.")
  }
  if ("boot.ci.lower" %in% cnames) {
    out0 <- c(out0,
              "- bCI.Lo, bCI.Up: Bootstrap confidence intervals.")
  }
  if ("boot.p" %in% cnames) {
    out0 <- c(out0,
              "- bp: Bootstrap p-values.")
  }
  attr(out0, "section_title") <- "Footnote:"
  attr(out0, "print_fun") <- "cat"
  attr(out0, "strwrap_args") <- list(exdent = 2)
  out0
}