Load climate predictions

pred_temp_do <- readRDS("data/predicted-DO-new.rds")

Choose settings spatial and temporal settings

region <- "West Coast Vancouver Island"
trend <- TRUE # IF TREND will automatically apply different thresholds for average change

# set years if not using an average trend
if (!trend) { 
start_time <- 2008
end_time <- 2011
} 

Choose what variables to base vectors on

# JUST TEMPERATURE
climate <- "temperature"
temp_threshold <- c(0.75)
optimal_temp <- 6.5
temp_trend_threshold <- c(0.25)

# # JUST DO
# climate <- "do"
# do_threshold <- c(0.25)
# optimal_do <- 2 
# do_trend_threshold <- c(0.1)  

# # BOTH
# climate <- "do and temperature"
# max_thresholds <- c(Inf, 0.75)
# min_thresholds <- c(0.25, Inf)
# optimal_value <- c(2, 6.5)
# # should be lower than used for raw changes
# min_trend_thresholds <- c(0.1, Inf)
# max_trend_thresholds <- c(Inf, 0.25)

Run all subsequent code...

library(dplyr)
library(ggplot2)
library(gfplot)
library(gfranges)
if (region == "Both odd year surveys") {
  survey <- c("SYN QCS", "SYN HS")
  model_ssid <- c(1, 3)
  ssid_string <- paste0(model_ssid, collapse = "n")
  trend_indices_do <- c(1, 1, 1, 2, 2)
  trend_indices_temp <- c(1, 1, 1, 2, 2)
  years <- NULL
}

if (region == "West Coast Vancouver Island") {
  survey <- c("SYN WCVI")
  model_ssid <- c(4)
  ssid_string <- paste0(model_ssid, collapse = "n")
  trend_indices_do <- c(1, 1, 1, 2, 2)
  trend_indices_temp <- c(1, 1, 1, 2, 2, 2)
  years <- NULL
}

if (region == "West Coast Haida Gwaii") {
  survey <- c("SYN WCHG")
  model_ssid <- c(16)
  ssid_string <- paste0(model_ssid, collapse = "n")
  trend_indices_do <- c(1, 1, 1, 2, 2)
  trend_indices_temp <- c(1, 1, 1, 2, 2)
  years <- NULL
}
# Set default cell sizes and scale
input_cell_size <- 2
scale_fac <- 1
delta_t_step <- 2
skip_time <- NULL

if (climate == "temperature") {
  variable_names <- c("temp")
  min_thresholds <- c(Inf)
  delta_threshold <- temp_threshold
  optimal_value <- optimal_temp

  if (trend) {
    max_thresholds <- temp_trend_threshold
    start_time <- 2008 # change if can add back in prior 2008 data
    end_time <- NULL
    delta_t_total <- 6 # change if can add back in prior 2008 data
    indices <- trend_indices_temp
  } else {
    max_thresholds <- temp_threshold
    indices <- c(1, 2)
    delta_t_total <- 2
  }
}

if (climate == "do") {
  variable_names <- c("do_est")
  max_thresholds <- c(Inf)
  delta_threshold <- -1*(do_threshold)
  optimal_value <- optimal_do

  if (trend) {
    min_thresholds <- do_trend_threshold
    start_time <- 2008
    skip_time <- 2016
    end_time <- NULL
    delta_t_total <- 6
    indices <- trend_indices_do
  } else {
    min_thresholds <- do_threshold
    indices <- c(1, 2)
    delta_t_total <- 2
  }
}

if (climate == "do and temperature") {
  variable_names <- c("do_est", "temp")

  if (trend) {
    min_thresholds <- min_trend_thresholds
    max_thresholds <- max_trend_thresholds
    start_time <- 2008
    skip_time <- 2016
    end_time <- NULL
    delta_t_total <- 6
    indices <- trend_indices_do
  } else {
    min_thresholds <- min_thresholds
    max_thresholds <- max_thresholds
    indices <- c(1, 2)
    delta_t_total <- 2
  }
}

min_string <- paste0(min_thresholds, collapse = "-")
max_string <- paste0(max_thresholds, collapse = "-")
optimal_string <- paste0(optimal_value, collapse = "-")
climate_string <- gsub(" ", "-", gsub("\\/", "-", tolower(climate)))

Calculate climate vectors

if (length(variable_names) > 1) {
  climate_predictions <- replicate(length(variable_names), pred_temp_do, simplify = FALSE)
} else {
  climate_predictions <- pred_temp_do
}

starttime <- Sys.time()
vocc <- make_vector_data(climate_predictions,
  variable_names = variable_names,
  ssid = model_ssid,
  start_time = start_time,
  end_time = end_time,
  skip_time = skip_time,
  input_cell_size = input_cell_size,
  scale_fac = scale_fac,
  delta_t_total = delta_t_total,
  delta_t_step = delta_t_step,
  indices = indices,
  min_thresholds = min_thresholds,
  max_thresholds = max_thresholds,
  round_fact = 10
)
endtime <- Sys.time()
time_vocc <- round(starttime - endtime)

vocc$var_1_min <- min_thresholds[1]
vocc$var_2_min <- min_thresholds[2]
vocc$var_1_max<- max_thresholds[1]
vocc$var_2_max<- max_thresholds[2]

#glimpse(vocc)

Plot climate vectors

if(climate== "temperature") {
vocc$diff <- vocc$units_per_decade/10*delta_t_total

vocc_plot <- plot_vocc(vocc,
  vec_aes = "distance",
  #max_vec_plotted = 100,
  min_vec_plotted = input_cell_size,
  NA_label = ".",
  fill_col = "diff",
  fill_label = paste("Change from\n", start_time, "to", start_time + delta_t_total, ""),
  raster_alpha = 1,
  white_zero = TRUE,
  vec_alpha = 0.25,
  axis_lables = FALSE,
  transform_col = no_trans
) + ggtitle(paste("Temperature within", max_thresholds, "°C"))
}

if (climate == "do") {
vocc$diff <- vocc$units_per_decade/10*delta_t_total

vocc_plot <- plot_vocc(vocc,
  vec_aes = "distance",
  #max_vec_plotted = 100,
  min_vec_plotted = input_cell_size,
  NA_label = ".",
  fill_col = "diff",
  fill_label = paste("Change from\n", start_time, "to", start_time + delta_t_total, ""),
  raster_alpha = 1,
  white_zero = TRUE,
  high_fill = "Steel Blue 4",
  low_fill = "Red 3",
  vec_alpha = 0.25,
  axis_lables = FALSE,
  transform_col = no_trans
) + ggtitle(paste("DO within", min_thresholds, "ml/L") )
}

if (climate == "do and temperature") {


vocc$temp_diff <- vocc$temp.units_per_decade/10*delta_t_total

vocc_plot_temp <- plot_vocc(vocc,
  vec_aes = "distance",
  #max_vec_plotted = 100,
  min_vec_plotted = input_cell_size,
  NA_label = ".",
  fill_col = "temp_diff",
  fill_label = paste("Change from\n", start_time, "to", start_time + delta_t_total, ""),
  raster_alpha = 1,
  white_zero = TRUE,
  vec_alpha = 0.25,
  axis_lables = FALSE,
  transform_col = no_trans
)+ ggtitle(paste("Temperature within", max_thresholds[2], "°C"))

vocc$do_diff <- vocc$do_est.units_per_decade/10*delta_t_total

vocc_plot_do <- plot_vocc(vocc,
  vec_aes = "distance",
  #max_vec_plotted = 100,
  min_vec_plotted = input_cell_size,
  NA_label = ".",
  fill_col = "do_diff",
  fill_label = paste("Change from\n", start_time, "to", start_time + delta_t_total, ""),
  raster_alpha = 1,
  white_zero = TRUE,
  high_fill = "Steel Blue 4",
  low_fill = "Red 3",
  vec_alpha = 0.25,
  axis_lables = FALSE,
  transform_col = no_trans
) + ggtitle(paste("DO within", min_thresholds[1], "ml/L") )

vocc_plot <- gridExtra::grid.arrange(
  vocc_plot_temp, vocc_plot_do,
  ncol = 2,
  top = grid::textGrob("VOCC: vectors based on two climate variables at once")
)
}

vocc_plot


pbs-assess/gfranges documentation built on Dec. 13, 2021, 4:50 p.m.