R/noise_simulations.R

In erikbsolbu/dynamicSAD:

# low-high-high ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 0, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 6, # species heterogeneity
                           sigma_s_sq = 3, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_lhh <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# high-high-high ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 2, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 6, # species heterogeneity
                           sigma_s_sq = 3, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_hhh <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# low-high-low ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 0, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 6, # species heterogeneity
                           sigma_s_sq = 0, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_lhl <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# low-low-high ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 0, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 0, # species heterogeneity
                           sigma_s_sq = 3, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_llh <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# low-low-low ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 0, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 0, # species heterogeneity
                           sigma_s_sq = 0, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_lll <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# high-low-low ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 2, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 0, # species heterogeneity
                           sigma_s_sq = 0, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_hll <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# high-high-low ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 2, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 6, # species heterogeneity
                           sigma_s_sq = 0, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_hhl <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))

# high-low-high ####
fish_boot_mt_input <- list(n_sp = 50, # number of species
                           n_year =  seq(0, 10, by = 1), # number of time points
                           n_loc = 3, # number of locations
                           n_obs = 1, # number of locations
                           mu = -0.5, # mean log abundance
                           # mu_sigma = fish_point_est_wide$sigma_mu_sq_p[2], # mean log abundance
                           mu_sigma = 2, # mean log abundance
                           # mu_sigma = 2, # mean log abundance
                           sigma_r_sq = 0, # species heterogeneity
                           sigma_s_sq = 3, # environmental variance
                           # sigma_s_sq = 0, # environmental variance
                           sigma_d_sq = 0, # sampling variance
                           delta = 0.2, # "strength of density regulation"
                           eta = 0,
                           pos_x = NULL,
                           pos_y = NULL) # 1 / "spatial scaling of noise"

fish_data_mt_boot <- f_sim_sad_fast(input = fish_boot_mt_input, time = "fixed", dependency = "time", var_mu = TRUE)

inter_fish_boot <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  dplyr::select(year_m, pos_m) %>%
  arrange(pos_m, year_m) %>%
  distinct() %>%
  mutate(inter = interaction(year_m, pos_m)) %>%
  group_by(pos_m) %>%
  mutate(min_int = inter[1],
         max_int = inter[length(inter)]) %>%
  ungroup() %>%
  dplyr::select(min_int, max_int) %>%
  distinct()

plot_hlh <- fish_data_mt_boot %>%
  filter(abundance > 0) %>%
  mutate(pos_m = rank(pos_m)) %>%
  ggplot(aes(x = interaction(year_m, pos_m), y = species)) +
  geom_raster(aes(fill = log(abundance)), hjust = 1) +
  geom_vline(data = inter_fish_boot, aes(xintercept = min_int), colour = "red") +
  scale_fill_viridis_c("log ab.") +
  theme(axis.text.x = element_text(hjust = 1, angle = 45))


library(cowplot)

plot_grid(plot_hhh, plot_hlh, plot_hhl, plot_hll,
          # labels = c("hhh", "hlh",
          #            "hhl", "hll"),
          labels = c("A", "B",
                     "C", "D"))

ggsave("logab_high_common_noise.pdf", width = 6, height = 9)

plot_grid(plot_hhh, plot_hlh,
          plot_lhh, plot_llh,
          # labels = c("hhh", "hlh",
          #            "lhh", "llh")
          labels = c("A", "B",
                     "C", "D"))

ggsave("logab_high_environmental_noise.pdf", width = 6, height = 9)