R/FilterMortByCondition.R

In bmcnellis/RSFIA: FIA Tree Mortality Analysis

#' @title Filter FIA_mortality* by COND table rules
#'
#' @description
#'
#' Some conditions needs to be filtered out - running this trims
#' down the dataset.
#'
#' @param dir Directory of *_COND.csv files
#' @param plot_type Either 'prev' for PREV_PLT_CN, or 'post' for PLT_CN
#' @param db_ver Database version
#'
#' @details
#'
#' Rules relevant to methods:
#'
#' Conditions included had to meet COND_STATUS_CD == 1, i.e. accessible forest land
#'
#' For plots with more than one forest type, sub-dominant (< sub_dom_cutoff)
#' forest types were excluded. If, after applying this filter, plots still had more
#' than one forest type, plots were excluded.
#'
#' For plots with multiple conditions, condition values were averaged.
#'
#'
#' @export
#' @examples FilterMortByCondition()
FilterMortByCondition <- function(dir = getwd(), data = NULL,
                                  plot_type = 'prev', db_ver = 7.2,
                                  sub_dom_cutoff = 0.75) {

  # Import data:
  # C:/Users/Brandon/Documents/docs/PHD/FIA/data
  plot_CNs <- integer()
  if (plot_type == 'prev') {
    if (length(data) > 0) {
      plot_CNs <- data$PREV_PLT_CN
    } else {
      plot_CNs <- FIA_mortality_with_explanatory$PREV_PLT_CN
    }
  } else if (plot_type == 'post') {
    if (length(data) > 0) {
      plot_CNs <- data$PLT_CN
    } else {
      plot_CNs <- FIA_mortality_with_explanatory$PLT_CN
    }
  } else {
    stop('bad plot_type input, see documentation')
  }
  if (plot_type == 'prev') {
    plot_CNs <- FIA_mortality_with_explanatory$PREV_PLT_CN
  } else if (plot_type == 'post') {
    plot_CNs <- FIA_mortality_with_explanatory$PLT_CN
  } else {
    stop('bad plot_type input, see documentation')
  }
  cond_df <- ImportConditions(dir = dir, filt_PLT_CN = plot_CNs)
  # Bulk condition filters:
  cond_df <- cond_df[which(cond_df$COND_STATUS_CD == 1), ]

  # Filter section:
  incl_tags <- c('PLT_CN', 'CN', 'FORTYPCD', 'STDAGE', 'STDSZCD', 'SITECLCD', 'STDORGCD',
                 'SLOPE', 'ASPECT', 'PHYSCLCD', 'ALSTKCD', 'CONDPROP_UNADJ',
                 'DSTRBCD1', 'DSTRBCD2', 'DSTRBCD3',
                 'DSTRBYR1', 'DSTRBYR2', 'DSTRBYR3',
                 'COND_STATUS_CD')
  for_typ <- data.frame(matrix(nrow = 0, ncol = length(incl_tags)))
  colnames(for_typ) <- incl_tags
  n_multi_con_drop <- integer()
  i0 <- length(unique(data$CN))

  cat('\n')
  message('Filter progress:')
  for (i in 1:i0) {
    cat('\r', format(i / i0 * 100, digits = 2, nsmall = 2), '%    ')
    i_CN <- unique(data$CN)[i]
    i_con <- cond_df[which(cond_df$PLT_CN == i_CN), incl_tags]
    #if (is.na(i_con$FORTYPCD)) browser()
    if (nrow(i_con) == 0) next
    if (nrow(i_con) > 2) {
      ag_for_typ <- aggregate(i_con$CONDPROP_UNADJ, by = list(i_con$FORTYPCD), FUN = sum)
      ag_prop <- ag_for_typ$x / sum(i_con$CONDPROP_UNADJ)
      if (any(ag_prop >= sub_dom_cutoff)) {
        i_con <- i_con[which(i_con$FORTYPCD == ag_for_typ$Group.1[which.max(ag_prop)]), ]
      }
      ag_for_typ <- aggregate(i_con$CONDPROP_UNADJ, by = list(i_con$FORTYPCD), FUN = sum)
      if (nrow(ag_for_typ) == 1) {
        ag_con <- apply(i_con, 2, function(x) {
          a <- mean(x, na.rm = T)
          b <- ifelse(is.nan(a), NA, a)
          c <- round(b, 2)
        })
        i_con[1, names(ag_con)] <- ag_con
        i_con <- i_con[1, ]
      }
    }
    if (nrow(i_con) > 2) {
      n_multi_con_drop <- append(n_multi_con_drop, i_CN)
      next
    }
    if (nrow(i_con) == 2) {
      i_con <- i_con[which.max(i_con$CONDPROP_UNADJ), ]
    }
    if (nrow(i_con) == 1) {
      if (!(i_con$COND_STATUS_CD %in% c(1, 2))) next
      for_typ <- plyr::rbind.fill(for_typ, i_con)
    }
    if (nrow(i_con) > 1) {
      stop('filtering error')
    }
  }
  colnames(for_typ)[which(colnames(for_typ) == 'CN')] <- 'COND_CN'

  cat('\n')
  message('Done!')
  cat('Sub-dominance forest type cutoff:', sub_dom_cutoff, '\n')
  cat('Number of multi-condition plots dropped:', length(n_multi_con_drop), '\n')
  perc_multi_drop <- round(length(n_multi_con_drop) / nrow(for_typ) * 100, 2)
  cat('Percent of multi-condition plots dropped:', perc_multi_drop, '\n')
  return(for_typ)
}