Merge branch '55-add-grain-dists-to-sim-result-handling' into 'master'

Resolve "add grain dists to sim result handling"

Closes #55

See merge request LDalby/almass!34

Grendel/ALMaSS_q28_goose_scenario_01

@@ -5,7 +5,7 @@
 #SBATCH --mem=240G
 #SBATCH --ntasks-per-node=1
 #SBATCH --cpus-per-task=4
-#SBATCH --time=48:00:00
+#SBATCH --time=168:00:00
 #SBATCH --exclusive
 #SBATCH --mem-per-cpu=60G
 #SBATCH --mail-type=ALL

Grendel/scenario_01_goose-numbers_file-dist.R

@@ -22,7 +22,7 @@ basedir <- "/home/ldalby/workspace/Goose/RunDirectory"
 pathtodirs <- path('/scratch', args)
 # Setup the directories
 npar <- 4  # Specifies the number of run directories
-years <- 10  # the number of seasons to run (goose sims run over the year boundary)
+years <- 100  # the number of seasons to run (goose sims run over the year boundary)
 basename <- "WD"  # The prefix to the directories
 # Make the directories and copy the files:
 for (i in 1:npar) {

R/01_plot-field-map.R

@@ -49,8 +49,8 @@ bckgr <- readRDS(here("data", "vejlerne-background.rds"))
 fields <- st_read(here("data", "vejlerne-fields.gpkg"))
 
 # Make base grid ----
-sf::st_make_grid(outline, n = c(36000/400, 24000/400))  %>% 
-  sf::st_sf(id = 1:5400, .) -> grd
+sf::st_make_grid(outline, n = c(36000/500, 24000/500))  %>% 
+  sf::st_sf(id = 1:3456, .) -> grd
 
 st_centroid(fields) -> field_centroids 
 
@@ -59,6 +59,7 @@ grd %>%
   drop_na() -> grd_id_polyref
 
 # Read results from scenario ----
+date_stamp <-  "2018-12-29"
 date_stamp <-  "2018-11-12"
 current <- path_home("Desktop", "goose", date_stamp)
 dir_ls(current, glob = "*.out") %>% 
@@ -70,39 +71,46 @@ res_file <- glue("scenario_results_{jobs}_{date_stamp}.fst")
 path_desktop("Goose", date_stamp, jobs, res_file) %>% 
   fst::read_fst() %>% 
   as_tibble() -> goose_numbers
-path_desktop("Goose", date_stamp, jobs, "scenario-grain-distributions.txt_1521869_2018-11-12.fst") %>% 
+path_desktop("Goose", date_stamp, jobs, "scenario-grain-distributions.txt_1521869_2018-11-12.fst") %>%
+# path_desktop("Goose", date_stamp, jobs, "scenario-grain-distributions.txt_2057427_2018-12-29.fst") %>% 
   fst::read_fst() %>% 
   as_tibble() -> grain_dists
 
+grain_dists %>% group_by(param, value) %>% count(grain_dist)
+
 left_join(goose_numbers, grain_dists, by = c("param", "value", "season")) -> goose_numbers
+goose_numbers %>% 
+  mutate(replicate = 1L) -> goose_numbers
 
+goose_numbers <- bind_rows(goose_numbers, goose_numbers_2)
 # Summarize field forage data ----  
 # Mean numbers per day
 goose_numbers %>% 
   mutate(day = day - 365*season) %>% 
-  select(param, value, season, day, pinkfoot, barnacle, greylag) %>% 
-  group_by(param, value, season, day) %>% 
+  select(replicate, grain_dist, param, value, season, day, pinkfoot, barnacle, greylag) %>% 
+  group_by(replicate, grain_dist, param, value, season, day) %>% 
   summarise_at(.vars = vars(pinkfoot, barnacle, greylag),
                .funs = funs(daily = sum)) %>% 
-  group_by(param, value, day) %>% 
+  group_by(grain_dist, param, value, day) %>% 
   summarise_at(.vars = vars(ends_with("daily")),
                .funs = funs(avg = mean,
                             min, max)) %>% 
   ungroup() %>% 
-  gather(key = species, value = numbers, -param, -value, -day) %>% 
+  gather(key = species, value = numbers, -grain_dist, -param, -value, -day) %>% 
   filter(numbers > 0) -> numbers_per_day
 
 numbers_per_day %>% 
   filter(value == 1) -> default
 left_join(numbers_per_day, default, 
-          by = c("param",  "day", "species"),
+          by = c("grain_dist", "param",  "day", "species"),
           suffix = c("", "_default")) %>% 
   mutate(standarized = numbers/numbers_default) -> numbers_per_day
 
-
+goose_sp <- "pinkfoot_daily_avg"
+goose_sp <- "greylag_daily_avg"
 numbers_per_day %>% 
   filter(param == "GOOSE_BN_STARTNO_SCALER",
-         species == "pinkfoot_daily_avg") %>% 
+         species == goose_sp) %>% 
   ggplot(aes(day, standarized)) +
   geom_line(aes(color = factor(value))) + 
   scale_color_viridis_d(name = "Barnacle scaler") +
@@ -111,20 +119,20 @@ numbers_per_day %>%
        subtitle = "Pinkfeet numbers as a function of barnacle numbers",
        caption = "Numbers are average of 10 consecutive seasons",
        y = "proportion of default numbers") +
-  facet_grid(~grain_dists)
+  facet_grid(~grain_dist)
   
 
-
 goose_numbers %>% 
-  select(param, value, season, pinkfoot, barnacle, greylag) %>% 
-  group_by(param, value, season) %>% 
+  select(replicate, grain_dist, param, value, season, pinkfoot, barnacle, greylag) %>% 
+  group_by(replicate, grain_dist, param, value, season) %>% 
   summarise(pinkfoot = sum(pinkfoot, na.rm = TRUE),
             barnacle = sum(barnacle, na.rm = TRUE),
             greylag = sum(greylag, na.rm = TRUE)) %>% 
-  ungroup() %>% 
-  gather(key = species, value = numbers, -param, -value, -season) -> numbers_per_season
-
-
+  group_by(grain_dist, param, value, season) %>%  
+  summarise(pinkfoot = mean(pinkfoot, na.rm = TRUE),
+            barnacle = mean(barnacle, na.rm = TRUE),
+            greylag = mean(greylag, na.rm = TRUE)) %>% 
+  gather(key = species, value = numbers, -grain_dist, -param, -value, -season) -> numbers_per_season
 
 
 numbers_per_season %>% 
@@ -132,7 +140,7 @@ numbers_per_season %>%
   ggplot(aes(value, numbers/1e6)) +
   geom_point(alpha = 0.5) +
   geom_smooth() + 
-  facet_wrap(~species, scales = "free_y") + 
+  facet_grid(species ~ grain_dist, scales = "free_y") + 
   hrbrthemes::theme_ipsum_rc(axis_title_size = 12) + 
   labs(title = "Goose days as a function of barnacle goose input numbers",
        y = expression("Goose days" %*% 1e-6),
@@ -142,16 +150,16 @@ numbers_per_season %>%
 
 # Average goose days per season
 goose_numbers %>% 
-  select(param, value, season, pinkfoot, barnacle, greylag, polyref) %>% 
-  group_by(param, value, season, polyref) %>% 
+  select(grain_dist, param, value, season, pinkfoot, barnacle, greylag, polyref) %>% 
+  group_by(grain_dist, param, value, season, polyref) %>% 
   summarise(pinkfoot = sum(pinkfoot, na.rm = TRUE),
             barnacle = sum(barnacle, na.rm = TRUE),
             greylag = sum(greylag, na.rm = TRUE)) %>% 
-  group_by(param, value, polyref) %>% 
+  group_by(grain_dist, param, value, polyref) %>% 
   summarise(pinkfoot = mean(pinkfoot, na.rm = TRUE),
             barnacle = mean(barnacle, na.rm = TRUE),
             greylag = mean(greylag, na.rm = TRUE)) %>% 
-  gather(key = species, value = numbers, -param, -value, -polyref) %>% 
+  gather(key = species, value = numbers, -grain_dist, -param, -value, -polyref) %>% 
   dplyr::filter(numbers > 0) -> season_goose_numbers
 
 
@@ -160,7 +168,7 @@ goose_numbers %>%
 # the same cell.
 left_join(grd_id_polyref, season_goose_numbers, by = "polyref") %>% 
   drop_na() %>% 
-  group_by(id, param, value, species) %>% 
+  group_by(id, grain_dist, param, value, species) %>% 
   summarise(goose_days = sum(numbers)) -> gdays
 
 readRDS(path_desktop("gdays.RDS")) -> gdays
@@ -169,7 +177,8 @@ readRDS(path_desktop("gdays.RDS")) -> gdays
 gdays %>% 
   ungroup() %>% 
   filter(param == "GOOSE_BN_STARTNO_SCALER",
-         species == "barnacle") %>% 
+         species == "barnacle",
+         grain_dist == "high") %>% 
   select(-param, -species) -> gdays_bn_scaler 
 
 ggmap(bckgr) + 
@@ -184,12 +193,13 @@ ggmap(bckgr) +
        y = "") +
   transition_manual(value) -> bn_map
 animate(bn_map)
-anim_save(filename = path_desktop("bn.gif"))
+anim_save(filename = path_desktop("bn-high.gif"))
 # Pinkfoot
 gdays %>% 
   ungroup() %>% 
   filter(param == "GOOSE_BN_STARTNO_SCALER",
-         species == "pinkfoot") %>% 
+         species == "pinkfoot",
+         grain_dist == "low") %>% 
   select(-param, -species) -> gdays_bn_scaler 
 
 ggmap(bckgr) + 
@@ -204,7 +214,7 @@ ggmap(bckgr) +
        caption = "Pinkfoot") +
   transition_manual(value) -> pf_map
 animate(pf_map)
-anim_save(filename = path_desktop("pf.gif"))
+anim_save(filename = path_desktop("pf-low.gif"))
 
 gdays_bn_scaler %>% 
   st_set_geometry(NULL) %>% 
@@ -244,7 +254,7 @@ ggmap(bckgr) +
 
 # Depletion
 forage %>% 
-  select(season, day, grain, maize) %>% 
+  select(grain_dist, season, day, grain, maize) %>% 
   group_by(season, day) %>% 
   summarise(grain = sum(grain),
             maize = sum(maize)) %>%