R/loocv_gear_temp.R

In afsc-gap-products/coldpool: AFSC/RACE Groundfish Assessment Program EBS and NBS temperature products

#' Compare bottom temperature interpolation methods using leave-one-out cross validation
#' 
#' Wrapper function for coldpool::loocv_2(), to run leave-one-out cross validation.
#' 
#' @param temp_data_path Filepath to temperature csv file generated by get_data()
#' @param proj_crs CRS for the project as a character vector.
#' @param interp_variable Name of the variable from the csv to be interpolated.
#' @export

loocv_gear_temp <- function(temp_data_path, 
                       proj_crs,
                       interp_variable = "gear_temperature") {
  
  if(!(class(temp_data_path) == "data.frame")) {
    temperature_df <- read.csv(file = temp_data_path,
                               stringsAsFactors = FALSE)
  }

  names(temperature_df) <- tolower(names(temperature_df))
  
  # Vector of years ----
  year_vec <- sort(unique(temperature_df$year))
  
  # Leave-one-out cross validation for each year ----
  # Sean Rohan's notes: 
  #  - loocv_2() doesn't model anisotropy. A redesign of loocv_2 using functions in geoR instead of gstat may allow a way to automatically estimate anisotropy.
  #  - Transformations (normal score, Box-Cox) may improve precision. Can investigate further if warranted.
  #  - For propagating uncertainty, open-source 2D empirical Bayesian kriging may be better... TBD.
  
  print("Starting cross validation using loocv_2()")
  for(i in 1:length(year_vec)) {
    gear_temp_raster <- coldpool::loocv_2(dat = filter(temperature_df, year == year_vec[i]),
                                             in.proj = "+proj=longlat",
                                             interp.proj = proj_crs,
                                             lon.col = "longitude",
                                             lat.col = "latitude",
                                             var.col = interp_variable,
                                             pre = paste0(interp_variable, "_", year_vec[i]),
                                             scale.vars = FALSE,
                                             center = FALSE,
                                             scale = FALSE)
  }
  
  # Generate vector of loocv files with "GEAR_TEMPERATURE" in the title.
  ifelse(!dir.exists(file.path(here::here("output"))), dir.create(file.path(here::here("output"))), FALSE)
  ifelse(!dir.exists(file.path(here::here("output", "loocv"))), dir.create(file.path(here::here("output", "loocv"))), FALSE)
  ifelse(!dir.exists(file.path(here::here("plots"))), dir.create(file.path(here::here("plots"))), FALSE)
  temp_dir <- dir(here::here("output", "loocv"), full.names = TRUE)
  temp_dir <- temp_dir[grep(interp_variable, temp_dir)]
  
  # Prediction error plots ----
  print("Making RSPE violin plots")
  for(i in seq(1, length(year_vec), by = 9)){
    plot_loocv_rmse(sel_paths = temp_dir[i:min(i+8, length(year_vec))],
                    y_lab = expression(RSPE~(degree*C)),
                    sel_var = interp_variable,
                    make_plot = TRUE,
                    by_cruise = FALSE,
                    suffix = paste0("_", match(i, seq(1, length(year_vec), by = 9))),
                    fig_res = 600)
  }
}