R/grsIn.R

In dcarver1/gapAnalysisR: An R package to calculate conservation indicators using species distribution models

#' @title Germplasm representativeness score estimation (In-situ conservation indicators).
#' @name grs_insitu
#' @description Performs an estimation of germplasm representativeness score for in-situ gap analysis (GRSin) using Khoury et al., (2019) methodology
#' This function uses a germplasm buffer raster file (e.g. CA50), a thresholded species distribution model, and a raster file of protected areas
#'  \deqn{GRSin = min(100,(Germplasm buffer area into protected area/species distribution area in protected areas)*100)}
#'
#' @param species A name species compiled using '_'  to call the raster files (SDM and germplasm buffer)
#'  from Workspace/parameters/inputs folder and occurrences files from Workspace/parameter/occurrences folder
#' @param Workspace A forder where the pipeline will be executed
#' @param  run_version The version of the analysis used (e.g 'v1')
#'
#' @return It returns a raster file with the species distribution model restricted to the protected areas raster file provided.
#' Also, this function returns a data frame file saved at gap_analysis folder with four columns:
#'
#' \tabular{lcc}{
#'  ID \tab Species name \cr
#'  SPP_AREA_km2 \tab Area occupied by the species using as input a SDM thresholded file in tiff format \cr
#'  SPP_WITHIN_PA_AREA_km2 \tab Area occupied by the germplasm accessions in a species distribution model \cr
#'  GRS \tab GRSex result \cr
#' }
#'
#' @examples grs_exsitu('Cucurbita_digitata',Workspace,'v1')
#'  \dontrun
#'  Workspace  <-  'E:/CIAT/workspace/Workspace_test/workspace'
#'  run_version  <- 'v1'
#'  species_list <- c('Cucurbita_cordata',
#'   'Cucurbita_digitata',
#'   'Cucurbita_foetidissima',
#'   'Cucurbita_palmata')
#'
#'   run_version <-'v1'
#
#'  lapply(1:length(species_list),function(i){
#'     species <- species_list[[i]]
#'     x <- grs_insitu(species,Workspace,run_version)
#'     print(paste0(species,' DONE!'))
#'  })
#'
#'@references
#'
#'Ramírez-Villegas, J., Khoury, C., Jarvis, A., Debouck, D. G., & Guarino, L. (2010).
#'A Gap Analysis Methodology for Collecting Crop Genepools: A Case Study with Phaseolus Beans.
#'PLOS ONE, 5(10), e13497. Retrieved from https://doi.org/10.1371/journal.pone.0013497
#'
#' Khoury, C. K., Amariles, D., Soto, J. S., Diaz, M. V., Sotelo, S., Sosa, C. C., … Jarvis, A. (2019).
#' Comprehensiveness of conservation of useful wild plants: An operational indicator for biodiversity
#' and sustainable development targets. Ecological Indicators. https://doi.org/10.1016/j.ecolind.2018.11.016
#'
#' @export

insitu_grs = function(species_list,occurrenceData,raster_list){

suppressMessages(require(rgdal))
suppressMessages(require(raster))
suppressMessages(require(tmap))
suppressMessages(require(fasterize))
suppressMessages(require(sf))

#importFrom("methods", "as")
#importFrom("stats", "complete.cases", "filter", "median")
#importFrom("utils", "data", "memory.limit", "read.csv", "write.csv")
  df <- data.frame(matrix(ncol=2, nrow = length(species_list)))
  colnames(df) <- c("species", "GRSin")
  # load in protect area raster
  proArea <- raster(system.file("data/protectedArea/wdpa_reclass.tif",
                                package = "gapAnalysisR"))
  # loop over species list
  for(i in 1:length(species_list)){
    # select threshold map for a given species
    for(j in 1:length(raster_list)){
      if(grepl(j, i, ignore.case = TRUE)){
        sdm <- raster_list[[j]]
      }
    }
    # determine the area of predicted presence of a species based on the threshold map
    sdm1 <- sdm
    proArea1 <- raster::crop(x = proArea,y = sdm1)
    sdm1[sdm1 == 0] <- NA
    cell_size <- raster::area(sdm1, na.rm=TRUE, weights=FALSE)
    cell_size <- cell_size[!is.na(cell_size)]
    thrshold_area <- length(cell_size)*median(cell_size)

    # mask the protected area Raster to the threshold map and calculate area
    proArea1[proArea1 == 0] <-NA
    proArea1 <- proArea1 * sdm1
    # calculate area
    cell_size <- raster::area(proArea1, na.rm=TRUE, weights=FALSE)
    cell_size <- cell_size[!is.na(cell_size)]
    protected_area <- length(cell_size)*median(cell_size)
    if(!is.na(protected_area)){
      # calculate GRSin
      grs <- min(c(100, protected_area/thrshold_area*100))
      df$species[i] <- as.character(species_list[i])
      df$GRSin[i] <- grs
    }else{
      df$species[i] <- as.character(species_list[i])
      df$GRSin[i] <- 0
    }
  }
  return(df)
}