Species in Atlantis for Rpath

The full species list from Atlantis, with group type and whether it is fished.

#localatlantisdir <- "/Users/sarah.gaichas/Documents/0_Data/ms-keyrun/simulated-data/atlantisoutput/NOBA_sacc_38"
localatlantisdir <- here("data-raw")

fwspp <- atlantisom::load_fgs(localatlantisdir, "nordic_groups_v04.csv") %>% 
  dplyr::filter(IsTurnedOn == 1) %>%
  dplyr::select(Code, Name, Long.Name, NumCohorts, isFished, InvertType) 

knitr::kable(fwspp)
Code Name Long.Name NumCohorts isFished InvertType
POB Polar_bear Polar Bear 10 0 MAMMAL
KWH Killer_whale Killer whale 10 0 MAMMAL
SWH Sperm_whale Sperm whale 8 0 MAMMAL
HWH Humpb_whale Humpback whale 10 0 MAMMAL
MWH Minke_whale Minke whale 10 1 MAMMAL
FWH Fin_whale Fin whale 10 0 MAMMAL
BES Beard_seal Bearded Seal 10 0 MAMMAL
HAS Harp_seal Harp Seal 10 0 MAMMAL
HOS Hood_seal Hooded Seal 10 0 MAMMAL
RIS Ring_seal Ringed Seal 10 0 MAMMAL
SBA Sea_b_arct Arctic seabirds 10 0 BIRD
SBB Sea_b_bor Boreal seabirds 10 0 BIRD
SHO Sharks_other Other sharks 10 0 SHARK
DEO Demersals_other Other demersals 10 1 FISH
PEL Pelagic_large Large pelagic fish 10 1 FISH
PES Pelagic_small Small pelagic fish 10 1 FISH
REO Redfish_other Other redfish 10 1 FISH
DEL Demersal_large Large demersals 10 1 FISH
FLA Flatfish_other Other flatfish 10 1 FISH
LRD Long_rough_dab Long rough dab 10 1 FISH
SSK Skates_rays Skates and rays 10 0 FISH
MES Mesop_fish Mesopelagic fish 10 1 FISH
GRH Green_halibut Greenland halibut 10 1 FISH
MAC Mackerel Mackerel 10 1 FISH
HAD Haddock Haddock 10 1 FISH
SAI Saithe Saithe 10 1 FISH
RED Redfish Redfish 10 1 FISH
BWH Blue_whiting Blue whiting 10 1 FISH
SSH Norwegian_ssh Norwegian spring spawning herring 10 1 FISH
NCO North_atl_cod Northeast Atlantic cod 10 1 FISH
PCO Polar_cod Polar cod 10 1 FISH
CAP Capelin Capelin 5 1 FISH
PWN Prawns Prawns 2 1 PWN
CEP Squid Squid 2 0 CEP
KCR King_crab Red king crab 1 1 LG_INF
SCR Snow_crab Snow crab 6 1 FISH
ZG Gel_zooplankton Gelatinous zooplankton 1 0 LG_ZOO
ZL Large_zooplankton Large zooplankton 1 0 LG_ZOO
ZM Medium_zooplankton Medium zooplankton 1 0 MED_ZOO
ZS Small_zooplankton Small zooplankton 1 0 SM_ZOO
DF Dinoflag Dinoflagellates 1 0 DINOFLAG
PS Small_phytop Small phytoplankton 1 0 SM_PHY
PL Large_phytop Large phytoplankton 1 0 LG_PHY
BC Predatory_ben Predatory benthos 1 0 LG_INF
BD Detrivore_ben Detrivore benthos 1 0 LG_INF
BFF Benthic_filterf Benthic filter feeders 1 0 SED_EP_FF
SPO Sponges Sponges 1 0 SED_EP_FF
COR Corals Corals 1 0 SED_EP_FF
PB Pel_bact Pelagic bacteria 1 0 PL_BACT
BB Ben_bact Sediment bacteria 1 0 SED_BACT
DR Ref_det Refractory detritus 1 0 REF_DET
DL Lab_det Labile detritus 1 0 LAB_DET
DC Carrion Carrion 1 0 CARRION

Total Catch

Already in mskeyrun (sarah_wgsamsim branch)? Yes.

Plot it. Turns out even though mammals can be caught, they are not in the catch output so probably not in this run. The catch we have added here in addition to the main mskeyrun data is Prawns, Redfish other, and Snow crabs (which appear to have 0 catch).

# make sure to install mskeyrun branch sarah_wgsamsim as not yet in main branch

# make a catch series function that can be split by fleet? this doesnt
# also note different time (days) from model timestep in all other output
plotC <- function(dat, truedat=NULL){
  
    ctbio <- dat %>% filter(variable=="catch")
    ctcv <- dat %>% filter(variable=="cv")
  
    ggplot() +
    geom_line(data=ctbio, aes(x=year,y=value, color="Catch biomass"), 
              alpha = 10/10) +
    {if(!is.null(truedat)) geom_line(data=truedat, aes(x=time/365,y=atoutput, color="True Catch"), alpha = 3/10)} + 
    theme_tufte() +
    theme(legend.position = "top") +
    xlab("model year") +
    ylab("tons") +
    labs(colour=dat$ModSim) +
    facet_wrap(~Name, scales="free") 
  
}

catchbio_ss <- mskeyrun::simCatchIndexFW #atlantisom::read_savedfisheries(d.name, 'Catch')

Fishery catch time series

# observed catch only
plotC(catchbio_ss)

Fishery catch subannual

catchbio_sub <- mskeyrun::simCatchIndexSubannualFW

# observed catch only
plotC(catchbio_sub)

Total Biomass

Update: working now for all groups.

Available for age structured fish groups, other age structured groups, and biomass pools:

# plot biomass time series facet wrapped by species
plotB <- function(dat, truedat=NULL){
  
    svbio <- dat %>% filter(variable=="biomass")
    svcv <- dat %>% filter(variable=="cv")
  
    ggplot() +
    geom_line(data=svbio, aes(x=year,y=value, color="Survey Biomass"), 
              alpha = 10/10) +
    {if(!is.null(truedat)) geom_line(data=truedat, aes(x=time/365,y=atoutput, color="True B"), alpha = 3/10)} + 
    theme_tufte() +
    theme(legend.position = "top") +
    xlab("model year") +
    ylab("tons") +
    labs(colour=dat$ModSim) +
    facet_wrap(~Name, scales="free") 
  
}

survObsBiom <- mskeyrun::simSurveyIndexFW #atlantisom::read_savedsurvs(d.name, 'survB')

Survey biomass index

# compare with true output (all timesteps)
# for(s in names(survObsBiom)){
#   cat("  \n##### ",  s,"  \n")
#   print(plotB(survObsBiom[[s]][[1]], omlist_ss$truetotbio_ss))
#   cat("  \n")
# }

# plots survey only
 for(s in unique(survObsBiom$survey)){
   cat("  \n#### ",  s,"  \n")
   print(plotB(survObsBiom %>%
                 filter((survey %in% s))))
   cat("  \n")
 }

BTS_fall_allbox_effic1_fw

BTS_fall_allbox_effic1

BTS_spring_allbox_effic1_fw

BTS_spring_allbox_effic1

Biomass pools added to workflow inatlantisom::run_truth, atlantisom::om_species, and atlantisom::om_index, new surveys generated in mskeyrun SimData vignette, new datasets generated with mskeyrun::create_sim_survey_index_fw.

Diets for all

Did I pull this? Just now yes.

How many prey categories across all diets?

preynames <- unique(c(mskeyrun::simSurveyDietcomp$prey, mskeyrun::simSurveyDietcompFW$prey))

preynames
##  [1] "Large_phytop"       "Large_zooplankton"  "Medium_zooplankton"
##  [4] "Capelin"            "Mesop_fish"         "Norwegian_ssh"     
##  [7] "Redfish"            "Squid"              "Prawns"            
## [10] "Small_zooplankton"  "Carrion"            "Detrivore_ben"     
## [13] "Gel_zooplankton"    "Skates_rays"        "Blue_whiting"      
## [16] "Mackerel"           "North_atl_cod"      "Redfish_other"     
## [19] "Ben_bact"           "Benthic_filterf"    "Lab_det"           
## [22] "Long_rough_dab"     "Pel_bact"           "Predatory_ben"     
## [25] "Haddock"            "King_crab"          "Polar_cod"         
## [28] "Flatfish_other"     "Pelagic_small"      "Demersal_large"    
## [31] "Green_halibut"      "Dinoflag"           "Demersals_other"   
## [34] "Small_phytop"       "Saithe"             "Ref_det"           
## [37] "Beard_seal"         "Harp_seal"          "Hood_seal"         
## [40] "Minke_whale"        "Pelagic_large"      "Ring_seal"         
## [43] "Fin_whale"          "Sea_b_arct"         "Sea_b_bor"         
## [46] "Sharks_other"

Plotting functions and colors for everyone

#  #http://medialab.github.io/iwanthue/
# 46 colors hard force vector colorblind friendly,sort by diff
preycol <- c("#647c00", 
             "#ff91cb", 
             "#114000", 
             "#6499ff",
             "#ffca82",
             "#1e1c5e",
             "#a2c032",
             "#820038",
             "#00b06a",
             "#4e005c",
             "#ffa83c",
             "#772c5c",
             "#9f4f00",
             "#78adff",
             "#ff6c80",
             "#003376",
             "#8e6a00",
             "#005ccb",
             "#00703a",
             "#a883ff",
             "#008c59",
             "#5f2192",
             "#b0e37f",
             "#6a0064",
             "#a7e66e",
             "#6a0022",
             "#68ae2f",
             "#d12e8c",
             "#7a9852",
             "#ca0a5c",
             "#e0d688",
             "#7d007a",
             "#3d9c28",
             "#9e5487",
             "#af9700",
             "#b391d6",
             "#005506",
             "#ae5acd",
             "#009639",
             "#0163b8",
             "#ff9065",
             "#fa98ff",
             "#692c00",
             "#ff8aa7",
             "#8f0012",
             "#b02713")
names(preycol) <- as.factor(preynames)

# going for more greyscale for unident categories, same website
unidcol <- c("#b8b8b2",
             "#302a1d",
             "#6b7069",
             "#1e3430")
names(unidcol) <- as.factor(c("Unid", "Unid_Fish",  "Unid_Invert", "Unid_Plankton"))

col <- c(preycol, unidcol)

# plot diet comp over time at age by species
plotdiet <- function(dat, compdat=NULL, namedat, namecomp=NULL){
  
  dat <- dat %>% add_column(run = namedat)
  if(!is.null(compdat)) compdat <- compdat %>% add_column(run = namecomp)
  
    ggplot() +
    geom_bar(data=dat, aes(year, value, fill=prey), stat = "identity") +
    {if(!is.null(compdat)) geom_bar(data=compdat, aes(year, value, fill=prey), stat = "identity")} + 
    theme_tufte() +
    theme(legend.position = "bottom") +
    xlab("year") +
    ylab("diet proportion") +
    facet_grid(agecl~run) + 
    scale_fill_manual(values=col) + 
    ggtitle(dat$Name)
  
}

# method for a single species diet, no comparisons
# plist = lapply(split(ms_diet, ms_diet$species), function(d) {
#   ggplot(d, aes(time.days/365, atoutput, fill=prey)) + 
#     geom_bar(stat = "identity") +
#     facet_wrap(species~agecl) +
#     xlab("year") +
#     ylab("diet proportion") +
#     theme_tufte() +
#     theme(legend.position="bottom")
# })

Visualize diet comparisons

Surveyed diet 11 species

Here we compare the surveyed diet comps:

preds <- unique(mskeyrun::simSurveyDietcomp$Name)

surveys <- unique(mskeyrun::simSurveyDietcomp$survey)

for(i in 1:length(preds)) {
  cat("  \n####",  as.character(preds[i]),"  \n")
  print(plotdiet(dat = mskeyrun::simSurveyDietcomp |> dplyr::filter(survey == surveys[1], Name %in% preds[i]), 
                 namedat = surveys[1], 
                 compdat = mskeyrun::simSurveyDietcomp |> dplyr::filter(survey == surveys[2], Name %in% preds[i]),
                 namecomp = surveys[2])) 
  cat("  \n")
}

Blue_whiting

Capelin

Green_halibut

Haddock

Long_rough_dab

Mackerel

North_atl_cod

Norwegian_ssh

Polar_cod

Redfish

Saithe

Surveyed diet remaining species

surveys <- unique(mskeyrun::simSurveyDietcompFW$survey)[c(1,3)] # has other surveys in it
preds <- mskeyrun::simSurveyDietcompFW |>
  dplyr::filter(survey %in% surveys) |>
  dplyr::select(Name) |>
  dplyr::distinct() |>
  as.vector() |>
  unname() |>
  unlist()


for(i in 1:length(preds)) {
  cat("  \n####",  as.character(preds[i]),"  \n")
  print(plotdiet(dat = mskeyrun::simSurveyDietcompFW |> dplyr::filter(survey == surveys[1], Name %in% preds[i]), 
                 namedat = surveys[1], 
                 compdat = mskeyrun::simSurveyDietcompFW |> dplyr::filter(survey == surveys[2], Name %in% preds[i]),
                 namecomp = surveys[2])) 
  cat("  \n")
}

Beard_seal

Benthic_filterf

Corals

Demersal_large

Demersals_other

Detrivore_ben

Fin_whale

Flatfish_other

Gel_zooplankton

Harp_seal

Hood_seal

Humpb_whale

Killer_whale

King_crab

Large_zooplankton

Medium_zooplankton

Mesop_fish

Minke_whale

Pelagic_large

Pelagic_small

Polar_bear

Prawns

Predatory_ben

Redfish_other

Ring_seal

Sea_b_arct

Sea_b_bor

Sharks_other

Skates_rays

Small_zooplankton

Sperm_whale

Sponges

Squid

Total Production

Can we get total production from Atlantis? What is in PROD.nc

Or do we need to sum the catch and consumption removals with the population growth in each year to get production

Total Consumption

Already have this from detailed diet check