In pbs-assess/gfranges: Groundfish Ranges
knitr::opts_chunk$set(
fig.keep = "all", fig.show = "asis", dpi = 600,
fig.width = 6, fig.height = 5,
echo = FALSE, warning = FALSE, message = FALSE
)
options(scipen = 999)
options(digits = 4)
library("gfranges")
library("dplyr")
library("ggplot2")
ggplot2::theme_set(gfplot::theme_pbs())
Set parameters
species <- params$species
region <- params$region
covs <- params$covs
climate_mod <- params$climate_mod
immature <- params$immature
scrambled_years <- params$scrambled_years
spp <- gsub(" ", "-", gsub("\\/", "-", tolower(species)))
paste("species =", species)
paste("region =", region)
paste("covs =", covs)
paste("scrambled =", scrambled_years)
Set cell sizes and scale
dist_intercept <- 0 # input_cell_size / 2
input_cell_size <- 2
scale_fac <- 2 # means that the raster is reprojected to _ X original grid (2 km)
if (region == "both odd year surveys") {
survey <- c("SYN QCS", "SYN HS")
model_ssid <- c(1, 3)
ssid_string <- paste0(model_ssid, collapse = "n")
start_year <- 2008
}
if (region == "West Coast Vancouver Island") {
survey <- c("SYN WCVI")
model_ssid <- c(4)
ssid_string <- paste0(model_ssid, collapse = "n")
start_year <- 2008
}
if (region == "West Coast Haida Gwaii") {
survey <- c("SYN WCHG")
model_ssid <- c(16)
ssid_string <- paste0(model_ssid, collapse = "n")
start_year <- 2008
}
rm(d_all)
if (immature) {
age <- "immature"
d_all <- readRDS(paste0(
"data/", spp,
"/predictions-", spp, covs, "-1n3n4n16-imm-biomass.rds" #-prior-FALSE.rds"
)) %>% filter(year >= start_year )
} else {
age <- "mature"
try({
d_all <- readRDS(paste0(
"data/", spp,
"/predictions-", spp, covs, "-1n3n4n16-total-biomass.rds"
)) %>% filter(year >= start_year )
})
try({
d_all <- readRDS(paste0(
"data/", spp,
"/predictions-", spp, covs, "-1n3n4n16-mature-biomass-new.rds"
)) %>% filter(year >= start_year )
})
}
# hist(d_all$log_effort)
SET TIME FRAME
multiyear <- TRUE
if (multiyear) {
# unique(d_all$year)
d <- d_all %>% filter(ssid %in% model_ssid)
length(unique(d$year))
years_sampled <- length(unique(d$year))
year_range <- (years_sampled - 1) * 2
# # alternate version of this model
# if (years_sampled == 5) indices <- c(1, 2, 2, 2, 2)
# if (years_sampled == 6) indices <- c(1, 2, 2, 2, 2, 2)
## USED for truncating distance-based velocities?
max_dist_all_years <- 10 * year_range
# should try running with gradient based on only later years?
if (years_sampled == 5) { indices <- c(1, 1, 1, 2, 2) }
if (years_sampled == 6) { indices <- c(1, 1, 1, 1, 2, 2) }
} else {
d <- d_all %>%
filter(year >= start_year) %>%
filter(year <= start_year + 3)
unique(d$year)
indices <- c(1, 2)
year_range <- 2
}
skip_time <- NULL
delta_t_step <- year_range
delta_t_total <- year_range
glimpse(d)
if (scrambled_years) {
d2 <- list()
# set.seed(9999)
# set.seed(12345)
# set.seed(98765)
original_years <- unique(d$year)
new_years <- sample(original_years, replace = F)
convert <- cbind(original_years, new_years)
convert
for (i in seq_len(length(original_years))) {
old_year <- as.integer(convert[[i, 1]])
fake_year <- as.integer(convert[[i, 2]])
d2[[i]] <- filter(d, year == old_year) %>%
mutate(year = fake_year)
# browser()
}
d <- do.call(rbind, d2)
convert
}
Load climate data
nd_all <- readRDS(paste0("data/predicted-DO-", climate_mod, "-roms.rds")) %>% select(X, Y, year, temp, do_est, roms_temp)
d <- left_join(d, nd_all, by = c("X", "Y", "year"))
glimpse(d)
Calculate climate velocity
Temperature
rm(gf2)
try({
gf2 <- readRDS(paste0("data/_climate_vocc/temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
})
if(!exists("gf2")){
gf2 <- vocc_gradient_calc(d, "temp",
scale_fac = scale_fac,
# all layers with max indice will be averaged for spatial gradient
# if default (NULL) will use all years
indices = indices,
# divisor = 1, #default of 10 gives per decade
quantile_cutoff = 0.05
)
gf2 <- gf2 %>%
mutate(
## these are only relevant if divisor = 1
# trend_per_year = trend,
# trend_per_decade = trend_per_year * 10,
trend_per_year = trend/10,
trend_per_decade = trend_per_year * 10,
velocity = velocity
) %>%
dplyr::mutate(X = x, Y = y)
# FULL DATASET
# mean(gf2$sd, na.rm = T) # [1] 0.3434
# min(gf2$sd, na.rm = T) # [1] 0.1266
# max(gf2$sd, na.rm = T) # [1] 1.12
# median(gf2$sd, na.rm = T) # [1] 0.287
#
# 2008 DATASET
# values currently depend on survey region
temp_mean_sd <- mean(gf2$sd, na.rm = T) # [1] 0.3274
temp_mean_sd
temp_median_sd <- median(gf2$sd, na.rm = T) # [1] 0.2882
temp_median_sd
temp_min_sd <- min(gf2$sd, na.rm = T) # [1] 0.1035
temp_min_sd*2
max(gf2$sd, na.rm = T) # [1] 1.248
# temp_sd_choice <- temp_min_sd*2
# temp_sd_choice <- temp_mean_sd
temp_sd_choice <- temp_median_sd
# lower_change_t <- temp_sd_choice
lower_change_t <- Inf
upper_change_t <- temp_sd_choice
rm(dvocc_t)
dvocc_t <- make_vector_data(d,
variable_names = c("temp"),
ssid = model_ssid,
input_cell_size = 2,
scale_fac = scale_fac,
min_dist = 0,
delta_t_total = years_sampled,
delta_t_step = 2,
indices = indices,
min_thresholds = lower_change_t,
max_thresholds = upper_change_t,
round_fact = 10
)
dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope)
dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10,
(max_dist_all_years/ year_range)*10, dvocc))
dvocc_t <- mutate(dvocc_t, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10,
-(max_dist_all_years / year_range)*10, dvocc))
dv_t <- dvocc_t %>% select(icell, x, y, dvocc) %>% unique()
# glimpse(dv_t)
# glimpse(gf2)
# gf2$dvocc <- dv_t$dvocc
gf2 <- left_join(gf2, dv_t)
saveRDS(gf2, file = paste0("data/_climate_vocc/temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# saveRDS(gf2, file = paste0("data/_climate_vocc/temp-dvocc-min2SD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# saveRDS(gf2, file = paste0("data/_climate_vocc/temp-dvocc-medSDasy-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
}
# mean(gf2$sd, na.rm = T) # [1] 0.3274
temp_min_sd <- round(min(gf2$sd, na.rm = T), digits = 3) #
# temp_mean_sd <- mean(gf2$sd, na.rm = T) #
# max(gf2$sd, na.rm = T) # [1] 1.248
temp_median_sd <- round(median(gf2$sd, na.rm = T), digits = 3) # [1] 0.2882
temp_mean_sd <- round(mean(gf2$sd, na.rm = T), digits = 3) #
print( paste("min SD =", round(temp_min_sd, digits = 2)))
print( paste("median SD =", round(temp_median_sd, digits = 2)))
print( paste("mean SD =", round(temp_mean_sd, digits = 2)))
print( paste("treshold chosen =", round(temp_mean_sd, digits = 2), "degrees C"))
plot_vocc(gf2,
raster_cell_size = 4,
coast = TRUE,
vec_aes = NULL,
NA_label = "NA",
fill_col = "dvocc",
fill_label = "temperature \nkm/decade",
raster_alpha = 1,
na_colour = "white",
mid_fill = "mistyrose1", grey_water = F,
# mid_fill = "lavenderblush1", grey_water = F,
low_fill = "royalblue4", # low_fill = "#5E4FA2",
high_fill = "Red 3",
raster_limits = c(0, 100),
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "distance-based temp velocity",
subtitle = paste0(min(d$year), "-", max(d$year), " mean temp SD = ", signif(temp_mean_sd, digits = 2), " for min threshold of change")
)
Annual temperature (ROMS)
rm(roms)
try({roms <- readRDS(paste0("data/_climate_vocc/annual-temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
})
if(!exists("roms")){
roms <- vocc_gradient_calc(d, "roms_temp",
scale_fac = scale_fac,
# all layers with max indice will be averaged for spatial gradient
# if default (NULL) will use all years
indices = indices,
# divisor = 1, #default of 10 gives per decade
quantile_cutoff = 0.05
)
roms <- roms %>%
mutate(
## these are only relevant if divisor = 1
# trend_per_year = trend,
# trend_per_decade = trend_per_year * 10,
trend_per_year = trend/10,
trend_per_decade = trend_per_year * 10,
velocity = velocity
) %>%
dplyr::mutate(X = x, Y = y)
# # FULL DATASET
# # mean(roms$sd, na.rm = T) # [1] 0.3434
# # min(roms$sd, na.rm = T) # [1] 0.1266
# # max(roms$sd, na.rm = T) # [1] 1.12
# # median(roms$sd, na.rm = T) # [1] 0.287
# #
# # 2008 DATASET
# mean(roms$sd, na.rm = T) # [1] 0.3274
# temp_min_sd <- min(roms$sd, na.rm = T) # [1] 0.1035
# max(roms$sd, na.rm = T) # [1] 1.248
# temp_median_sd <- median(roms$sd, na.rm = T) # [1] 0.2882
#
# temp_sd_choice <- temp_min_sd
# lower_change_t <- temp_sd_choice*2
# upper_change_t <- temp_sd_choice*2
#
# dvocc_t <- make_vector_data(d,
# variable_names = c("temp"),
# ssid = model_ssid,
# input_cell_size = 2,
# scale_fac = scale_fac,
# min_dist = 0,
# delta_t_total = years_sampled,
# delta_t_step = 2,
# indices = indices,
# min_thresholds = lower_change_t,
# max_thresholds = upper_change_t,
# round_fact = 10
# )
#
# dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope)
# dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10,
# (max_dist_all_years/ year_range)*10, dvocc))
# dvocc_t <- mutate(dvocc_t, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10,
# -(max_dist_all_years / year_range)*10, dvocc))
#
# dv_t <- dvocc_t %>% select(icell, x, y, dvocc) %>% unique()
#
# # glimpse(dv_t)
# # glimpse(roms)
#
# # roms$dvocc <- dv_t$dvocc
#
# roms <- left_join(roms, dv_t)
saveRDS(roms, file = paste0("data/_climate_vocc/annual-temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# mine <- gf2 %>% select(x, y, temp_trend = trend, temp_vel = velocity, temp_grad= gradient)
# other <- roms %>% select(x, y, roms_temp_trend = trend, roms_temp_vel = velocity, roms_temp_grad= gradient)
#
# check <- left_join(mine, other)
#
# ggplot(check, aes(roms_temp_trend, temp_trend)) +
# geom_point(alpha=0.2) +
# coord_cartesian(xlim = c(-1,2.5), ylim = c(-1,2.5)) +
# geom_abline(intercept = 0, slope = 1) +
# gfplot::theme_pbs()
}
# ggplot(check, aes(roms_temp_vel, temp_vel)) +
# geom_point(alpha=0.2) +
# coord_cartesian(xlim = c(-20,80), ylim = c(-20,80)) +
# geom_abline(intercept = 0, slope = 1) +
# gfplot::theme_pbs()
Dissolved Oxygen
rm(gf0)
try({
gf0 <- readRDS(paste0("data/_climate_vocc/DO-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
})
if(!exists("gf0")){
gf0 <- vocc_gradient_calc(d, "do_est",
scale_fac = scale_fac,
# all layers with max indice will be averaged for spatial gradient
# if default (NULL) will use all years
indices = indices,
# divisor = 1, #default of 10 gives per decade
quantile_cutoff = 0.05
)
gf0 <- gf0 %>%
mutate(
## these are only relevant if divisor = 1
# trend_per_year = trend,
# trend_per_decade = trend_per_year * 10,
trend_per_year = trend/10,
trend_per_decade = trend_per_year * 10,
velocity = velocity
) %>%
dplyr::mutate(X = x, Y = y) %>%
gfplot:::utm2ll(., utm_zone = 9)
# 2008 DATASET
# hist(gf0$sd)
# mean(gf0$sd, na.rm = T)
# 0.2828
do_min_sd <- min(gf0$sd, na.rm = T)
# [1] 0.04122
max(gf0$sd, na.rm = T)
# [1] 1.531
do_median_sd <- median(gf0$sd, na.rm = T)
# [1] 0.2343
do_mean_sd <- mean(gf0$sd, na.rm = T)
do_min_sd
do_median_sd
do_mean_sd
# do_sd_choice <- do_min_sd*2
do_sd_choice <- do_median_sd
lower_change_d <- c(do_sd_choice)
upper_change_d <- c(do_sd_choice)
# variable_names <- c("temp")
# input_cell_size <- 2
## USED for truncating velocities?
max_dist_all_years <- 10 * year_range
dvocc_d <- make_vector_data(d,
variable_names = c("do_est"),
ssid = model_ssid,
# start_time = 2013,
# end_time = 2015,
# skip_time = skip_time,
input_cell_size = 2,
scale_fac = scale_fac,
min_dist = 0,
delta_t_total = years_sampled,
delta_t_step = 2,
indices = indices,
min_thresholds = lower_change_d,
max_thresholds = upper_change_d,
round_fact = 10
)
# glimpse(dvocc_d)
dvocc_d <- dvocc_d %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope)
dvocc_d <- dvocc_d %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10,
(max_dist_all_years/ year_range)*10, dvocc))
dvocc_d <- mutate(dvocc_d, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10,
-(max_dist_all_years / year_range)*10, dvocc))
dv_d <- dvocc_d %>% select(icell, x, y, dvocc) %>% unique()
# glimpse(dv_b)
# glimpse(gf0)
# gf0$dvocc <- dv_d$dvocc
gf0 <- left_join(gf0, dv_d)
lower_change_td <- c(Inf, do_sd_choice)
upper_change_td <- c(temp_sd_choice, do_sd_choice)
# variable_names <- c("temp")
# input_cell_size <- 2
## USED for truncating velocities?
max_dist_all_years <- 10 * year_range
dvocc_td <- make_vector_data(list(d, d),
variable_names = c("temp", "do_est"),
ssid = model_ssid,
# start_time = 2013,
# end_time = 2015,
# skip_time = skip_time,
input_cell_size = 2,
scale_fac = scale_fac,
min_dist = 0,
delta_t_total = years_sampled,
delta_t_step = 2,
indices = indices,
min_thresholds = lower_change_td,
max_thresholds = upper_change_td,
round_fact = 10
)
# glimpse(dvocc_d)
dvocc_td <- dvocc_td %>% mutate(dvocc = if_else(do_est.slope < 0, -km_per_decade, km_per_decade), trend = do_est.slope)
dvocc_td <- dvocc_td %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10,
(max_dist_all_years/ year_range)*10, dvocc))
dvocc_td <- mutate(dvocc_td, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10,
-(max_dist_all_years / year_range)*10, dvocc))
dv_td <- dvocc_td %>% select(icell, x, y, dvocc_both = dvocc) %>% unique()
# glimpse(dv_b)
# glimpse(gf0)
# gf0$dvocc_both <- dv_td$dvocc
gf0 <- left_join(gf0, dv_td)
#
saveRDS(gf0, file = paste0("data/_climate_vocc/DO-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-min2SD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-meanSD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
# saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-medSD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds"))
}
do_min_sd <- round(min(gf0$sd, na.rm = T), digits = 3)
do_mean_sd <- round(mean(gf0$sd, na.rm = T), digits = 3)
do_median_sd <- round(median(gf0$sd, na.rm = T), digits = 3)
print( paste("min SD =", do_min_sd))
print( paste("mean SD =", do_mean_sd))
print( paste("median SD =", do_median_sd))
do_sd_choice <- do_mean_sd
print( paste("treshold chosen =", signif(do_sd_choice, digits = 2), "ml/L"))
plot_vocc(gf0,
raster_cell_size = 4,
coast = TRUE,
vec_aes = NULL,
NA_label = "NA",
fill_col = "dvocc",
fill_label = "km/decade",
raster_alpha = 1,
na_colour = "white",
mid_fill = "honeydew", grey_water = F,
# mid_fill = "lightyellow", grey_water = F,
# mid_fill = "azure", grey_water = F,
high_fill = "gold",
low_fill = "darkcyan", # "lightseagreen",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "distance-based DO velocity",
subtitle = paste0(min(d$year), "-", max(d$year), " DO SD = ", signif(do_sd_choice, 2), " for min threshold of DO change")
)
plot_vocc(gf0,
raster_cell_size = 4,
coast = TRUE,
vec_aes = NULL,
NA_label = "NA",
fill_col = "dvocc_both",
fill_label = "km/decade",
raster_alpha = 1,
na_colour = "white",
mid_fill = "honeydew", grey_water = F,
# mid_fill = "lightyellow", grey_water = F,
# mid_fill = "azure", grey_water = F,
high_fill = "gold",
low_fill = "darkcyan", # "lightseagreen",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "distance-based temp and DO velocity"
)
BIOTIC VECTORS
gf3 <- vocc_gradient_calc(d, "est",
scale_fac = scale_fac,
indices = indices,
# divisor = 1, #default of 10 gives per decade
log_space = TRUE,
quantile_cutoff = 0.025
)
gf3 <- gf3 %>%
mutate(
CV_formula = sqrt(exp(cv^2) - 1),
CV = sd / exp(mean),
mean_biomass = exp(mean) * 10000,
## these are only relevant if divisor = 1
# trend_per_year = trend,
# trend_per_decade = trend_per_year * 10,
## if default of divisor = 10
trend_per_year = exp(trend/10), # percent change not raw
trend_per_decade = trend_per_year * 10,
velocity = velocity
) %>%
dplyr::mutate(X = x, Y = y)
# hist((gf3$trend))
# hist((gf3$sd))
b_mean_sd <-mean(gf3$sd, na.rm = T)
# 0.2828
b_min_sd <- min(gf3$sd, na.rm = T)
# [1] 0.04122
max(gf3$sd, na.rm = T)
# [1] 1.531
b_median_sd <- median(gf3$sd, na.rm = T)
# [1] 0.2343
b_choice <- 0.1 # in log space roughly equal 10% change
lower_change_t <- c(b_choice)
upper_change_t <- c(b_choice)
# variable_names <- c("temp")
# input_cell_size <- 2
## USED for truncating velocities?
max_dist_all_years <- 10 * year_range
bio5perc <- sum(gf3$mean_biomass, na.rm = TRUE) * 0.05
s <- sort(gf3$mean_biomass)
bio_sum <- cumsum(s)
lower_density_threshold <- s[which(bio_sum >= bio5perc)[1]]
dvocc_b <- make_vector_data(d,
variable_names = c("est"),
ssid = model_ssid,
# start_time = 2013,
# end_time = 2015,
# skip_time = skip_time,
input_cell_size = 2,
scale_fac = scale_fac,
min_dist = 0,
delta_t_total = years_sampled,
delta_t_step = 2,
indices = indices,
min_thresholds = lower_change_t,
max_thresholds = upper_change_t,
round_fact = 10
)
# glimpse(dvocc_b)
dvocc_b <- dvocc_b %>%
mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend=slope)
dvocc_b <- dvocc_b %>%
mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10,
(max_dist_all_years/ year_range)*10, dvocc))
dvocc_b <- dvocc_b %>%
mutate(dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10,
-(max_dist_all_years / year_range)*10, dvocc))
dv_b <- dvocc_b %>% select(icell, x, y, dvocc) %>% unique()
# glimpse(dv_b)
# glimpse(gf3)
# gf3$dvocc <- dv_b$dvocc
gf3 <- left_join(gf3, dv_b)
gf3t <- filter(gf3, mean_biomass > lower_density_threshold)
gradbv <- plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "dvocc",
fill_label = "km/decade",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
# raster_limits = c(raster_min, raster_max),
na_colour = "black",
# white_zero = TRUE,
mid_fill = "lightyellow",
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = no_trans # fourth_root_power # sqrt #
) + labs(
title = paste0(species, " (", age, ") distance-based biotic velocity"),
subtitle = paste0(min(d$year), "-", max(d$year), " threshold for change = ", b_choice, "")
)
gradbv
gf3t$squashed_velocity <- collapse_outliers(gf3t$velocity, c(0.005, 0.995))
grad3 <- gfranges::plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "squashed_velocity",
fill_label = "biotic \nvelocity",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "black",
# white_zero = TRUE,
mid_fill = "lightyellow",
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt # no_trans #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("Biotic velocity ", min(d$year), "-", max(d$year))
)
grad3
# gf3t$squashed_velocity <- collapse_outliers(gf3t$velocity, c(0.005, 0.995))
grad3 <- gfranges::plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "gradient",
fill_label = "biotic \ngradient",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "black",
# white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt # no_trans #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("Biotic gradient for final year")
)
grad3
#raster_max <- quantile(gf3$mean_biomass, 0.995, na.rm = T)
gf3$mean_biomass2 <- collapse_outliers(gf3$mean_biomass, c(0, 0.995))
grad4 <- gfranges::plot_vocc(gf3,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "mean_biomass2",
fill_label = "mean",
raster_alpha = 1,
# raster_limits = c(0, raster_max),
na_colour = "white",
# white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("Mean biomass ", min(d$year), "-", max(d$year))
)
grad4
grad5 <- gfranges::plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "trend",
fill_label = "biomass \ntrend",
raster_alpha = 1,
# raster_limits = c(0, 10),
na_colour = "black",
white_zero = TRUE,
mid_fill = "lightyellow",
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("Biomass trend ", min(d$year), "-", max(d$year))
)
grad5
grad6 <- gfranges::plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "cv",
fill_label = "sd in log space",
raster_alpha = 1,
# raster_limits = c(0, 7),
na_colour = "black",
# white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("Biotic variability (cv from vocc_gradient_calc) ", min(d$year), "-", max(d$year))
)
grad6
raster_max <- quantile(gf3t$CV, 0.90, na.rm = TRUE)
grad7 <- gfranges::plot_vocc(gf3t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "CV",
fill_label = "CV",
raster_alpha = 1,
raster_limits = c(0, raster_max),
na_colour = "black",
# white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = paste0(species, " (", age, ")"), #
subtitle = paste0("CV in estimated biomass ", min(d$year), "-", max(d$year))
)
grad7
Save VOCCs
biotic_values <- gf3 %>%
rename(
biotic_vel = velocity,
biotic_dvocc = dvocc,
biotic_trend = trend, # currently per decade
biotic_CV = CV,
sd_est = cv,
biotic_grad = gradient
) %>%
mutate(
sd_biomass = sd * 10000,
mean_biomass = exp(mean) * 10000
) %>%
select(x, y, biotic_vel, biotic_dvocc, biotic_trend, biotic_grad, biotic_CV, sd_est, sd_biomass, mean_biomass) %>%
mutate(start_year = !!start_year)
DO_values <- gf0 %>%
rename(
DO_vel = velocity,
DO_dvocc = dvocc,
DO_trend = trend_per_decade,
DO_grad = gradient,
mean_DO = mean
) %>%
select(x, y, DO_vel, DO_dvocc, DO_trend, DO_grad, mean_DO, dvocc_both)
temp_values <- gf2 %>%
rename(
temp_vel = velocity,
temp_dvocc = dvocc,
temp_trend = trend_per_decade,
temp_grad = gradient,
mean_temp = mean
) %>%
select(x, y, temp_vel, temp_dvocc, temp_trend, temp_grad, mean_temp)
fished <- readRDS("data/fishing-effort-change.rds")
vocc <- left_join(biotic_values, fished) # add fishing effort
vocc <- left_join(vocc, DO_values) # add DO
vocc <- left_join(vocc, temp_values) # add temp
vocc <- mutate(vocc, species = !!species, method = "both")
write.csv(vocc, file = paste0(
"data/_optimized_", age, "_new/", age, "-w-do-nulls-", spp, "-", ssid_string, covs, "-", min(d$year), "-", max(d$year), ".csv"
))
vocc <- read.csv(paste0(
"data/_optimized_", age, "_new/", age, "-w-do-nulls-", spp, "-", ssid_string, covs, "-", min(d$year), "-", max(d$year), ".csv"
))
vocc_t <- filter(vocc, mean_biomass > lower_density_threshold)
plot_vocc(vocc_t,
raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "temp_dvocc",
fill_label = "temperature \nkm/decade",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
# viridis_option = "C",
mid_fill = "lightyellow",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "distance-based temp velocity",
subtitle = paste0(min(d$year), "-", max(d$year), " min temp SD = ", temp_mean_sd, " for 95% min threshold of change")
)
vocc_t$squashed_temp <- collapse_outliers(vocc_t$temp_vel, c(0.005, 0.995))
plot_vocc(vocc_t, raster_cell_size = 4,
# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "squashed_temp",
fill_label = "temperature \nvelocity",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
mid_fill = "lightyellow",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "gradient-based temp velocity",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "temp_trend",
fill_label = "temperature \ntrend",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "grey",
white_zero = TRUE,
mid_fill = "lightyellow",
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = no_trans # log10 # fourth_root_power # sqrt #
) + labs(
title = "temp trend",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "mean_temp",
fill_label = "mean \ntemperature",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "grey",
white_zero = TRUE,
mid_fill = "lightyellow",
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = no_trans # log10 # fourth_root_power # sqrt #
) + labs(
title = "mean temp",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "DO_dvocc",
fill_label = "km/decade",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
# viridis_option = "C",
mid_fill = "lightyellow",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "distance-based DO and temp velocity",
subtitle = paste0(min(d$year), "-", max(d$year), " DO SD = ", signif(do_sd_choice, 2), " x 2 for threshold for do change")
)
vocc_t$squashed_do <- collapse_outliers(vocc_t$DO_vel, c(0.005, 0.995))
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "squashed_do",
fill_label = "DO \nvelocity",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
mid_fill = "lightyellow",
white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt #
) + labs(
title = "DO gradient-based velocity",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "DO_trend",
fill_label = "DO \ntrend",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
mid_fill = "lightyellow",
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = no_trans # log10 # fourth_root_power # sqrt #
) + labs(
title = "DO trend",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme = "black",
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
NA_label = "NA",
fill_col = "mean_DO",
fill_label = "mean \nDO",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
white_zero = TRUE,
high_fill = "Steel Blue 4",
low_fill = "Red 3",
mid_fill = "lightyellow",
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = no_trans # log10 # fourth_root_power # sqrt #
) + labs(
title = "mean DO",
subtitle = paste0(min(d$year), "-", max(d$year))
)
vocc_t$squashed_effort <- collapse_outliers(vocc_t$mean_effort, c(0.005, 0.995))
plot_vocc(vocc_t,
raster_cell_size = 4,# theme_black = TRUE,
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
# NA_label = "NA",
fill_col = "squashed_effort",
fill_label = "mean effort",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
na_colour = "white",
# high_fill = "Steel Blue 4",
# low_fill = "Red 3",
# white_zero = TRUE,
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = fourth_root_power #
) + labs(
title = "Fishing presssure",
subtitle = paste0(min(d$year), "-", max(d$year))
)
plot_vocc(vocc_t,
raster_cell_size = 4,# theme_black = TRUE,
coast = TRUE,
vec_aes = NULL,
# grad_vec_aes = "velocity",
# vec_lwd_range = c(0.5, 0.5),
# NA_label = "NA",
fill_col = "fishing_trend",
fill_label = "trend",
raster_alpha = 1,
# raster_limits = c(raster_min, raster_max),
# high_fill = "Steel Blue 4",
# low_fill = "Red 3",
# white_zero = TRUE,
mid_fill = "lightyellow",
vec_alpha = 0.35,
axis_lables = FALSE,
transform_col = sqrt # no_trans #
) + labs(
title = "Fishing trend",
subtitle = paste0(min(d$year), "-", max(d$year))
)
pbs-assess/gfranges documentation built on Dec. 13, 2021, 4:50 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
knitr::opts_chunk$set( fig.keep = "all", fig.show = "asis", dpi = 600, fig.width = 6, fig.height = 5, echo = FALSE, warning = FALSE, message = FALSE ) options(scipen = 999) options(digits = 4)
library("gfranges") library("dplyr") library("ggplot2") ggplot2::theme_set(gfplot::theme_pbs())
Set parameters
species <- params$species region <- params$region covs <- params$covs climate_mod <- params$climate_mod immature <- params$immature scrambled_years <- params$scrambled_years spp <- gsub(" ", "-", gsub("\\/", "-", tolower(species))) paste("species =", species) paste("region =", region) paste("covs =", covs) paste("scrambled =", scrambled_years)
Set cell sizes and scale
dist_intercept <- 0 # input_cell_size / 2 input_cell_size <- 2 scale_fac <- 2 # means that the raster is reprojected to _ X original grid (2 km)
if (region == "both odd year surveys") { survey <- c("SYN QCS", "SYN HS") model_ssid <- c(1, 3) ssid_string <- paste0(model_ssid, collapse = "n") start_year <- 2008 } if (region == "West Coast Vancouver Island") { survey <- c("SYN WCVI") model_ssid <- c(4) ssid_string <- paste0(model_ssid, collapse = "n") start_year <- 2008 } if (region == "West Coast Haida Gwaii") { survey <- c("SYN WCHG") model_ssid <- c(16) ssid_string <- paste0(model_ssid, collapse = "n") start_year <- 2008 }
rm(d_all) if (immature) { age <- "immature" d_all <- readRDS(paste0( "data/", spp, "/predictions-", spp, covs, "-1n3n4n16-imm-biomass.rds" #-prior-FALSE.rds" )) %>% filter(year >= start_year ) } else { age <- "mature" try({ d_all <- readRDS(paste0( "data/", spp, "/predictions-", spp, covs, "-1n3n4n16-total-biomass.rds" )) %>% filter(year >= start_year ) }) try({ d_all <- readRDS(paste0( "data/", spp, "/predictions-", spp, covs, "-1n3n4n16-mature-biomass-new.rds" )) %>% filter(year >= start_year ) }) } # hist(d_all$log_effort)
SET TIME FRAME
multiyear <- TRUE if (multiyear) { # unique(d_all$year) d <- d_all %>% filter(ssid %in% model_ssid) length(unique(d$year)) years_sampled <- length(unique(d$year)) year_range <- (years_sampled - 1) * 2 # # alternate version of this model # if (years_sampled == 5) indices <- c(1, 2, 2, 2, 2) # if (years_sampled == 6) indices <- c(1, 2, 2, 2, 2, 2) ## USED for truncating distance-based velocities? max_dist_all_years <- 10 * year_range # should try running with gradient based on only later years? if (years_sampled == 5) { indices <- c(1, 1, 1, 2, 2) } if (years_sampled == 6) { indices <- c(1, 1, 1, 1, 2, 2) } } else { d <- d_all %>% filter(year >= start_year) %>% filter(year <= start_year + 3) unique(d$year) indices <- c(1, 2) year_range <- 2 } skip_time <- NULL delta_t_step <- year_range delta_t_total <- year_range glimpse(d) if (scrambled_years) { d2 <- list() # set.seed(9999) # set.seed(12345) # set.seed(98765) original_years <- unique(d$year) new_years <- sample(original_years, replace = F) convert <- cbind(original_years, new_years) convert for (i in seq_len(length(original_years))) { old_year <- as.integer(convert[[i, 1]]) fake_year <- as.integer(convert[[i, 2]]) d2[[i]] <- filter(d, year == old_year) %>% mutate(year = fake_year) # browser() } d <- do.call(rbind, d2) convert }
Load climate data
nd_all <- readRDS(paste0("data/predicted-DO-", climate_mod, "-roms.rds")) %>% select(X, Y, year, temp, do_est, roms_temp) d <- left_join(d, nd_all, by = c("X", "Y", "year")) glimpse(d)
Calculate climate velocity
Temperature
rm(gf2) try({ gf2 <- readRDS(paste0("data/_climate_vocc/temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) }) if(!exists("gf2")){ gf2 <- vocc_gradient_calc(d, "temp", scale_fac = scale_fac, # all layers with max indice will be averaged for spatial gradient # if default (NULL) will use all years indices = indices, # divisor = 1, #default of 10 gives per decade quantile_cutoff = 0.05 ) gf2 <- gf2 %>% mutate( ## these are only relevant if divisor = 1 # trend_per_year = trend, # trend_per_decade = trend_per_year * 10, trend_per_year = trend/10, trend_per_decade = trend_per_year * 10, velocity = velocity ) %>% dplyr::mutate(X = x, Y = y) # FULL DATASET # mean(gf2$sd, na.rm = T) # [1] 0.3434 # min(gf2$sd, na.rm = T) # [1] 0.1266 # max(gf2$sd, na.rm = T) # [1] 1.12 # median(gf2$sd, na.rm = T) # [1] 0.287 # # 2008 DATASET # values currently depend on survey region temp_mean_sd <- mean(gf2$sd, na.rm = T) # [1] 0.3274 temp_mean_sd temp_median_sd <- median(gf2$sd, na.rm = T) # [1] 0.2882 temp_median_sd temp_min_sd <- min(gf2$sd, na.rm = T) # [1] 0.1035 temp_min_sd*2 max(gf2$sd, na.rm = T) # [1] 1.248 # temp_sd_choice <- temp_min_sd*2 # temp_sd_choice <- temp_mean_sd temp_sd_choice <- temp_median_sd # lower_change_t <- temp_sd_choice lower_change_t <- Inf upper_change_t <- temp_sd_choice rm(dvocc_t) dvocc_t <- make_vector_data(d, variable_names = c("temp"), ssid = model_ssid, input_cell_size = 2, scale_fac = scale_fac, min_dist = 0, delta_t_total = years_sampled, delta_t_step = 2, indices = indices, min_thresholds = lower_change_t, max_thresholds = upper_change_t, round_fact = 10 ) dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope) dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10, (max_dist_all_years/ year_range)*10, dvocc)) dvocc_t <- mutate(dvocc_t, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10, -(max_dist_all_years / year_range)*10, dvocc)) dv_t <- dvocc_t %>% select(icell, x, y, dvocc) %>% unique() # glimpse(dv_t) # glimpse(gf2) # gf2$dvocc <- dv_t$dvocc gf2 <- left_join(gf2, dv_t) saveRDS(gf2, file = paste0("data/_climate_vocc/temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # saveRDS(gf2, file = paste0("data/_climate_vocc/temp-dvocc-min2SD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # saveRDS(gf2, file = paste0("data/_climate_vocc/temp-dvocc-medSDasy-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) } # mean(gf2$sd, na.rm = T) # [1] 0.3274 temp_min_sd <- round(min(gf2$sd, na.rm = T), digits = 3) # # temp_mean_sd <- mean(gf2$sd, na.rm = T) # # max(gf2$sd, na.rm = T) # [1] 1.248 temp_median_sd <- round(median(gf2$sd, na.rm = T), digits = 3) # [1] 0.2882 temp_mean_sd <- round(mean(gf2$sd, na.rm = T), digits = 3) # print( paste("min SD =", round(temp_min_sd, digits = 2))) print( paste("median SD =", round(temp_median_sd, digits = 2))) print( paste("mean SD =", round(temp_mean_sd, digits = 2))) print( paste("treshold chosen =", round(temp_mean_sd, digits = 2), "degrees C"))
plot_vocc(gf2, raster_cell_size = 4, coast = TRUE, vec_aes = NULL, NA_label = "NA", fill_col = "dvocc", fill_label = "temperature \nkm/decade", raster_alpha = 1, na_colour = "white", mid_fill = "mistyrose1", grey_water = F, # mid_fill = "lavenderblush1", grey_water = F, low_fill = "royalblue4", # low_fill = "#5E4FA2", high_fill = "Red 3", raster_limits = c(0, 100), white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "distance-based temp velocity", subtitle = paste0(min(d$year), "-", max(d$year), " mean temp SD = ", signif(temp_mean_sd, digits = 2), " for min threshold of change") )
Annual temperature (ROMS)
rm(roms) try({roms <- readRDS(paste0("data/_climate_vocc/annual-temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) }) if(!exists("roms")){ roms <- vocc_gradient_calc(d, "roms_temp", scale_fac = scale_fac, # all layers with max indice will be averaged for spatial gradient # if default (NULL) will use all years indices = indices, # divisor = 1, #default of 10 gives per decade quantile_cutoff = 0.05 ) roms <- roms %>% mutate( ## these are only relevant if divisor = 1 # trend_per_year = trend, # trend_per_decade = trend_per_year * 10, trend_per_year = trend/10, trend_per_decade = trend_per_year * 10, velocity = velocity ) %>% dplyr::mutate(X = x, Y = y) # # FULL DATASET # # mean(roms$sd, na.rm = T) # [1] 0.3434 # # min(roms$sd, na.rm = T) # [1] 0.1266 # # max(roms$sd, na.rm = T) # [1] 1.12 # # median(roms$sd, na.rm = T) # [1] 0.287 # # # # 2008 DATASET # mean(roms$sd, na.rm = T) # [1] 0.3274 # temp_min_sd <- min(roms$sd, na.rm = T) # [1] 0.1035 # max(roms$sd, na.rm = T) # [1] 1.248 # temp_median_sd <- median(roms$sd, na.rm = T) # [1] 0.2882 # # temp_sd_choice <- temp_min_sd # lower_change_t <- temp_sd_choice*2 # upper_change_t <- temp_sd_choice*2 # # dvocc_t <- make_vector_data(d, # variable_names = c("temp"), # ssid = model_ssid, # input_cell_size = 2, # scale_fac = scale_fac, # min_dist = 0, # delta_t_total = years_sampled, # delta_t_step = 2, # indices = indices, # min_thresholds = lower_change_t, # max_thresholds = upper_change_t, # round_fact = 10 # ) # # dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope) # dvocc_t <- dvocc_t %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10, # (max_dist_all_years/ year_range)*10, dvocc)) # dvocc_t <- mutate(dvocc_t, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10, # -(max_dist_all_years / year_range)*10, dvocc)) # # dv_t <- dvocc_t %>% select(icell, x, y, dvocc) %>% unique() # # # glimpse(dv_t) # # glimpse(roms) # # # roms$dvocc <- dv_t$dvocc # # roms <- left_join(roms, dv_t) saveRDS(roms, file = paste0("data/_climate_vocc/annual-temp-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # mine <- gf2 %>% select(x, y, temp_trend = trend, temp_vel = velocity, temp_grad= gradient) # other <- roms %>% select(x, y, roms_temp_trend = trend, roms_temp_vel = velocity, roms_temp_grad= gradient) # # check <- left_join(mine, other) # # ggplot(check, aes(roms_temp_trend, temp_trend)) + # geom_point(alpha=0.2) + # coord_cartesian(xlim = c(-1,2.5), ylim = c(-1,2.5)) + # geom_abline(intercept = 0, slope = 1) + # gfplot::theme_pbs() }
# ggplot(check, aes(roms_temp_vel, temp_vel)) + # geom_point(alpha=0.2) + # coord_cartesian(xlim = c(-20,80), ylim = c(-20,80)) + # geom_abline(intercept = 0, slope = 1) + # gfplot::theme_pbs()
Dissolved Oxygen
rm(gf0) try({ gf0 <- readRDS(paste0("data/_climate_vocc/DO-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) }) if(!exists("gf0")){ gf0 <- vocc_gradient_calc(d, "do_est", scale_fac = scale_fac, # all layers with max indice will be averaged for spatial gradient # if default (NULL) will use all years indices = indices, # divisor = 1, #default of 10 gives per decade quantile_cutoff = 0.05 ) gf0 <- gf0 %>% mutate( ## these are only relevant if divisor = 1 # trend_per_year = trend, # trend_per_decade = trend_per_year * 10, trend_per_year = trend/10, trend_per_decade = trend_per_year * 10, velocity = velocity ) %>% dplyr::mutate(X = x, Y = y) %>% gfplot:::utm2ll(., utm_zone = 9) # 2008 DATASET # hist(gf0$sd) # mean(gf0$sd, na.rm = T) # 0.2828 do_min_sd <- min(gf0$sd, na.rm = T) # [1] 0.04122 max(gf0$sd, na.rm = T) # [1] 1.531 do_median_sd <- median(gf0$sd, na.rm = T) # [1] 0.2343 do_mean_sd <- mean(gf0$sd, na.rm = T) do_min_sd do_median_sd do_mean_sd # do_sd_choice <- do_min_sd*2 do_sd_choice <- do_median_sd lower_change_d <- c(do_sd_choice) upper_change_d <- c(do_sd_choice) # variable_names <- c("temp") # input_cell_size <- 2 ## USED for truncating velocities? max_dist_all_years <- 10 * year_range dvocc_d <- make_vector_data(d, variable_names = c("do_est"), ssid = model_ssid, # start_time = 2013, # end_time = 2015, # skip_time = skip_time, input_cell_size = 2, scale_fac = scale_fac, min_dist = 0, delta_t_total = years_sampled, delta_t_step = 2, indices = indices, min_thresholds = lower_change_d, max_thresholds = upper_change_d, round_fact = 10 ) # glimpse(dvocc_d) dvocc_d <- dvocc_d %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend = slope) dvocc_d <- dvocc_d %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10, (max_dist_all_years/ year_range)*10, dvocc)) dvocc_d <- mutate(dvocc_d, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10, -(max_dist_all_years / year_range)*10, dvocc)) dv_d <- dvocc_d %>% select(icell, x, y, dvocc) %>% unique() # glimpse(dv_b) # glimpse(gf0) # gf0$dvocc <- dv_d$dvocc gf0 <- left_join(gf0, dv_d) lower_change_td <- c(Inf, do_sd_choice) upper_change_td <- c(temp_sd_choice, do_sd_choice) # variable_names <- c("temp") # input_cell_size <- 2 ## USED for truncating velocities? max_dist_all_years <- 10 * year_range dvocc_td <- make_vector_data(list(d, d), variable_names = c("temp", "do_est"), ssid = model_ssid, # start_time = 2013, # end_time = 2015, # skip_time = skip_time, input_cell_size = 2, scale_fac = scale_fac, min_dist = 0, delta_t_total = years_sampled, delta_t_step = 2, indices = indices, min_thresholds = lower_change_td, max_thresholds = upper_change_td, round_fact = 10 ) # glimpse(dvocc_d) dvocc_td <- dvocc_td %>% mutate(dvocc = if_else(do_est.slope < 0, -km_per_decade, km_per_decade), trend = do_est.slope) dvocc_td <- dvocc_td %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10, (max_dist_all_years/ year_range)*10, dvocc)) dvocc_td <- mutate(dvocc_td, dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10, -(max_dist_all_years / year_range)*10, dvocc)) dv_td <- dvocc_td %>% select(icell, x, y, dvocc_both = dvocc) %>% unique() # glimpse(dv_b) # glimpse(gf0) # gf0$dvocc_both <- dv_td$dvocc gf0 <- left_join(gf0, dv_td) # saveRDS(gf0, file = paste0("data/_climate_vocc/DO-vocc-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-min2SD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-meanSD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) # saveRDS(gf0, file = paste0("data/_climate_vocc/DO-dvocc-medSD-", climate_mod, "-scale-", scale_fac, "-ssid", ssid_string, "-years-", years_sampled, ".rds")) } do_min_sd <- round(min(gf0$sd, na.rm = T), digits = 3) do_mean_sd <- round(mean(gf0$sd, na.rm = T), digits = 3) do_median_sd <- round(median(gf0$sd, na.rm = T), digits = 3) print( paste("min SD =", do_min_sd)) print( paste("mean SD =", do_mean_sd)) print( paste("median SD =", do_median_sd)) do_sd_choice <- do_mean_sd print( paste("treshold chosen =", signif(do_sd_choice, digits = 2), "ml/L"))
plot_vocc(gf0, raster_cell_size = 4, coast = TRUE, vec_aes = NULL, NA_label = "NA", fill_col = "dvocc", fill_label = "km/decade", raster_alpha = 1, na_colour = "white", mid_fill = "honeydew", grey_water = F, # mid_fill = "lightyellow", grey_water = F, # mid_fill = "azure", grey_water = F, high_fill = "gold", low_fill = "darkcyan", # "lightseagreen", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "distance-based DO velocity", subtitle = paste0(min(d$year), "-", max(d$year), " DO SD = ", signif(do_sd_choice, 2), " for min threshold of DO change") )
plot_vocc(gf0, raster_cell_size = 4, coast = TRUE, vec_aes = NULL, NA_label = "NA", fill_col = "dvocc_both", fill_label = "km/decade", raster_alpha = 1, na_colour = "white", mid_fill = "honeydew", grey_water = F, # mid_fill = "lightyellow", grey_water = F, # mid_fill = "azure", grey_water = F, high_fill = "gold", low_fill = "darkcyan", # "lightseagreen", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "distance-based temp and DO velocity" )
BIOTIC VECTORS
gf3 <- vocc_gradient_calc(d, "est", scale_fac = scale_fac, indices = indices, # divisor = 1, #default of 10 gives per decade log_space = TRUE, quantile_cutoff = 0.025 ) gf3 <- gf3 %>% mutate( CV_formula = sqrt(exp(cv^2) - 1), CV = sd / exp(mean), mean_biomass = exp(mean) * 10000, ## these are only relevant if divisor = 1 # trend_per_year = trend, # trend_per_decade = trend_per_year * 10, ## if default of divisor = 10 trend_per_year = exp(trend/10), # percent change not raw trend_per_decade = trend_per_year * 10, velocity = velocity ) %>% dplyr::mutate(X = x, Y = y) # hist((gf3$trend)) # hist((gf3$sd)) b_mean_sd <-mean(gf3$sd, na.rm = T) # 0.2828 b_min_sd <- min(gf3$sd, na.rm = T) # [1] 0.04122 max(gf3$sd, na.rm = T) # [1] 1.531 b_median_sd <- median(gf3$sd, na.rm = T) # [1] 0.2343 b_choice <- 0.1 # in log space roughly equal 10% change lower_change_t <- c(b_choice) upper_change_t <- c(b_choice) # variable_names <- c("temp") # input_cell_size <- 2 ## USED for truncating velocities? max_dist_all_years <- 10 * year_range bio5perc <- sum(gf3$mean_biomass, na.rm = TRUE) * 0.05 s <- sort(gf3$mean_biomass) bio_sum <- cumsum(s) lower_density_threshold <- s[which(bio_sum >= bio5perc)[1]]
dvocc_b <- make_vector_data(d, variable_names = c("est"), ssid = model_ssid, # start_time = 2013, # end_time = 2015, # skip_time = skip_time, input_cell_size = 2, scale_fac = scale_fac, min_dist = 0, delta_t_total = years_sampled, delta_t_step = 2, indices = indices, min_thresholds = lower_change_t, max_thresholds = upper_change_t, round_fact = 10 ) # glimpse(dvocc_b) dvocc_b <- dvocc_b %>% mutate(dvocc = if_else(slope < 0, -km_per_decade, km_per_decade), trend=slope) dvocc_b <- dvocc_b %>% mutate(dvocc = if_else(dvocc > (max_dist_all_years/ year_range)*10, (max_dist_all_years/ year_range)*10, dvocc)) dvocc_b <- dvocc_b %>% mutate(dvocc = if_else(dvocc < -(max_dist_all_years / year_range)*10, -(max_dist_all_years / year_range)*10, dvocc)) dv_b <- dvocc_b %>% select(icell, x, y, dvocc) %>% unique() # glimpse(dv_b) # glimpse(gf3) # gf3$dvocc <- dv_b$dvocc gf3 <- left_join(gf3, dv_b) gf3t <- filter(gf3, mean_biomass > lower_density_threshold)
gradbv <- plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "dvocc", fill_label = "km/decade", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), # raster_limits = c(raster_min, raster_max), na_colour = "black", # white_zero = TRUE, mid_fill = "lightyellow", high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = no_trans # fourth_root_power # sqrt # ) + labs( title = paste0(species, " (", age, ") distance-based biotic velocity"), subtitle = paste0(min(d$year), "-", max(d$year), " threshold for change = ", b_choice, "") ) gradbv
gf3t$squashed_velocity <- collapse_outliers(gf3t$velocity, c(0.005, 0.995)) grad3 <- gfranges::plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "squashed_velocity", fill_label = "biotic \nvelocity", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "black", # white_zero = TRUE, mid_fill = "lightyellow", high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # no_trans # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("Biotic velocity ", min(d$year), "-", max(d$year)) ) grad3
# gf3t$squashed_velocity <- collapse_outliers(gf3t$velocity, c(0.005, 0.995)) grad3 <- gfranges::plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "gradient", fill_label = "biotic \ngradient", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "black", # white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # no_trans # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("Biotic gradient for final year") ) grad3
#raster_max <- quantile(gf3$mean_biomass, 0.995, na.rm = T) gf3$mean_biomass2 <- collapse_outliers(gf3$mean_biomass, c(0, 0.995)) grad4 <- gfranges::plot_vocc(gf3, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "mean_biomass2", fill_label = "mean", raster_alpha = 1, # raster_limits = c(0, raster_max), na_colour = "white", # white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("Mean biomass ", min(d$year), "-", max(d$year)) ) grad4
grad5 <- gfranges::plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "trend", fill_label = "biomass \ntrend", raster_alpha = 1, # raster_limits = c(0, 10), na_colour = "black", white_zero = TRUE, mid_fill = "lightyellow", high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("Biomass trend ", min(d$year), "-", max(d$year)) ) grad5
grad6 <- gfranges::plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "cv", fill_label = "sd in log space", raster_alpha = 1, # raster_limits = c(0, 7), na_colour = "black", # white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("Biotic variability (cv from vocc_gradient_calc) ", min(d$year), "-", max(d$year)) ) grad6
raster_max <- quantile(gf3t$CV, 0.90, na.rm = TRUE) grad7 <- gfranges::plot_vocc(gf3t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "CV", fill_label = "CV", raster_alpha = 1, raster_limits = c(0, raster_max), na_colour = "black", # white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", axis_lables = FALSE, transform_col = sqrt # ) + labs( title = paste0(species, " (", age, ")"), # subtitle = paste0("CV in estimated biomass ", min(d$year), "-", max(d$year)) ) grad7
Save VOCCs
biotic_values <- gf3 %>% rename( biotic_vel = velocity, biotic_dvocc = dvocc, biotic_trend = trend, # currently per decade biotic_CV = CV, sd_est = cv, biotic_grad = gradient ) %>% mutate( sd_biomass = sd * 10000, mean_biomass = exp(mean) * 10000 ) %>% select(x, y, biotic_vel, biotic_dvocc, biotic_trend, biotic_grad, biotic_CV, sd_est, sd_biomass, mean_biomass) %>% mutate(start_year = !!start_year) DO_values <- gf0 %>% rename( DO_vel = velocity, DO_dvocc = dvocc, DO_trend = trend_per_decade, DO_grad = gradient, mean_DO = mean ) %>% select(x, y, DO_vel, DO_dvocc, DO_trend, DO_grad, mean_DO, dvocc_both) temp_values <- gf2 %>% rename( temp_vel = velocity, temp_dvocc = dvocc, temp_trend = trend_per_decade, temp_grad = gradient, mean_temp = mean ) %>% select(x, y, temp_vel, temp_dvocc, temp_trend, temp_grad, mean_temp) fished <- readRDS("data/fishing-effort-change.rds") vocc <- left_join(biotic_values, fished) # add fishing effort vocc <- left_join(vocc, DO_values) # add DO vocc <- left_join(vocc, temp_values) # add temp vocc <- mutate(vocc, species = !!species, method = "both") write.csv(vocc, file = paste0( "data/_optimized_", age, "_new/", age, "-w-do-nulls-", spp, "-", ssid_string, covs, "-", min(d$year), "-", max(d$year), ".csv" ))
vocc <- read.csv(paste0( "data/_optimized_", age, "_new/", age, "-w-do-nulls-", spp, "-", ssid_string, covs, "-", min(d$year), "-", max(d$year), ".csv" )) vocc_t <- filter(vocc, mean_biomass > lower_density_threshold)
plot_vocc(vocc_t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "temp_dvocc", fill_label = "temperature \nkm/decade", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", # viridis_option = "C", mid_fill = "lightyellow", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "distance-based temp velocity", subtitle = paste0(min(d$year), "-", max(d$year), " min temp SD = ", temp_mean_sd, " for 95% min threshold of change") )
vocc_t$squashed_temp <- collapse_outliers(vocc_t$temp_vel, c(0.005, 0.995)) plot_vocc(vocc_t, raster_cell_size = 4, # theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "squashed_temp", fill_label = "temperature \nvelocity", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", mid_fill = "lightyellow", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "gradient-based temp velocity", subtitle = paste0(min(d$year), "-", max(d$year)) )
plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "temp_trend", fill_label = "temperature \ntrend", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "grey", white_zero = TRUE, mid_fill = "lightyellow", vec_alpha = 0.35, axis_lables = FALSE, transform_col = no_trans # log10 # fourth_root_power # sqrt # ) + labs( title = "temp trend", subtitle = paste0(min(d$year), "-", max(d$year)) )
plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "mean_temp", fill_label = "mean \ntemperature", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "grey", white_zero = TRUE, mid_fill = "lightyellow", vec_alpha = 0.35, axis_lables = FALSE, transform_col = no_trans # log10 # fourth_root_power # sqrt # ) + labs( title = "mean temp", subtitle = paste0(min(d$year), "-", max(d$year)) )
plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "DO_dvocc", fill_label = "km/decade", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", # viridis_option = "C", mid_fill = "lightyellow", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "distance-based DO and temp velocity", subtitle = paste0(min(d$year), "-", max(d$year), " DO SD = ", signif(do_sd_choice, 2), " x 2 for threshold for do change") )
vocc_t$squashed_do <- collapse_outliers(vocc_t$DO_vel, c(0.005, 0.995)) plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "squashed_do", fill_label = "DO \nvelocity", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", mid_fill = "lightyellow", white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # ) + labs( title = "DO gradient-based velocity", subtitle = paste0(min(d$year), "-", max(d$year)) )
plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "DO_trend", fill_label = "DO \ntrend", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", mid_fill = "lightyellow", vec_alpha = 0.35, axis_lables = FALSE, transform_col = no_trans # log10 # fourth_root_power # sqrt # ) + labs( title = "DO trend", subtitle = paste0(min(d$year), "-", max(d$year)) )
plot_vocc(vocc_t, raster_cell_size = 4,# theme = "black", coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), NA_label = "NA", fill_col = "mean_DO", fill_label = "mean \nDO", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", white_zero = TRUE, high_fill = "Steel Blue 4", low_fill = "Red 3", mid_fill = "lightyellow", vec_alpha = 0.35, axis_lables = FALSE, transform_col = no_trans # log10 # fourth_root_power # sqrt # ) + labs( title = "mean DO", subtitle = paste0(min(d$year), "-", max(d$year)) )
vocc_t$squashed_effort <- collapse_outliers(vocc_t$mean_effort, c(0.005, 0.995)) plot_vocc(vocc_t, raster_cell_size = 4,# theme_black = TRUE, coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), # NA_label = "NA", fill_col = "squashed_effort", fill_label = "mean effort", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), na_colour = "white", # high_fill = "Steel Blue 4", # low_fill = "Red 3", # white_zero = TRUE, vec_alpha = 0.35, axis_lables = FALSE, transform_col = fourth_root_power # ) + labs( title = "Fishing presssure", subtitle = paste0(min(d$year), "-", max(d$year)) ) plot_vocc(vocc_t, raster_cell_size = 4,# theme_black = TRUE, coast = TRUE, vec_aes = NULL, # grad_vec_aes = "velocity", # vec_lwd_range = c(0.5, 0.5), # NA_label = "NA", fill_col = "fishing_trend", fill_label = "trend", raster_alpha = 1, # raster_limits = c(raster_min, raster_max), # high_fill = "Steel Blue 4", # low_fill = "Red 3", # white_zero = TRUE, mid_fill = "lightyellow", vec_alpha = 0.35, axis_lables = FALSE, transform_col = sqrt # no_trans # ) + labs( title = "Fishing trend", subtitle = paste0(min(d$year), "-", max(d$year)) )
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.