In noaa-garfo/discaRd: Cochran bycatch estimation

knitr::opts_chunk$set(echo=FALSE, warning = FALSE, 
                                            message = FALSE, cache = FALSE,
                                            progress = TRUE, verbose = FALSE, comment = F
                                            , error = FALSE, dev = 'png', dpi = 200)

options(knitr.kable.NA = '')

Background

Exploration of discard estimates from november year module

Steps

pull all merged trips from CAMS_OBS_CATCH
Set stratification variables for all trips

# FULL Stratification variables

stratvars = c('SPECIES_STOCK'
              ,'CAMS_GEAR_GROUP'
                            , 'MESHGROUP'
                          , 'TRIPCATEGORY'
                          , 'ACCESSAREA')

- SPECIES_STOCK is taken from CAMS support table `MAPS.CAMS_STATAREA_STOCK`
- CAMS_GEAR_GROUP is derived from a support table (`MAPS.CAMS_GEARCODE_STRATA`)
- MESHGROUP is hardcoded for all trips according to decisions made by the mesh subgroup (see summary when available)
- SECTID comes from a CAMS matching table (`MAPS.MATCH_MULT_SECTID`)

Run discaRd
- there are sub-passes for year t and year t-1

# Assumed Stratification variables

stratvars = c('SPECIES_STOCK'
              , 'CAMS_GEAR_GROUP',
              ,'MESHGROUP')

The discaRd functions allow for an assumed rate to be calculated. This assumed rate is relative to the stratification used in the functions. Here, the stratification is coarsened to 'SPECIES_STOCK', 'CAMS_GEAR_GROUP'and 'MESHGROUP'.
A transition rate is calculated between year t and year t-1. This rate determines how much, if any, information is used from previous years.
A broad stock/gear stratification is applied to subtrips that do not meet the requirements for the full or assumed stratification. This broad stock/gear stratification utilizes observed coverage from both the focal year and the previous year as it includes gear types with low observer coverage.
Rates and DISCARD_SOURCE (in parentheses) are assigned for each trip according to:
(I) in season rate; >= 5 trips in Full Stratification
(T) Transition in season rate; < 5 trips in Full Stratification, year t, AND >= 5 trips in year t-1
(GM) Assumed rate. This is the rate when there were <5 trips in season and <5 in previous season.
(G) Broad stock/gear rate <5 trips in season and <5 trips in previous season
(N) No coverage in focal season and previous season; discard rate is not possible to estimate.
(O) Observed values used from observed rate trips; discard rate is NOT USED.
CV calculations are available for (I), (T), (GM) and (G).
Discard pounds per trip are calculated according to

    mutate(DISC_MORT_RATIO = coalesce(DISC_MORT_RATIO, 1)) %>%
    mutate(DISCARD = case_when(!is.na(LINK1) ~ DISC_MORT_RATIO*OBS_DISCARD
                                                         , is.na(LINK1) ~ DISC_MORT_RATIO*COAL_RATE*LIVE_POUNDS)

            # COAL_RATE is the final discard rate used. It is 'coalesced' from the (I), (A) and (G) rates

By assigning SPECIES_STOCK as a stratification variable, the computation time is reduced. Each subtrip may only occur in a single statistical area so it should never cross stock boundaries.

Once the full table (CAMS_OBS_CATCH) is loaded, each species takes ~X seconds to process on the server.

Output tables are produced for each species. These can easily be recombined. Preliminary species tables have been shared on MAPS and CAMS_GARFO.

    MAPS.CAMS_DISCARD_GOLDENTILEFISH_19
    CAMS_GARFO.CAMS_DISCARD_GOLDENTILEFISH_19

    MAPS.CAMS_DISCARD_BLUELINETILEFISH_19
    CAMS_GARFO.CAMS_DISCARD_BLUELINETILEFISH_19

Diagnostic Plots/Tables

# library(odbc)
# library(dplyr, warn.conflicts = FALSE)
# # library(dbplyr)
library(ggplot2)
# # library(config)
# library(stringr)
# library(discaRd)
# library(knitr)
# library(kableExtra)
devtools::load_all()
options(scipen = 999)

discard <- read.csv("~/discaRd/CAMS/MODULES/NOVEMBER/November_comparisons.csv")

'%ni%' <- Negate("%in%")

discard %>% filter(Source == 'Federal') %>%  pivot_wider(id_cols = !Source, names_from = Method,values_from =  c(Discards, kall)) %>% mutate(diff = Discards_CAMS-Discards_QM, perc_diff = 100 * (Discards_QM - Discards_CAMS)/(Discards_CAMS)) %>% dplyr::select(Species, 
Discards_CAMS, Discards_QM, diff, perc_diff) %>%dplyr::mutate(across(where(is.numeric), round, 2))%>% 
        DT::datatable(caption = 'Comparison of total federal discards in FY 2019 by species and source')

discard <- read.csv("~/discaRd/CAMS/MODULES/NOVEMBER/November_comparisons.csv")

'%ni%' <- Negate("%in%")

discard %>% filter(Source %in% c("State")) %>% dplyr::select(Species, Discards) %>%
    DT::datatable(caption = 'Total discards by species for state trips')

 dw_maps <- config::get(config = "maps", file = "~/config_group.yml")

# Connect to database - move this to config file in the future - quick addition for server
  connectString <- paste(
    "(DESCRIPTION=",
    "(ADDRESS=(PROTOCOL=tcp)(HOST=", dw_maps$host, ")(PORT=", dw_maps$port, "))",
    "(CONNECT_DATA=(SERVICE_NAME=",dw_maps$svc, ")))",
    sep = ""
  )

# Connect to oracle each loop in case of timeouts
    bcon <- ROracle::dbConnect(
      drv = ROracle::Oracle(),
      username = dw_maps$uid,
      password = dw_maps$pwd,
      dbname = connectString
    )

db_goldentilefish = tbl(bcon, sql("select * from MAPS.CAMS_DISCARD_TILEFISHGOLDEN_2019")) %>% 
    collect

db_goldentilefish$SPECIES <- "GOLDENTILEFISH"

db_bluelinetilefish = tbl(bcon, sql("select * from MAPS.CAMS_DISCARD_TILEFISHBLUELINE_2019")) %>% 
    collect

db_bluelinetilefish$SPECIES <- "BLUELINETILEFISH"

db_final <- rbind(db_goldentilefish, db_bluelinetilefish)

strat_goldentilefish <- db_goldentilefish %>%   
   dplyr::group_by(SPECIES, STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                   dplyr::select(
                              SPECIES
                              ,STRATA_FULL
                             # STRATA 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD
              ,DISCARD_SOURCE
              ,FED_OR_STATE)

strat_bluelinetilefish <- db_bluelinetilefish %>%   
    dplyr::group_by(SPECIES, STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                   dplyr::select(
                              SPECIES
                              ,STRATA_FULL
                             # STRATA 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD
              ,DISCARD_SOURCE
              ,FED_OR_STATE)


strat_final <- rbind(strat_goldentilefish, strat_bluelinetilefish)


  strat_final %>% 
    filter(FED_OR_STATE=='FED') %>% 
#   dplyr::group_by(SPECIES, STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>% 
    #slice(1) %>% 
    ggplot()+
    geom_histogram(aes(x = CAMS_DISCARD_RATE, fill = DISCARD_SOURCE), bins=60, position = "dodge") +
    facet_wrap(~SPECIES)+
    theme_light()

strat_final %>% 
    filter(DISCARD_SOURCE != 'O' & FED_OR_STATE=='STATE') %>% 
    #dplyr::group_by(SPECIES, STRATA, DISCARD_SOURCE, FED_OR_STATE) %>% 
    #slice(1) %>% 
    ggplot()+
    geom_histogram(aes(x = CAMS_DISCARD_RATE, fill = DISCARD_SOURCE), bins=60, position = "dodge") +
    facet_wrap(~SPECIES)+
    theme_light()

db_final %>% 
    # filter(DISCARD_SOURCE != 'O')&
  filter(FED_OR_STATE=='FED') %>% 
    group_by(SPECIES, SPECIES_ITIS, DISCARD_SOURCE,SPECIES_STOCK) %>% 
    dplyr::summarise(DSUM = sum(CAMS_DISCARD, na.rm = T)) %>% 
    # slice(1) %>% 
    ggplot()+
    geom_bar(aes(x = SPECIES_STOCK, y = DSUM, fill = DISCARD_SOURCE), stat = 'identity', position = 'dodge')+
 facet_wrap(~SPECIES, scales = 'free')+
    theme_light()

db_final %>% 
    # filter(DISCARD_SOURCE != 'O')&
  filter(FED_OR_STATE=='STATE') %>% 
    group_by(SPECIES, SPECIES_ITIS, DISCARD_SOURCE,SPECIES_STOCK) %>% 
    dplyr::summarise(DSUM = sum(CAMS_DISCARD, na.rm = T)) %>% 
    # slice(1) %>% 
    ggplot()+
    geom_bar(aes(x = SPECIES_STOCK, y = DSUM, fill = DISCARD_SOURCE), stat = 'identity', position = 'dodge')+
 facet_wrap(~SPECIES, scales = 'free')+
    theme_light()

db_goldentilefish %>%   
  filter(DISCARD_SOURCE != 'O' & FED_OR_STATE=="FED") %>%
  group_by(STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                      # , N_OBS_TRIPS_F = mean(N_OBS_TRIPS_F)) %>% 
   ungroup() %>% dplyr::select(STRATA_FULL,
                                                                                                                                          DISCARD_SOURCE 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD) %>%
             # ,CV) 
             dplyr::arrange(desc(DISCARD)) %>%# top_n(10) %>%
        DT::datatable(caption = 'Discard rates for federal trips by strata for Goldentilefish')

db_goldentilefish %>%   
  filter(DISCARD_SOURCE != 'O' & FED_OR_STATE=="STATE") %>%
  group_by(STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                      # , N_OBS_TRIPS_F = mean(N_OBS_TRIPS_F)) %>% 
   ungroup() %>% dplyr::select(STRATA_FULL,
                                                                                                                                          DISCARD_SOURCE 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD) %>%
             # ,CV) 
             dplyr::arrange(desc(DISCARD)) %>%# top_n(10) %>%
        DT::datatable(caption = 'Discard rates for state trips by strata for goldentilefish')

db_bluelinetilefish %>%     
   filter(DISCARD_SOURCE != 'O' & FED_OR_STATE=="FED") %>%
  group_by(STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                      # , N_OBS_TRIPS_F = mean(N_OBS_TRIPS_F)) %>% 
   ungroup() %>% dplyr::select(STRATA_FULL,
                                                                                                                                          DISCARD_SOURCE 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD) %>%
             # ,CV) 
             dplyr::arrange(desc(DISCARD)) %>%# top_n(10) %>%
        DT::datatable(caption = 'Discard rates by Strata for Bluelinetilefish')

db_bluelinetilefish %>%     
   filter(DISCARD_SOURCE != 'O' & FED_OR_STATE=="STATE") %>%
  group_by(STRATA_FULL, DISCARD_SOURCE, FED_OR_STATE) %>%
  dplyr::summarise(SUBTRIP_KALL = sum(SUBTRIP_KALL)
                       , CAMS_DISCARD_RATE = round(mean(CAMS_DISCARD_RATE), digits=5)
                       , DISC_MORT_RATIO = mean(DISC_MORT_RATIO)
                      # , CV = mean(CV, na.rm=TRUE)
                       , DISCARD = round(sum(CAMS_DISCARD))) %>%
                      # , N_OBS_TRIPS_F = mean(N_OBS_TRIPS_F)) %>% 
   ungroup() %>% dplyr::select(STRATA_FULL,
                                                                                                                                          DISCARD_SOURCE 
             # ,N_OBS_TRIPS_F
              ,CAMS_DISCARD_RATE
              ,SUBTRIP_KALL
              ,DISC_MORT_RATIO
              ,DISCARD) %>%
             # ,CV) 
             dplyr::arrange(desc(DISCARD)) %>%# top_n(10) %>%
        DT::datatable(caption = 'Discard rates by Strata for Bluelinetilefish')

db_goldentilefish %>% filter(FED_OR_STATE == 'FED') %>% group_by(DISCARD_SOURCE, CAMS_GEAR_GROUP, SPECIES_STOCK) %>% dplyr::summarise(D = sum(CAMS_DISCARD), K = sum(SUBTRIP_KALL), Drate = mean(CAMS_DISCARD_RATE), Dmort = mean(DISC_MORT_RATIO)) %>% arrange(desc(D))%>% top_n(10) %>% dplyr::mutate(across(where(is.numeric), round, 4))%>% mutate(across(ends_with("D"), round, 0))%>%
        DT::datatable(caption = 'Discard totals for federal trips by stock, gear and discard source for Goldentilefish')

db_goldentilefish %>% filter(FED_OR_STATE == 'STATE') %>% group_by(DISCARD_SOURCE, CAMS_GEAR_GROUP, SPECIES_STOCK) %>% dplyr::summarise(D = sum(CAMS_DISCARD), K = sum(SUBTRIP_KALL), Drate = mean(CAMS_DISCARD_RATE), Dmort = mean(DISC_MORT_RATIO)) %>% arrange(desc(D))%>% top_n(10) %>% dplyr::mutate(across(where(is.numeric), round, 4))%>% mutate(across(ends_with("D"), round, 0))%>%
        DT::datatable(caption = 'Discard totals for state trips by stock, gear and discard source for Goldentilefish')

db_bluelinetilefish %>% filter(FED_OR_STATE == 'FED') %>% group_by(DISCARD_SOURCE, CAMS_GEAR_GROUP, SPECIES_STOCK) %>% dplyr::summarise(D = sum(CAMS_DISCARD), K = sum(SUBTRIP_KALL), Drate = mean(CAMS_DISCARD_RATE), Dmort = mean(DISC_MORT_RATIO)) %>% arrange(desc(D))%>% top_n(10) %>% dplyr::mutate(across(where(is.numeric), round, 4))%>% mutate(across(ends_with("D"), round, 0))%>%
        DT::datatable(caption = 'Discard totals for federal trips by stock, gear and discard source for Bluelinetilefish')

db_bluelinetilefish %>% filter(FED_OR_STATE == 'STATE') %>% group_by(DISCARD_SOURCE, CAMS_GEAR_GROUP, SPECIES_STOCK) %>% dplyr::summarise(D = sum(CAMS_DISCARD), K = sum(SUBTRIP_KALL), Drate = mean(CAMS_DISCARD_RATE), Dmort = mean(DISC_MORT_RATIO)) %>% arrange(desc(D))%>% top_n(10) %>% dplyr::mutate(across(where(is.numeric), round, 4 ))%>% mutate(across(ends_with("D"), round, 0))%>%
        DT::datatable(caption = 'Discard totals for state trips by stock, gear and discard source for Bluelinetilefish')