R/CalcPairwiseWilcox.R

In bmcnellis/RSFIA: FIA Tree Mortality Analysis

#' @title Calculates pairwise wilcox for mort_df variables
#' @description TODO: this
#'
#' @export
CalcPairwiseWilcox <- function(xvar, yvar, mort_df = NULL, out_dir = getwd(),
                               alpha = c(0.05, 0.01, 0.001),
                               write = F, tpose = F, method = 'num_pval',
                               dunn = F, raw = F) {
  # Prep
  require(stats)
  od <- getwd()
  on.exit(setwd(od))
  if (write) setwd(out_dir)
  if (is.null(mort_df)) {
    mdf <- FIA_mortality_with_explanatory
  }
  if (dunn) {
    # Passes options that make wilcox.test behave as the Dunn test
    cor0 <- F
    ex <- F
  } else {
    cor0 <- T
    ex <- NULL
  }

  # Calcs
  p <- pairwise.wilcox.test(mdf[[yvar]], g = mdf[[xvar]], correct = cor0, exact = ex)[[3]]
  if (raw) pl <- p

  # Formatting
  if (tpose) p <- TransposeUnevenMatrix(p, fill = F)
  if (method == 'star_notation') {
    pw <- p
    pc <- as.character(p)
    pc <- as.character(p)
    pc[which(p >= alpha[1])] <- 'NS'
    pc <- ifelse(is.na(p), '', pc)
    pc <- ifelse(p < alpha[1], '*', pc)
    pc <- ifelse(p < alpha[2], '**', pc)
    pc <- ifelse(p < alpha[3], '***', pc)
    p <- matrix(data = pc, nrow = nrow(pw))
    row.names(p) <- row.names(pw)
    colnames(p) <- colnames(pw)
  } else if (method == 'logical_sig') {
    p <- ifelse(p > alpha[1], 0, 1)
    p <- as.logical(p)
  } else if (method == 'num_pval') {
    p <- round(p, 3)
  } else {
    stop('allowed methods are star_notation, logical_sig, num_pval')
  }

  if (xvar == 'section') {
    p <- data.frame(p, stringsAsFactors = F)
    cc <- ClelandEcoregions::ScaleUpClelandName(row.names(p), 'section', belt = 'M332D')
    p <- cbind(cc, row.names(p), p)
    row.names(p) <- NULL
    colnames(p)[1:2] <- c('province', 'section')
  }

  # Return
  ftag <- paste0('pwise_wilxcox_', yvar, '_', xvar, '.csv')
  if (write) {
    message('Writing to:')
    cat(out_dir, '/', ftag, '\n')
  }
  ftag <- ifelse(dunn, paste0('dunn_', ftag), paste0('wilcox', ftag))
  if (raw) {
    if (write) write.csv(pl, paste0('long_', ftag))
    invisible(pl)
  } else {
    if (write) write.csv(p, ftag)
    invisible(p)
  }
}
#' @describeIn CalcPairwiseWilcox Wrapper for significance group calculations
#' @family package_utilities
#' @export
CalcSectionWilcox <- function(x) {
  # CalcPairwiseWilcox was written with different formatting in mind, this fixes
  wm <- CalcPairwiseWilcox(xvar = x, yvar = 'mort_rate', method = 'star_notation')
  if (x == 'section') {
    row.names(wm) <- wm[, 2]
    wm <- wm[, -c(1, 2)]
    colnames(wm) <- gsub('\\.', ' ', colnames(wm))
    wm <- as.matrix(wm)
  }
  #wm <- TransposeUnevenMatrix(wm, fill = T)
  # T = groups are not significantly different, i.e. should be grouped
  #wm <- ifelse(wm < 0.05, F, T)
  #diag(wm) <- T
  return(wm)
}
#' @describeIn CalcPairwiseWilcox Wrapper for some experimental stuff
#' @family package_utilities
#' @export
CalcWilcoxGroups <- function() {
  wm <- CalcPairwiseWilcox(xvar = 'section', yvar = 'mort_rate')
  row.names(wm) <- wm[, 2]
  wm <- wm[, -c(1, 2)]
  colnames(wm) <- gsub('\\.', ' ', colnames(wm))
  wm <- as.matrix(wm)
  wm <- TransposeUnevenMatrix(wm, fill = T)
  # T = groups are not significantly different, i.e. should be grouped
  wm <- ifelse(wm < 0.05, F, T)
  diag(wm) <- T
  z <- GroupLogicalMatrix(wm)
  z
}