R/RSFIA_stats.R

In bmcnellis/RSFIA: FIA Tree Mortality Analysis

#' @title Stat analysis functions for RSFIA data
#'
#' @description
#'
#' These functions print stats/return tables/etc for RSFIA data
#' @name RSFIA_stats
NULL
#' @describeIn RSFIA_stats Prep function for running RSFIA_stats functions
PrepStatEnvir <- function() {
  require(RSFIA)
  require(ClelandEcoregions)
  data('FIA_mortality_with_explanatory')
  mort_df <- FIA_mortality_with_explanatory
  return(mort_df)
}
#' @describeIn RSFIA_stats Calculates summary statistics
CalcSummary <- function(x, y = NULL) {
  if (inherits(x, 'formula')) {
    mod <- lm(x)
  } else {
    mod <- lm(y ~ x)
  }
  message('\nSummary:')
  print(summary(mod))
  message('\nNormality:')
  if (length(mod$resid) < 5001) {
    print(shapiro.test(mod$resid))
  } else {
    cat('Sample size too large for shapiro.test\n')
  }
  message('\nAOV:')
  print(anova(mod))
}
#' @describeIn RSFIA_stats Mortality affects at ecosystem level
CalcEcosystemMort <- function() {
  # Prep function environment:
  mort_df <- RSFIA::PrepStatEnvir()
  mean_NA <- function(x) suppressWarnings(mean(x, na.rm = T))
  sd_NA <- function(x) suppressWarnings(sd(x, na.rm = T))
  message('Ecosystem code table:\n')
  print(table(mort_df$Cleland_province))
  prov_name <- RSFIA::KeyClelandCode(mort_df$Cleland_province, lvl = 'province')
  message('\nEcosystem name table:\n')
  print(table(prov_name))
  mort_df <- data.frame(mort_df, prov_name, stringsAsFactors = F)
  state <- KeyStateCode(mort_df$STATECD, db_ver = 7.2)
  # ANOVAs
  browser()
  # Overall mort:
  mort_agg <- aggregate(mort_df, by = list(prov_name), FUN = mean_NA)
  mort_agg <- lapply(mort_agg, function(x) ifelse(is.nan(x), NA, x))
  mort_agg <- data.frame(mort_agg, stringsAsFactors = F)
  colnames(mort_agg)[1] <- 'Cleland_province'
  mort_agg <- mort_agg[, sapply(mort_agg, function(x) !all(is.na(x)))]
  browser()
  invisible()
}