CohoFRAMModelValidation.Rmd

---
title: "Coho FRAM Model Validation and Mixed Stock Model Updating<br><br>"
author:
- Angelika Hagen-Breaux^[Washington Department of Fish and Wildlife, hagenafh@dfw.wa.gov]
- Carrie Cook-Tabor^[U.S. Fish and Wildlife Service, carrie_cook-tabor@fws.gov],
- Dan Auerbach^[Washington Department of Fish and Wildlife, Daniel.Auerbach@dfw.wa.gov]
date: "`r Sys.Date()`"
output:
  bookdown::html_document2:
    fig_caption: yes
    theme: cerulean
    toc: yes
    toc_depth: 3
    toc_float: yes
    code_folding: hide
  html_document:
    toc: yes
    toc_depth: '3'
    df_print: paged
  word_document: null
always_allow_html: yes
editor_options:
  chunk_output_type: console
  markdown:
    wrap: 72
subtitle: Southern Fund Project 2019-SP-3A<br><br><br>
---

<style type="text/css">
h1.title {
  text-align: center;
}
h3.subtitle {
  text-align: center;
}
h4.author { /* Header 4 - and the author and data headers use this too  */
  text-align: center;
}
h4.date { /* Header 4 - and the author and data headers use this too  */
  text-align: center;
}
</style>


```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE, results = FALSE, warning = FALSE, message = FALSE, fig.width = 9, fig.height = 9)

library("tidyverse")
library("patchwork")
library("gt")

# install.packages("devtools")
# devtools::install_github("BSchamberger/RDCOMClient")
# devtools::install_github("FRAMverse/framr")

library(framr)
library(rstatix)
library(ggpubr)
library(broom)
library(modelr)
library(reactable)

pal = c("#465177", "#E4C22B", "#965127")

table_dir<-"C:/Data/CoTC/SouthernFund/2019/2019SouthernFundProject/DataSF19/"

mdb<-"C:/data/CoTC/PostseasonRuns/PostDB/PSC_CoTC_PostSeason_CohoFRAMDB_thru2020_021622.mdb"  

d_fishmap<-read.csv(file.path(table_dir,"FishMapAug26_2022.csv"))
d_uniquefishery<-read.csv(file.path(table_dir,"MR_FishMapUnique.csv"))
d_uniquefishery2<-read.csv(file.path(table_dir,"MR_FishMapUnique2.csv"))
d_tSlu<-read.csv(file.path(table_dir,"ts_lu.csv"))
d_fishscalar<-framr::read_coho_fish_sclr(mdb,runs = 34:42) |> 
    select(RunYear, FisheryID, TimeStep, FisheryFlag, MarkReleaseRate)
d_exclude_fish<-read.csv(file.path(table_dir,"MR_ExcludeFish.csv"))
t_unclipped<-read.csv(file.path(table_dir,"MR_unclipped.csv"))
t_prmu<-read.csv(file.path(table_dir,"PR_MU.csv"))
d_stockmap<-read.csv(file.path(table_dir,"ER_Stock_CWT07112022.csv"))
d_ERtrollsum<-read.csv(file.path(table_dir,"ER_Ttrollsum.csv"))
d_top10<-read.csv(file.path(table_dir,"Fish43top10.csv"))

#needed to combine T/NT net pairs. Do not combine T/NT troll (different MSF regs)
Tpair<-c(80,82,87,96,101,109,111,119,121,123,130,132,137,139,141,143,145,153,155,157,159)

#load DA objects for ER_comp
list2env(readRDS("T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/ER_comp_fram_rmis_obj.rds"), .GlobalEnv)

```

```{r d_samp}
d_samp <- read_csv(file.path(table_dir, "MR Sampling Data.csv")) 

# get FishIDs for sport fisheries for later use
sampsportID <- unique(d_samp[d_samp$Gear=="Sport",3]) |>
  filter(!is.na(FisheryID))
sort(sampsportID$FisheryID)

#exclude some non-FRAM areas without FishID
d_samp <- d_samp|>
  filter(!is.na(FisheryID)) |> 
  select(RunYear = Year, FisheryID, TimeStep = TS, #Flag,
         MK = `Mrkd Landed`,
         UK = `UM Landed`,
         MR = `Mrkd Released`,
         UR = `UM Released`
  ) |> 
  mutate(
    FisheryID = if_else(FisheryID %in% Tpair, FisheryID+1, FisheryID) 
  ) |> #count(Flag)
  group_by(RunYear, FisheryID, TimeStep) |> #, Flag
  summarise(across(MK:UR, sum), .groups = "drop")  #count(Flag) 

# only calculate landed mark rate if sum of marked and unmarked >20
# calculate encounter mark rate if sum of all encounters >20
d_samp <- d_samp |>
  mutate(
    mr_kept = MK/(MK+UK),
    mr_enc = (MK+MR)/(MK+UK+MR+UR)
  ) |> 
  dplyr::rename_at(vars(MK:mr_enc), ~paste0(.,"_samp")) |> 
  filter(!is.na(mr_enc_samp))|>
  mutate(Source = "Sampling")

#use Encounters for MSF Sport and Landed for NS 
d_samp <- left_join(
  d_samp,
  d_fishscalar |> select(RunYear,FisheryID,TimeStep,FisheryFlag),
  by=c("RunYear","FisheryID","TimeStep"))

d_samp <- d_samp |>
  mutate(Type = ifelse(FisheryFlag>2,"Encounters","Landed"))|>
  mutate(
    TotalLanded = MK_samp+UK_samp,
    TotalEncounters = MK_samp+UK_samp+MR_samp+UR_samp,
    MrkdEncounters = MK_samp + MR_samp
  )|>
  mutate(
    Total_samp = ifelse(Type=="Encounters",TotalEncounters,TotalLanded),
    Marked_samp = ifelse(Type=="Encounters",MrkdEncounters,MK_samp),
    MR_sampl = ifelse(Type=="Encounters",mr_enc_samp,mr_kept_samp)
  )

d_samp <- d_samp[,c(1:3,16:18,10,11,12)] |>
  filter(!is.na(MR_sampl)) |>
  filter(!is.na(FisheryFlag)) |>
  filter(!FisheryFlag == 0) |>
  filter(!FisheryFlag == 28)

#eliminate MSF troll fisheries Type= Encounters FisheryID = 34,35,38,42. Sampling for these fisheries only records landed catch
troll <- c(34,35,38,42)
d_samp <- d_samp |>
  filter(!FisheryID %in% troll ) |>
  filter(!Total_samp<20)

```

```{r d_rmisraw}
# data pulled from RMIS CatchSample May 17, 2022 with only a species (Coho) and years (2010 to 2020) filter
d_rmisraw<-read.csv(file.path(table_dir,"MR_RMIS2010_20CohoMrkSample.csv"))
# fisheries that remained unassigned after the first round were added to 
# the fishmap file based on best fit, with some exceptions that remained unmapped
# such as 1M2, 1M3, 1M4, Oregon Net etc.

#Clean and summarize data
#maps CatchSample data to FRAM fisheries
d_rmisraw2<-d_rmisraw%>%
  add_column (FisheryID = 0)
# add blanks to end of recovery location code
d_rmisraw2$catch_location_code<-paste(d_rmisraw2$catch_location_code,"    ")

for (i in 1:nrow(d_fishmap)) {
  numchars<-d_fishmap[i,4]
  gear<-d_fishmap[i,3]
  FishID2<-d_fishmap[i,1]
  FishName2<-d_fishmap[i,2]
  locc<-substr(d_fishmap[i,5], 1, numchars)
  i=i+1
#find all the records in rmis_recrel with matching fishery and location codes
# update FishID and FishName fields to FishID and FishName
# assign fishery 50-59 to escapement  
  if (gear == 5){
    d_rmisraw2<-within(d_rmisraw2,FisheryID[floor(fishery/10) ==gear ]<-FishID2)
  } else {
     d_rmisraw2<-within(d_rmisraw2,FisheryID[floor(fishery/10) ==gear & substr(catch_location_code,start=1,stop=numchars) ==locc]<-FishID2)
 }  
}
#assign unassigned records with a freshwater sampling location to escapement
d_rmisraw2<- within(d_rmisraw2, FisheryID[FisheryID==0 & substr(catch_location_code,2,2) == "F"]<-300) 

#add Treaty ocean troll designation
NT_ocean_troll<-c(35,38,42)
d_rmisraw2$FisheryID = ifelse(d_rmisraw2$fishery ==15 & d_rmisraw2$FisheryID %in% NT_ocean_troll, d_rmisraw2$FisheryID + 1,d_rmisraw2$FisheryID)


#add Fishery Name
SoFtroll<-c(18,20,22)

d_rmisraw2<-left_join(d_rmisraw2,d_uniquefishery %>% select(FisheryID,FisheryName), by = c("FisheryID"))|>
  # eliminate Newport,Coos Bay, Tillamook troll, very low mark rate based on only one or two years of sampling with low sample rate
filter(!FisheryID %in% SoFtroll)

for (i in 1:nrow(d_rmisraw2)) {
  if (d_rmisraw2$FisheryID[i] %in% Tpair){
    d_rmisraw2$FisheryID[i]<- d_rmisraw2$FisheryID[i]+1
  } else {
    d_rmisraw2$FisheryID[i]<- d_rmisraw2$FisheryID[i]
 }  
}

d_junk<-d_rmisraw2|>
  group_by(FisheryID)|>
  summarize(n = n())
 

# eliminate records with no sampling to later compute a new sampling expansion
d_rmis_no_samp<-d_rmisraw2[,c(3,13,14,15,20,21)]|>
  filter(number_sampled>0,
         number_caught>0,
         sample_type!=5)|>
  select(-sample_type)|>
  dplyr::rename("RunYear"=catch_year)|>
  group_by(RunYear,FisheryID,FisheryName)|>
  summarise(across(c(number_caught,number_sampled),sum))


#for use in comparing FRAM and RMIS Catches see Data Preparation
rmis_landed<-d_rmisraw2[,c(3,14,15,20,21)]|>
  filter(!is.na(number_caught))|>
dplyr::rename("Total_rmis"=number_caught,
       "Sampled_rmis"=number_sampled,
       "RunYear"=catch_year)|>
  group_by(RunYear,FisheryID,FisheryName)|>
  summarise(across(c(Total_rmis,Sampled_rmis),sum),.groups="drop")


 # add time steps
d_rmis3<-left_join(d_rmisraw2,d_tSlu,by=c('period_type'='PeriodType','period'='Period'))|>
  #clean data
#Eliminate records where mark rate is na
filter(!is.na(mark_rate))|>
#eliminate records where number caught is na
filter(!is.na(number_caught))|>
#eliminate records where number caught is zero
filter(!number_caught==0) |>
#eliminate records where number sampled is na
filter(!is.na(number_sampled))|>
#eliminate records where number sampled is zero
filter(!number_sampled==0)|>
#calculate the number marked
mutate(number_marked = number_caught * mark_rate) |> 
#eliminate records where adclip_selective_fishery = M (mixed) 
filter(!adclip_selective_fishery=="M")|>
#eliminate Skagit FW sport. Only one year of sampling in the Cascade which appears to have a different mark rate (higher) than the Skagit River fishery modeled in FRAM. Eliminate Quinault net. Only one record that has a mark sample. This record has a mark rate of zero. This result is highly improbable.
filter(!FisheryID %in% c(108,63))


#sum number marked by fishery, year, and time step
d_rmis<-d_rmis3 |>
  group_by(catch_year, FisheryID, TimeStep) |> 
summarise(across(c(number_caught,number_marked,number_sampled), sum), .groups = "drop")|>
  mutate(MR_rmis = number_marked/number_caught)|>
#eliminate records where number sampled < 20
filter(!number_sampled<20)
d_rmis<-d_rmis[,c(1:5,7)]|>
  mutate(Source = "RMIS")|>
  mutate(Type = "Landed")
colnames(d_rmis)[1]<- "RunYear"
colnames(d_rmis)[4]<- "Total_rmis"
colnames(d_rmis)[5]<- "Marked_rmis"
colnames(d_rmis)[2]<-"FisheryID"

```

```{r d_fram}
d_fram <- left_join(
  framr::read_coho_mort(mdb, runs = 34:44) |> 
    select(RunYear, StockID, FisheryID, TimeStep, LandedCatch, MSFLandedCatch, MSFNonRetention) |> 
    mutate(RunYear = as.numeric(RunYear))
  ,
  d_fishscalar
  ,
  #d_fram_mort, d_fram_fs,
  by = c("RunYear","FisheryID","TimeStep")
) |> #filter(FisheryID %in% 96:97) |> count(FisheryID, FisheryName)
  mutate(
    MarkStatus = if_else(StockID %% 2 == 0, "M", "UM"),
    FisheryID = ifelse(FisheryID %in% Tpair, FisheryID+1, FisheryID),
    Kept = LandedCatch + MSFLandedCatch,
    Enc = LandedCatch + MSFLandedCatch + if_else(MarkReleaseRate>0,MSFNonRetention/MarkReleaseRate,0)
  ) 
#re-flag 1s to 2s.In post-season runs fisheries should all be flagged as 2. This will avaoid potential duplicates for T/NT fishing pairs
  d_fram[which(d_fram$FisheryFlag == 1),]$FisheryFlag=2
 
d_fram<-d_fram|>
  group_by(RunYear, FisheryID,MarkStatus, TimeStep,FisheryFlag)|>
  summarize(across(Kept:Enc,sum),.groups = "drop")|>
  pivot_wider(names_from = MarkStatus, values_from = c(Kept,Enc))|>
  mutate(
    mr_kept = Kept_M/(Kept_M+Kept_UM),
    mr_enc = Enc_M/(Enc_M+Enc_UM)
  ) |>
  dplyr::rename_at(vars(Kept_M:mr_enc), ~paste0(.,"_fram"))  
  
# for landed catch comparison with RMIS
fram_landed<-d_fram[,c(1:6)]|>
  mutate(Total_fram =Kept_M_fram + Kept_UM_fram )|>
  group_by(RunYear,FisheryID)|>
summarise(across(c(Kept_M_fram,Total_fram),sum),.groups="drop")
     
d_fram<-d_fram|>     
  filter(!is.na(mr_enc_fram))|>
  filter(!FisheryFlag>8)|> #exclude mixed MSF/NS fisheries
  mutate(Source = "FRAM")

d_fram<-d_fram[,c(1:4,7,8,10,11)]|>
  mutate(Type = "Encounters",
         Total = Enc_M_fram + Enc_UM_fram
         )
  colnames(d_fram)[5]<-"Marked"
  colnames(d_fram)[7]<-"MR"
d_fram<-d_fram[,c(1:4,10,5,7:9)]

```

```{r joined_MR_objs}
 #nothing in sampling dataset that isn't in fram dataset
setdiff(unique(d_samp$FisheryID), unique(d_fram$FisheryID))
#but lots of fram fisheries for which we don't have sampling obs...
sort(setdiff(unique(d_fram$FisheryID), unique(d_samp$FisheryID)))

#join to sampling if MSF, join to RMIS if non-selective

# d_samp<-d_samp|>
#   group_by(RunYear,FisheryID,TimeStep,Source,Type)|>
#   summarise(across(c(Total_samp,Marked_samp,MR_sampl),sum),.groups="drop")

d_fram_joined1 <- d_fram|>
  filter(FisheryFlag ==8)|>
  inner_join(
    d_samp[,-c(8)],
    by=c("RunYear","FisheryID","TimeStep")
)

colnames(d_fram_joined1)[10]<-"Total"
colnames(d_fram_joined1)[11]<-"Marked"
colnames(d_fram_joined1)[12]<-"MR"


d_fram_joined2<-d_fram|>
  filter(!FisheryFlag==8)|>
  inner_join(
    d_rmis,
    by=c("RunYear","FisheryID","TimeStep")
)
colnames(d_fram_joined2)[10]<-"Total"
colnames(d_fram_joined2)[11]<-"Marked"
colnames(d_fram_joined2)[12]<-"MR"

mr_fram_samp_rmis <-rbind(d_fram_joined1,d_fram_joined2)

colnames(mr_fram_samp_rmis)[10]<-"Tot_samp"
colnames(mr_fram_samp_rmis)[11]<-"Marked_samp"
colnames(mr_fram_samp_rmis)[12]<-"MR_samp"

mr_ts_fram_samp_rmis <-left_join(mr_fram_samp_rmis,d_uniquefishery %>% select(FisheryID,FisheryName), by = c("FisheryID"))


#sum over time steps before computing annual rates. It's okay to ignore Flags and Type (based on data QAQC). Conduct analyses on annual mark rates, because many FRAM fisheries are not reported in the time step they actually occur; i.e., PS fw fisheries are modeled in time 5 even if most catch occurs in time 4.

d_fram_yr<-d_fram|>
  group_by(RunYear,FisheryID,FisheryFlag,Source,Type)|>
  summarise(across(c(Total,Marked,MR),sum),.groups="drop")|>
  mutate(MR=Marked/Total)|>
  left_join(d_uniquefishery|>
              select(FisheryID,FisheryName), by = c("FisheryID"))

d_samp_yr<-d_samp|>
  filter(Type != "Landed")|>
  group_by(RunYear,FisheryID,Source,Type)|>
  summarise(across(c(Total_samp,Marked_samp,MR_sampl),sum),.groups="drop")|>
  mutate(MR_sampl=Marked_samp/Total_samp)

#include RMIS time steps for non-selective fisheries that are modeled in a different time step in FRAM

d_rmis_yr<-d_rmis|>
  left_join(d_fram|>
              select(RunYear,FisheryID,TimeStep,FisheryFlag), by = c("RunYear","FisheryID","TimeStep"))|>
rowid_to_column("uid")

d_rmis_yr<-d_rmis_yr|>
  nest_by(FisheryID, RunYear) |>
  #same as: group_by(FisheryID, RunYear) |> nest() |> rowwise() |>
  filter(
    nrow(data) >= 2 &
      nrow(data |> filter(is.na(FisheryFlag) & MR_rmis < 0.95) > 0) &
        nrow(data |> filter(FisheryFlag == 2) > 0)
    ) |>
  unnest(data) |> ungroup() |>
  split(f = ~FisheryID + RunYear, drop = T) |>
  map_df(
    function(d){
      ts_na = d |> filter(is.na(FisheryFlag) & MR_rmis < 0.95) |> pluck("TimeStep")
      
      if(length(ts_na) > 1) {
        map_df(
          ts_na,
          function(ts){
            ff_pre = d |> filter(TimeStep == ts - 1) |> pluck("FisheryFlag", .default = NA)
            ff_pst = d |> filter(TimeStep == ts + 1) |> pluck("FisheryFlag", .default = NA)
            if(!all(is.na(c(ff_pre, ff_pst)))){
              if(ff_pre == 2 | ff_pst == 2){
                d |> filter(TimeStep == ts) |> select(uid) |> mutate(flag = "3")
              }
            } 
          })
        
      } else {
        ff_pre = d |> filter(TimeStep == ts_na - 1) |> pluck("FisheryFlag", .default = NA)
        ff_pst = d |> filter(TimeStep == ts_na + 1) |> pluck("FisheryFlag", .default = NA)
        if(!all(is.na(c(ff_pre, ff_pst)))){
          if(ff_pre == 2 | ff_pst == 2){
            d |> filter(TimeStep == ts_na) |> select(uid) |> mutate(flag = "3")
          }
        }
      }
    }
  ) %>%
  left_join(d_rmis_yr, y = ., by = "uid")|>
  filter(FisheryFlag==2 | flag==3)|>
  #eliminate non-selective fisheries with MR >0.98% These are likely MSFs placed in the wrong time step due to RMIS management week structure
filter(MR_rmis<0.98)

# for (i in 1:nrow(d_rmis_yr)) {
#   Tstep<-d_rmis_yr$TimeStep[i]
#   runyr<- d_rmis_yr$RunYear[i] #d_rmis_yr[i,"RunYear"]
#   fid<-d_rmis_yr$FisheryID[i]
#   fflag<-d_rmis_yr$FisheryFlag[i]
#   mrkrate<-d_rmis_yr$MR_rmis[i]
#   if (is.na(fflag) & mrkrate<0.95) {
#     other_ts <- d_rmis_yr[ 
#       d_rmis_yr$RunYear == runyr &
#       d_rmis_yr$FisheryID == fid & 
#       d_rmis_yr$FisheryFlag == 2 &
#       (d_rmis_yr$TimeStep == Tstep + 1 | d_rmis_yr$TimeStep == Tstep - 1 ) ,]
#     if(nrow(other_ts) > 0) {
#       d_rmis_yr$FisheryFlag[i]<-3
#     }
#   }
# }



d_rmis_yr<-d_rmis_yr[,c(2:9)]|>
   group_by(RunYear,FisheryID,Source,Type)|>
  summarise(across(c(Total_rmis,Marked_rmis,MR_rmis),sum),.groups="drop")|>
  mutate(MR_rmis=Marked_rmis/Total_rmis)

mr_yr_fram_join_samp<-d_fram_yr|>
  filter(FisheryFlag ==8)|>
  inner_join(
    d_samp_yr,
    by=c("RunYear","FisheryID")
)
mr_yr_fram_join_rmis<-d_fram_yr|>
filter(FisheryFlag !=8)|>
  inner_join(
    d_rmis_yr,
    by=c("RunYear","FisheryID")
  )|>
  dplyr::rename(Total_samp = Total_rmis)|>
  dplyr::rename(Marked_samp = Marked_rmis)|>
  dplyr::rename(MR_sampl=MR_rmis)

mr_yr_fram_samp_rmis<-rbind(mr_yr_fram_join_samp,mr_yr_fram_join_rmis)|>
  mutate(RelPctDiff = (MR_sampl-MR)/MR_sampl,
         AbsPctDiff = (MR_sampl-MR))
  

# compare landed RMIS to landed sampling for non-selective fisheries
       
# mr_rmis_samp<-left_join(mr_rmis_samp,d_uniquefishery %>% select(FRAM_Fish_ID,FRAM_Fish_Name), by = c("FisheryID"="FRAM_Fish_ID"))
mr_rmis_samp<-d_samp|>
  #filter(Type=="Landed")|>
  inner_join(
  d_rmis,
  by=c("RunYear","FisheryID","TimeStep")
)|>
  left_join(
    d_uniquefishery,
    by=c("FisheryID")
  )
# #add fishery flag
# mr_rmis_samp<-left_join(mr_rmis_samp, d_fram%>%                          select(RunYear,FisheryID,TimeStep,FisheryFlag),
#                         by=c("RunYear","FisheryID","TimeStep"))
#eliminate MSF
mr_rmis_samp<- mr_rmis_samp|>
  filter(FisheryFlag<=2)


mr_fram_samp_long <- left_join(mr_fram_samp_rmis,d_uniquefishery %>% select(FisheryID,FisheryName), by = c("FisheryID"))|>
  select(RunYear, FisheryID, FisheryName, TimeStep, starts_with("MR")) |> 
  pivot_longer(names_to = "type", values_to = "val", cols = starts_with("MR"))
colnames(mr_fram_samp_long)[3]<-"FisheryName"
mr_fram_samp_long$type[mr_fram_samp_long$type!="MR_samp"]<-"MR_fram"

mr_rmis_samp_long <- mr_rmis_samp |> 
  select(RunYear, FisheryID, FisheryName, TimeStep, MR_sampl ,MR_rmis) |>
  pivot_longer(names_to = "type", values_to = "val", cols =  c("MR_sampl","MR_rmis"))
colnames(mr_rmis_samp_long)[3]<-"FisheryName"

```

\newpage

# Overview

## Acknowledgements

Financial support for this project was provided by the Southern Endowment Fund, with support from the Washington Department of Fish and Wildlife, U.S. Fish and Wildlife Service, and members of the Coho Technical Committee. <br>
The authors also wish to thank Jeremiah Shrovnal, Derek Dapp, Marlene Bellman, and Jon Carey for their help and feedback.

\newpage

## Glossary

<a id='M'></a>**Adipose Mark**<br>
A missing adipose fin on salmon as a result of clipping by a hatchery program. Commonly, a visual identifier of a hatchery fish that can be retained in a [mark selective fishery](#MSF). The use of "marked" and "clipped" is used interchangeably within the document.

<a id='BP'></a>**Base Period (BP)** <br>
A range of years from which CWT data are used to estimate exploitation rates and other parameters through a process of cohort reconstruction. Resulting base period reference parameters are used to populate the FRAM model to predict annual stock/fishery specific impacts.

<a id='BkFRAM'></a>**Backwards FRAM (BkFRAM)** <br>
A [FRAM](#FRAM) utility that reconstructs starting Cohorts by adding fisheries mortalities and escapements for the purpose of producing post-season model runs (Hagen-Breaux 2018). FRAM iteratively adjusts stock recruit scalars (a surrogate for [starting cohorts](#SC)) and runs FRAM forward until the resulting escapements match the target escapements.

<a id='CRC'></a>**Catch Record Cards (CRC)** <br>
Catch Record Cards are used to estimate the recreational catch of salmon, steelhead, sturgeon, halibut and Puget Sound Dungeness crab. The CRC system houses the recreational catch estimates.

<a id='CWT'></a>**Coded-Wire Tag (CWT)** <br>
Coded micro-wire implanted in juvenile salmon prior to release. When recovered, the binary or numeric code on the tag identifies the tag group, and thus provides information such as location and timing of release, special hatchery treatments, etc.

<a id='CoTC'></a>**Coho Technical Committee (CoTC)** <br>
Technical committee formed by the [Pacific Salmon Commission](#PSC) providing the panels with scientific data to make decisions affecting Southern Coho salmon.

<a id='ER'></a>**Exploitation Rate (ER)**<br>
The mortality (landed or total) of a stock or stock aggregate in a fishery(s) divided by the stock's abundance (fishery mortality + escapement).

<a id='ET'></a>**Extreme Terminal**<br>
Marine area directly adjacent to salmon native streams with a very high proportion of catch made up of the "local stock" due to river mouth proximity, i.e., Elliott Bay is the extreme terminal area for Green River Chinook and Coho.

<a id='FRAM'></a>**Fishery Regulation Assessment Model (FRAM)** <br>
[Fishery simulation model](https://framverse.github.io/fram_doc/index.html) used by the Pacific Salmon Commission (Coho only), Pacific Fishery Management Council, Washington Department of Fish and Wildlife, and Puget Sound Treaty Tribes to evaluate the effects of fisheries on Chinook and Coho salmon and to evaluate management objectives.

<a id='fulljoin'></a>**Full Join** <br>
A data join returning all the rows from the joined tables, whether they are matched or not.

<a id='MR'></a>**Mark Rate** <br>
The mark rate of a fishery is computed as [adipose marked](#M) landed catch or encounters divided by total landed catch or total encounters in the fishery.

<a id='MSF'></a>**Mark-Selective Fishery (MSF)** <br>
A fishery in which only marked (adipose fin clipped) fish can be legally retained.

<a id='MSM'></a>**Mixed Stock Model (MSM)** <br>
Main computer program used to build the current Coho [base period](#BP). Model inputs included [coded-wire tags](#CWT) and estimated catches and escapements of Coho coast-wide from catch years 1986-1992. Model output includes base period [exploitation rates](#ER) utilized in contemporary [FRAM](#FRAM) runs.

<a id='NOF'></a>**North of Falcon (NOF)**<br>
Co-manager, annual, pre-season Puget Sound salmon fishery planning process; referencing the region North of Cape Falcon, Oregon, which marks the southern border of active management for Washington salmon stocks.
 
<a id='PSC'></a>**Pacific Salmon Commission [(PSC)](https://www.psc.org/)** <br>
A commission composed of Canadian and United States representatives tasked with developing and upholding the [Pacific Salmon Treaty](#PST).

<a id='PFMC'></a>**Pacific Fishery Management Council [(PFMC)](https://www.pcouncil.org/)** <br>
The Pacific Fishery Management Council manages fisheries for approximately 119 species of salmon, groundfish, coastal pelagic species (sardines, anchovies, and mackerel), and highly migratory species (tunas, sharks, and swordfish) on the West Coast of the United States.

<a id='PST'></a>**Pacific Salmon Treaty [(PST)](https://www.psc.org/about-us/history-purpose/pacific-salmon-treaty/)** <br>
Treaty between the Governments of Canada and the United States of America concerning Pacific Salmon.

<a id='PTerm'></a>**Preterminal Fisheries**<br>
Fisheries on salmon in marine preterminal areas. Preterminal areas are distant to natal streams and consist of predominately mixed stocks. 

<a id='RMIS'></a>**Regional Mark Information System (RMIS)** <br>
[RMIS](https://www.rmpc.org/) is an online data repository, storing Pacific Coast-wide anadromous salmon [adipose marking](#M) information, [coded-wire-tag](#CWT) releases and recoveries, recovery locations, and sampling information.

<a id='SOF'></a> **South of Falcon (SOF)**<br>
Referencing the region South of Cape Falcon, Oregon, which marks the southern border of active management for Washington State salmon stocks.

<a id='SC'></a> **Starting Cohort**<br>
Initial adult Coho run sizes prior to subtracting natural         and fishing mortality at the beginning of the run year.

<a id='Term'></a>**Terminal Fisheries**<br>
Fisheries on adult salmon in marine terminal areas. Terminal areas are close to the natal stream and consist of predominately local stocks. 

<a id='TS'></a>**Time Steps**<br>
Coho [FRAM](#FRAM) has 5 time steps.<br>
- Time 1: January – June<br>
- Time 2: July<br>
- Time 3: August<br>
- Time 4: September<br>
- Time 5: October-December

<a id='VR'></a>**Validation Runs** <br>
Post-season FRAM model runs populated with observed mortalities and escapements, reflecting the "best information". [Exploitation rates](#ER) resulting from these runs are compared to pre-season exploitation rates to assess whether management objectives were achieved. 


---

```{js, echo=FALSE}
// Place footnotes after glossary section
$(document).ready(function() {
  $('.footnotes ol').appendTo('#glossary');
  $('.footnotes').remove();
});
```
\newpage

## List of Figures
[Figure 4.1. Exploitation Rates of Key Wild Stocks Resulting from Varying the Abundance of Test Stocks](#ERWild)<br>

[Figure 4.2. Exploitation Rate Differences Resulting from Varying the Abundances of Test Stocks from High to Low](#ERWildDiff)<br>

[Figure 5.1. Differences between Post-season Coho FRAM Modeled Mark Rates and Sampling Observations, 2010-18](#MRDiff)<br>

[Figure 5.2. 2010-18 Average Fishery mark Rates from FRAM and Sampling Observations](#AvgMR)<br>

[Figure 5.3. Wilcoxon p-values by Fishery](#Wilcox1)<br>

[Figure 5.4. 2010-18 Fishery Mark Rates from FRAM and Sampling by Year](#MRbyYr)<br>

[Figure 6.1. Recovery Counts by Time Period and Type](#RecCount)<br>

[Figure 6.2. Annual CWT Recoveries by Release Location in Fisheries and Escapement](#AnnualRecCount)<br>

[Figure 6.3.1. Skagit Pseudo ERs and Proportional Fishery Rankings](#Skagit1)<br>

[Figure 6.3.2. Tulalip Pseudo ERs and Proportional Fishery Rankings](#Tulalip1)<br>

[Figure 6.3.3. Quilcene Pseudo ERs and Proportional Fishery Rankings](#Quilcene1)<br>

[Figure 6.3.4. George Adams Pseudo ERs and Proportional Fishery Rankings](#GeorgeAdams1)<br>

[Figure 6.3.5. Nisqually Pseudo ERs and Proportional Fishery Rankings](#Nisqually1)<br>

[Figure 6.3.6. Puyallup Pseudo ERs and Proportional Fishery Rankings](#Puyallup1)<br>

[Figure 6.3.7. Green Pseudo ERs and Proportional Fishery Rankings](#Green1)<br>

[Figure 6.3.8. Quillayute Pseudo ERs and Proportional Fishery Rankings](#Quillayute1)<br>

[Figure 6.3.9. Queets Pseudo ERs and Proportional Fishery Rankings](#Queets1)<br>

[Figure 6.3.10. Quinault Pseudo ERs and Proportional Fishery Rankings](#Quinault1)<br>

[Figure 6.3.11. Chehalis Pseudo ERs and Proportional Fishery Rankings](#Chehalis1)<br>

[Figure 6.3.12. Willapa Pseudo ERs and Proportional Fishery Rankings](#Willapa1)<br>

[Figure 6.3.13. Columbia Earlies Pseudo ERs and Proportional Fishery Rankings](#ColumbiaEarlies1)<br>

[Figure 6.3.14. Columbia Lates Pseudo ERs and Proportional Fishery Rankings](#ColumbiaLates1)<br>

[Figure 6.3.15. Strait of Georgia Vancouver Island Pseudo ERs and Proportional Fishery Rankings](#SGVI1)<br>

[Figure 6.3.16. Fraser Lower Pseudo ERs and Proportional Fishery Rankings](#FraserLower1)<br>

[Figure 6.3.17. Fraser Upper Pseudo ERs and Proportional Fishery Rankings](#FraserUpper1)<br>


\newpage

## List of Tables
[Table 3.1. Unique Release Codes Mapped](#Unique)<br>

[Table 3.2. Columbia River Stock Assignments and Numbers Released](#ColRStk)<br>

[Table 3.3. Stock CWT 1998](#StkCWT)<br>

[Table 3.4. Fisheries Excluded from CWT Analysis and Main Reasons for Exclusion](#ExcludeFish)<br>

[Table 4.1. List of Response Stock Aggregates](#ResponseStk)<br>

[Table 4.2. List of Test Stock Aggregates](#TestStk)<br>

[Table 4.3. 2010-2018 Starting Cohorts of Coho Test Stocks](#StartingCohorts)<br>

[Table 4.4. List of Model Runs, Fisheries, and Time Steps with Catch Reductions to Avoid Negative Escapements](#NegEsc)<br>

[Table 4.5. Exploitation Rates of Key Wild Response Stocks Resulting from Varying Abundances of Test Stocks](#ERKey)<br>

[Table 4.6. Average Exploitation Rates Differences of Key Wild Stocks Resulting from Varying the Abundances of Test Stocks ](#ERDiff)<br>

[Table 4.7. Top 5 Fisheries Responsible for Total Mortality Differences Caused by Varying the Abundances of Test Stocks](#Top5)<br>

[Table 5.1. Sources for Mark Rate Comparisons by Fishery Type](#MRSource)<br>

[Table 5.2. Number of Fisheries Evaluated with Wilcoxon and Number of Fisheries with Significant Mark Rate Differences in Each Area](#Wilcox2)<br>

[Table 5.3. 'Good', 'Moderate', 'Poor' Agreement between FRAM and Sampling Mark Rates](#GMP)<br>

[Table 5.4. Mean Mark Rate Differences by Fishery Area](#MeanMRDiff)<br>

[Table 5.5. Comparison of FRAM and RMIS Catches by Mark Status](#MarkStatus)<br>

[Table 5.6. Comparison of FRAM and RMIS Marked Catches in the Area 4/4B Treaty Troll Fishery](#MarkedCatch)<br>

[Table 6.1. Number of CWT Recoveries for Assessed Stocks, 2010-1019](#FocalCWT)<br>


\newpage

## Executive Summary

Evaluations of Coho salmon (*Oncorhynchus kisutch*) exploitation rates ([ERs](#ER)) in fisheries along the Canadian and U.S. West Coast rely on the Fishery Regulation Assessment Model [(FRAM)](#FRAM). Coded-wire-tag [(CWT)](#CWT) recoveries from 1986-1992, also referred to as the [base-period](#BP), were used to develop the stock-fishery-time period specific exploitation rates that act as core model parameters driving FRAM calculations of the impacts associated with proposed and observed fisheries and returns.

The main goal of this project was to assess whether [FRAM](#FRAM) outputs based on decades old parameterization reflect current patterns of stock distribution within mixed stock fisheries. This objective was addressed by (1) comparison of post-season FRAM ["pseudo-ER"](#Pseudo) estimates against rates derived directly from contemporary [CWT](#CWT) recoveries, and (2) through comparison of fishery mark-rates produced by [FRAM](#FRAM) against reported values in those same fisheries. This project also aimed to use sensitivity analyses to study how changes in the abundance of [FRAM](#FRAM) stock aggregates affect the resulting modeled [exploitation rates](#ER) of key co-migrating stocks. Depending on these assessments, the final project goal was to investigate potential pathways to updating the existing base-period. 

**[Exploitation Rate](#ER) Analysis**<br>

[CWT](#CWT) data queried from the Regional Mark Information System [(RMIS)](#RMIS) were used for comparison of [exploitation rates](#ER) against post-season [FRAM](#FRAM). Error checking these data and determining their mapping to FRAM stocks and fisheries formed major elements of the project. In the course of mapping RMIS recoveries to FRAM fisheries we discovered that the **Southeast Alaska sport fishery is entirely missing from FRAM**. The updated look-up tables that reference [RMIS](#RMIS) recovery locations to FRAM fisheries and [RMIS](#RMIS) release information to FRAM stocks can form a valuable starting point to any future work seeking to build on the results presented here.

Initial exploration of these data narrowed fisheries to mainly those in mixed stock areas with adequate sampling and to selection of several focal stocks with sufficient [CWT](#CWT) representation in these fisheries. Most freshwater and many [extreme terminal](#ET) fisheries were excluded due to sampling limitations and the presumption that catches consisted of predominantly local stocks. [Pseudo-exploitation rates](#Pseudo) were then calculated by dividing [FRAM](#FRAM) landed catches and [RMIS](#RMIS) [CWT](#CWT) expanded recoveries by the sum of these fisheries plus escapement. In addition, the rank abundances of fishery-specific mortality proportions were calculated to provide a complementary depiction of among-fishery patterns unaffected by issues related to monitoring and sampling escapement.

**It was initially hypothesized that values of these measures from Coho [FRAM](#FRAM) [validation runs](#VR) for the fishing years 2010-2019 would correspond poorly with those calculated directly from [RMIS](#RMIS) recoveries during the same time span. However, for the selected stocks and fisheries, [FRAM](#FRAM) output exhibited better than expected agreement in terms of the relative ranking of impacts among pre-terminal fisheries and generally modest differences in the calculated [pseduo-ERs](#Pseudo).** Most larger differences between the two datasets were plausibly explained by sampling limitations in fisheries, sampling limitations in escapement, or both. Potential errors in the post-season [FRAM](#FRAM) datasets of annual catch and escapement also likely affected individual stocks, fisheries, and years, indicating the need to **review and update these critical inputs** to Coho FRAM "validation" runs as high priority next step for model improvements.

**[Mark Rate](#MR) Analysis** <br>

External mark status information (i.e., the presence or absence of an adipose fin) is collected by professional fishery samplers and reported to [RMIS](#RMIS). For the [mark rate](#MR) analysis [adipose marked](#M) statistics from [RMIS](#RMIS)' Catch/Sample database were summarized by [FRAM](#FRAM) fishery and time period and compared to the corresponding FRAM output. [Mark rates](#MR) were calculated for each fishery as the number of marked Coho encountered and/or landed divided by all Coho encountered and/or landed. In [FRAM](#FRAM), the fishery-wide [mark rate](#MR) results from aggregation across the stocks that are modeled as caught or encountered, and divergence of this composite rate from observations may therefore indicate errors in the modeled relative abundances and/or stock composition of the fishery.

This project also compared [RMIS](#RMIS) [mark rates](#MR) to mark rates directly reported by the sampling programs. As [RMIS](#RMIS) mark rates are derived from these data, we found, perhaps unsurprisingly, that [RMIS](#RMIS) and sampling mark rates closely match. Upon examination, the small differences in mark rates were caused by slight differences in [time step](#TS) definitions.

**Agreement between [FRAM](#FRAM) [mark rates](#MR) and sampling mark rates was worse than expected in many fisheries.** In particular, fishery mark rates in [FRAM](#FRAM) were consistently higher than sampled rates in several Washington coastal and ocean fisheries, and these differences warrant [further investigation](#casestudy). We were unable to conclusively determine whether [mark rate](#MR) differences were caused by inaccurate stock distributions, incorrect abundance inputs, or other data errors. [FRAM](#FRAM) and [RMIS](#RMIS) data have various shared and distinct sources of inaccuracy. In some cases, preliminary findings suggest that missing wild escapements may be responsible for some higher than observed FRAM mark rates (see [Case Study](#casestudy)).

**Sensitivity Analysis**<br>

Stock aggregate abundances were artificially adjusted from the starting values in several [validation runs](#VR) selected to capture qualitatively distinct "good" and "bad" return years. One-at-a-time increases and decreases were evaluated for the effects on modeled [exploitation rates](#ER) of other model "response" stocks, providing an indication of the relative importance of among-stock interactions in the model.

**With the exception of the Lower Fraser aggregate, the Canadian Coho stocks included in the [Pacific Salmon Treaty](#PST) have a relatively minor effect on the [FRAM](#FRAM)-modeled [exploitation rates](#ER) of U.S. Puget Sound and Coastal stocks. However, several Puget Sound stocks can significantly affect both Canadian and other Washington stocks.** In particular, changes to the Nooksack/Samish, Hood Canal, and South Puget Sound stock aggregates were influential for other modeled stock units. While improved forecast accuracy may benefit all model stocks, investments to increase the performance of preseason projections could yield the broadest returns if they were targeted to these aggregates (Nooksack/Samish, Hood Canal, South Puget Sound, and Lower Fraser). Similarly, dedicating additional resources to refined escapement estimation for these model units may tend to have greater influence on the overall accuracy of model results.

**[Updating the Mixed Stock Model](#MSM)**<br>

The [Mixed Stock Model (MSM)](#MSM) computer program was developed to estimate stock composition in [pre-terminal](#Pterm) mixed stock fisheries using [coded-wire-tag](#CWT) recovery data. Output from the [MSM](#MSM) was combined with [terminal](#Term) fishery data before producing cohort reconstructions and calculating the [base period](#BP) [exploitation rates](#ER) that form the foundation of Coho [FRAM](#FRAM).

This project investigated steps associated with combining contemporary [CWT](#CWT) recoveries with older [base period](#BP) CWTs to update base period [exploitation rates](#ER). **We identified several potential challenges and considerations, including a lack of CWT sampling in some contemporary fisheries as well as changes to marking and tagging practices that ultimately prevented a model update. We believe that updating the MSM would require an intense, concerted effort by regional experts to address the many data shortcomings.**

To calculate [base period](#BP) [exploitation rates](#ER) for contemporary tags, modelers must develop estimates of cohort size (i.e., the ER denominator). Because cohort size combines escapement and mortality, it is not possible to accurately estimate cohort size while excluding poorly sampled fisheries. Incorrect cohort sizes influence base period exploitation rate calculations for all fishery-[time steps](#TS) in the model, and therefore, even fisheries that are sampled well will require correct cohort sizes to estimate base exploitation rates.

However, to facilitate future [base period](#BP) development and updates, the following products were developed:

 - A description of challenges and potential solutions related to base period changes
 - Development of [Mixed Stock Model](#MSM) documentation
 - Updated instructions to map [RMIS](#RMIS) recoveries to [FRAM](#FRAM) stocks and fisheries

\newpage

# Introduction

The Coho [FRAM](#FRAM) is a forward projecting computer application that uses estimates of stock abundances and fishery mortalities scaled to a [base period](#BP) average.  This model is the primary tool used to evaluate performance of fisheries regimes by the Pacific Fishery Management Council [(PFMC)](#PFMC) and the parties to the Pacific Salmon Treaty [(PST)](#PST). Improvement, validation, and updating of this model are high priority activities needed to ensure accurate estimates of fishery impacts for stocks ranging from Alaska to California. [FRAM](#FRAM) is based upon historical (catch years 1986 through 1992) coded-wire-tag [(CWT)](#CWT) data to assess stock specific [exploitation rates](#ER). Members of the Coho Technical Committee [(CoTC)](#CoTC) and Coho Workgroup have grown increasingly concerned about the ability of this 30-year-old [base period](#BP) to accurately reflect current stock distribution patterns. To re-establish confidence in the model and assure that changes in stock distribution and fishery patterns are accurately captured, contemporary assessments are needed. An updated assessment of stock migratory patterns is particularly relevant given concerns about changing environmental conditions.

Unfortunately, [CWT](#CWT) recoveries have been in steep decline for the past two decades ([see recovery patterns since 1986](#RP)), with too few recoveries to complete full-scale coast-wide cohort reconstructions. Updating stock distribution patterns, which are modeled through the use of static [base period](#BP) [exploitation rates](#ER), is extremely data intense and time consuming, with the current [base period](#BP) requiring over 10 years of development time. Rather than attempting to create a brand new [base period](#BP), a project that previously failed due to its large scope and sparsity of data, this project partitioned the tasks into four smaller, achievable sub-components.

1. A sensitivity analysis evaluating which stocks have the greatest influence on the [exploitation rates](#ER) of key stocks of concern.
2. A comparison of fishery [adipose mark](#M) rates from sampling to [FRAM](#FRAM) estimated [mark rates](#MR).
3. A comparison of contemporary [exploitation rates](#ER) to [base period](#BP) derived exploitation rates
4. Investigation of potential pathways to updating the existing base-period.<br>

**The R computer code for this project as well as supporting data files are located on [GitHub](https://github.com/PSC-CoTC/coho_CWT_analysis). Many data tables can also be downloaded directly from the tables in the [Appendix](#A).**
```{r}
# focal_cols_recs <- c("recovery_id", "tag_code", "run_year", "recovery_location_code", "recovery_location_name", "fishery", "gear", "estimated_number")
# 
# rmis_recs<-purrr::map_dfr(
#   list.files(dir_proj, full.names =T)
#   ,
#   ~readr::read_csv(.x,col_select = all_of(focal_cols_recs))
# )
# #write.csv(rmis_recs, "C:/Data/CoTC/SouthernFund/2019/ERAnalysis/rmis_recov_86_20.txt")
# 
# # get unique tag codes to query RMIS releases
# 
# unique.tc<- unique(rmis_recs$tag_code) #12937 unique codes
# # split into two groups, because RMIS does not accept more than 10000 records
# top<-length(unique.tc)
# unique.tc1<-unique.tc[1:6000]
# unique.tc2<-unique.tc[6001:top]
# setwd("C:/Data/CoTC/SouthernFund/2019/ERAnalysis/")
# 
# write(unique.tc1, "utc1.txt")
# write(unique.tc2, "utc2.txt")
# 
# focal_cols_rel <- c("tag_code_or_release_id", "brood_year", "release_location_name", "hatchery_location_name", "stock_location_name", "release_location_state", "release_location_rmis_region", "release_location_rmis_basin", "cwt_1st_mark", "cwt_1st_mark_count", "cwt_2nd_mark", "cwt_2nd_mark_count")
# 
# rmis_rel<-map_df(
#   c(
#     "https://www.rmis.org/reports/CSV6577.txt",
#     "https://www.rmis.org/reports/CSV7053.txt"
#   ),
#   #~read.csv(.x)
#   ~readr::read_csv(.x,col_select = all_of(focal_cols_rel))
# )
# 
# #add area aggregation column to releases rmis_rel$Region
# rmis_rel$Region<-rmis_rel$release_location_state
# attach(rmis_rel)
# rmis_rel$Region<-case_when(is.na(release_location_state) ~ "PS",
#                            release_location_state=="BC" ~ "BC",
#                            release_location_state=="AK" ~ "AK",
#                            
#                            release_location_rmis_region  =="CECR" | release_location_rmis_region  =="CRGN"|
#                              release_location_rmis_region  =="LOCR" | release_location_rmis_region  =="UPCR"~ "Col R",
#                            
#                            release_location_rmis_region  =="HOOD" | release_location_rmis_region  =="JUAN"|
#                              release_location_rmis_region  =="MPS" | release_location_rmis_region  =="NOWA"|
#                              release_location_rmis_region  =="NPS" | release_location_rmis_region  =="SKAG"|
#                              release_location_rmis_region  =="SPS" ~ "PS",
#                            
#                            release_location_rmis_region  =="GRAY" | release_location_rmis_region  =="NWC"|
#                              release_location_rmis_region  =="WAGN" | release_location_rmis_region  =="WILP" ~ "Coastal",
#                            
#                            release_location_rmis_region  =="SOOR" | release_location_rmis_region  =="KLTR"|
#                              release_location_rmis_region  =="NOCA"  ~ "SONCC",
#                            release_location_rmis_region  =="NOOR" | release_location_rmis_region  =="ORGN"
#                            ~ "OR,CA,SNAK",
#                            release_location_rmis_region  =="CECA" | release_location_rmis_region  =="SAFA"
#                            ~ "OR,CA,SNAK",
#                            release_location_rmis_region  =="SNAK" ~ "OR,CA,SNAK"
# )
# 
# #write.csv(rmis_rel, "RMIS_RelforUTC.txt")
# # join release and recovery table
# rmis_recs2<-left_join(rmis_recs,rmis_rel, by = c("tag_code" = "tag_code_or_release_id"))
# # eliminate unmarked
# rmis_recs3<-rmis_recs2[,c(-18,-19)]
# rmis_recs4<-rmis_recs3[rmis_recs3$cwt_1st_mark>4999,]
# rmis_recs4<-rmis_recs4[!is.na(rmis_recs4$run_year),]
# 
# # #plot by region and year
# # rmis_recs4|>
# #   count(run_year, Region)|>
# #   ggplot(aes(run_year, n,fill = Region)) +
# #   geom_col(show.legend = F) + 
# #   facet_wrap(~Region, ncol=1)
# # map to fisheries
# 
# fish_agg<-read.csv("fishaggmap_no_fw.txt")   #freshwater not mapped in file and hence excluded"    
# # loop through each row in fish_agg
# # find all the records in rmis_recs4 that have a matching location and gear code
# # assign FishID to these records
# rmis_recrel<-rmis_recs4%>%
#   add_column (FishID = 0, FishName = 0)
# 
# 
# x<-nrow(fish_agg)
# 
# for(i in 1:x) {
#   # if (i > 31) {
#   #   browser()
#   #   }
#   locc<-fish_agg[i,4]
#   gearc<-fish_agg[i,3]
#   FishID2<-fish_agg[i,1]
#   FishName2<-fish_agg[i,2]
#   i=i+1
# #find all the records in rmis_recrel with matching gear and location codes
# # update FishID and FishName fields to FishID and FishName
#   l<-nchar(locc)
#   if (gearc == 5){
#     rmis_recrel<-within(rmis_recrel,FishID[floor(fishery/10) ==gearc ]<-FishID2)
#   } else {
#     rmis_recrel<-within(rmis_recrel,FishID[floor(fishery/10) ==gearc & substr(recovery_location_code,start=1,stop=l) ==locc]<-FishID2)
#   }  
# }
# 
# #rm(list=ls())
# #populate FishName field through look-up in fish_agg table
# rmis_recrel$FishName <-fish_agg[match(rmis_recrel$FishID,fish_agg$Fish_ID),2]
# rmis_recrel_fish<-rmis_recrel[!rmis_recrel$FishID==18,]
# rmis_recrel_esc<<-rmis_recrel[rmis_recrel$FishID==18,]
# #recrel_test<-rmis_recrel[1:20000,c(4,6,8,19,20)]
# #write.csv(recrel_test,"RecrelTest.txt")
# 
# #plot number of fishery recoveries by year for each region as stacked fishery bar graph
# #plot number of escapement recoveries by year for each region
# 
# a<-rmis_recrel_fish|>
#    count(run_year,Region,FishName)|>
#    filter(!is.na(FishName))  
# plot1<-ggplot(a,aes(run_year,n,fill = FishName)) + geom_col() + 
#    ggtitle("Marine Fishery CWTs by Stock Origin and Year")+
#    labs(y="CWT_Fisheries")+
#    facet_wrap(~Region, ncol = 1)+
#    theme(legend.position="none")
# 
# b<-rmis_recrel_esc|>
#   count(run_year,Region)
# plot2<-ggplot(b,aes(run_year,n, fill = "darkred")) + geom_col() + 
#     theme(legend.position = "none")+
#     ggtitle("Escapement CWTs by Stock Origin and Year")+
#     labs(y="CWT_Escapement")+
#     facet_wrap(~Region, ncol = 1)
# #arrange escapement and fisheries plots side-by-side
# require(gridExtra)
# plot3<-grid.arrange(plot1, plot2, ncol=2)
# ggsave("cwtsovertime.png", plot=plot3, height=7, width=10)
```

\newpage

# Data Preparation

## Assigning Coded-Wire-Tagged Groups to [FRAM](#FRAM) Model Stocks {#SM}

This section describes the steps taken to identify and map available Coho salmon [CWTs](#CWT) to the stocks in [FRAM](#FRAM). All available Coho salmon release and recovery data were downloaded from the Regional Mark Processing Center [(RMIS)](https://www.rmis.org/cgi-bin/queryfrm.mpl?Table=catch_sample&Version=4.1) in July of 2020 and uploaded into a Microsoft Access Database called 'MappingDatabase.mdb'. A total of 17,268 CWT groups were downloaded and ranged from brood year 1969 to 2019. After the raw data were loaded, the following steps were taken to create the stock assignments:

1. Mark Status was determined for each [CWT](#CWT) group released. Release and recovery data were examined to determine the intended mark status of each of the CWT release groups. For release groups identified as being marked, the estimated number of CWTs recovered from fisheries were summarized by 'recorded_mark' and compared to the mark type listed under 'cwt_1st_mark'. The CWT was flagged for further review if more than 20% of the fishery recoveries were identified as having a different mark than what was listed in the release data. Final review of all of these flagged CWTs yielded a total of 15 tag codes with the intended mark type differing from the mark type listed under 'cwt_1st_mark' (see Appendix Table [8.4](#WrongClip) for the list of these tags). A table called 'CWT_1st_mark_fixes' was created in the Access database and a query titled 'Step1- Fix CWT_1st_Mark value' was used to update the release data that was downloaded from [RMIS](#RMIS) into the database to correct the 'cwt_1st_mark'.

2. Releases were grouped and assigned to Production Regions (PRs) and Management Units (MUs). To streamline assigning release groups to PRs and MUs, a list of unique combinations of the following data was generated and saved in a table titled '[RMIS](#RMIS)_UniqueReleaseCodes' within the Access database: 
 - `release_location_state`
 - `release_location_rmis_region`
 - `release_location_rmis_basin`
 - `release_location_name`
 - `hatchery_location_name`
 - `stock_location_name`<br>
This table, consisting of a total of 2,058 unique combinations, was then copied to a table titled 'RMIS_UniqueReleaseCodesMapped' and a column containing the assigned PR and another containing the MU were added.  Production Regions and Management Units were manually assigned to each unique record. An example of a unique combination within the 'RMIS_UniqueReleaseCodesMapped' table that was mapped to a PR and an MU is below. See Appendix [Table 8.6](#StkMap) for a list of mapped [CWT](#CWT) codes and Appendix [Table 8.5](#CoPRMU) for a list of Coho Production Regions and Management Units[.]{#Unique}
```{r t41, results ='asis'}
tibble(
  `release location state` = "AK",
  `release location rmis region` = "ALSR",
  `release location rmis basin` = "ALSRG",
  `hatchery location name` = "NA",
  `stock location_name` = "AKWE R",
  `release_location name` = "AKWE R 182-40",
  `PR` = "Northern Alaska Outside",
  `MU` = "Alaska Northern Outside Hat/Wild",
) |> 
  gt::gt() |> 
  gt::tab_header(title = "Table 3.1. Unique Release Codes Mapped")
```
<br>
3. Assign MUs to the Columbia River releases. The Columbia River Production Region consists of two types of Coho Salmon: South or Early returning stocks (also known as A-Type or fall type run) and North or Late returning stocks (also known as B-Type or late-fall run). While there is a field in [RMIS](#RMIS) that provides the run type, this data was missing for many of the releases. Over 1,500 [CWT](#CWT) releases within the Columbia River Basin did not have a run type associated with them. Because of this, a summary of these tag codes was provided to Eric Kinne and Jillian Cady of WDFW to review. They provided MU assignments for the CWTs with blank MUs and confirmed all the preliminary assignments for the remaining records. This resulted in the following assignments summarized by MU, and an Access query updated all the Columbia River CWTs with these newly assigned MUs[.]{#ColRStk}
```{r t42, results ='asis'}
tibble(
  `Management Unit` = c("Columbia River Early Hatchery", "Columbia River Late Hatchery", "Lower Col R Oregon Wild","Wash Early Wild","Wash Late Wild","Youngs Bay Hatchery", "Grand Total"),
  `Number of CWTs Released` = c(1430,594,9,114,25,836,3008)
) |> 
  gt::gt() |> 
  gt::cols_label(.list = list(`Management Unit` = gt::md("**Management Unit**"), 
                              `Number of CWTs Released` = gt::md("**Number of CWTs Released**")))|>
  gt::tab_style(
    style = gt::cell_text(weight = "bold"),  
    locations = gt::cells_body(columns = everything(), rows=7))|>                        
  gt::tab_header(title = "Table 3.2. Columbia River Stock Assignments and Numbers Released")
```
<br>
4. Assigned PRs and MUs were used to map [coded wire tag](#CWT) codes to [FRAM](#FRAM) stocks through a look-up table with the following fields[:]{#StkCWT}
<br>
```{r t43, results ='asis'}
tibble(
  `Run ID` = 9,
  `PR Number` = 16,
  `PR MU Number` = 0,
  `CWT Code` = 011044,
  `Number Released` = 1000,
  `Release name` = "BURNT HILL CR (OPSR)",
  `Hatchery name` = "ORE-PAC SALMON RANCH",
  `Stock Name` = "NESTUCCA R PRIVATE",
  `brood year` = 1995
) |> 
  gt::gt() |> 
  gt::tab_header(title = "Table 3.3. Stock CWT 1998")
```
<br>
The look-up file to map [RMIS](#RMIS) recoveries to [FRAM](#FRAM) stocks is included on the [GitHub site](https://github.com/PSC-CoTC/coho_CWT_analysis/blob/main/lu_cct_MSM_Stock_CWT_20220913.csv) for this project. 

## Assigning Coded-Wire-Tagged Recoveries to [FRAM](#FRAM) Fisheries {#FM}

Coho [FRAM](#FRAM) contains [198 fisheries](#F) for complete West-Coast-wide Coho fishery coverage. [FRAM](#FRAM) fisheries are assigned by geographic areas and gear type. Generally, geographic fishing areas overlap with important management boundaries. Fisheries in the same geographic area are usually divided into troll, net, and sport fisheries. Frequently, treaty (first nations) and non-treaty (all citizen) fisheries are also separated into discrete fisheries. An assignment of [RMIS](#RMIS) [CWT](#CWT) recovery information to [FRAM](#FRAM) fisheries is usually possible by examining [RMIS](#RMIS)' `fishery` and `recovery_location_code` fields. [RMIS](#RMIS) `fishery` coding specifies the gear type, while the `recovery_location_code` designates a geographic area. [RMIS](#RMIS) contains 90 unique `fishery` codes and approximately 5,000 unique `recovery_location_codes`, potentially designating tens of thousands of unique fisheries. A 19-character location code strip is used by all parties to the [Pacific Salmon Treaty](#PST) to geographically identify locations (see [Lapi et al. 1989](#Lapi1)). Assignments of [RMIS](#RMIS) recoveries to [FRAM](#FRAM) fisheries are accomplished through a mix of a look-up on 'fishery' and 'recovery_location_code' fields and some coded logic. The coded logic relies on the hierarchical nature of the `recovery_location_code`, where the location code string narrows down the geographic area going from left to right in the coded sequence, e.g., the first character specifies the state or country, the second freshwater or marine, the third the management region, etc. When developing the look-up, the analyst determines the minimum number of `recovery_location_code` characters needed to make a distinct [FRAM](#FRAM) location assignment. A similar approach is taken to make [FRAM](#FRAM) gear assignments from the [RMIS](#RMIS) `fishery` field, e.g., a `fishery` code between 40 and 49 identifies a sport fishery. Recoveries with a `fishery` code between 50 and 59 are automatically assigned to escapement. Thus, the final fishery look-up table contains only 537 records. This is possible because the R code developed to map to FRAM fisheries truncates the `recovery_location_code` sequence to just the characters needed to make a unique [FRAM](#FRAM) fishery assignment, significantly reducing the number of possible permutations. The look-up file is included on the [GitHub site](https://github.com/PSC-CoTC/coho_CWT_analysis/blob/main/lu_ahb_FishMap_20220826.csv) for this project. 

## Data Issues

Despite the dedication of staff involved in data collection and management across organizations, the sheer number of records involved in constructing and populating a coast-wide model of salmon harvest means that problems will arise. The following briefly describes some of the major issues related to data acquisition and preparation that were encountered during this project. Some concerns cannot be remedied (e.g., past year fisheries cannot be retrospectively sampled), but further collective review efforts may improve the quality of key [FRAM](#FRAM) datasets and thereby any inferences based on model runs. In particular, we observed a need to ensure recent time-series of post-season runs reflect current catch and escapement estimates, some of which appear to have changed since being originally compiled.

Reported catch values are essential inputs to post-season [FRAM](#FRAM) and are critical to computing [CWT](#CWT) expansions. This project was impacted by inaccurate, uncertain, and missing values in catch reporting and accounting, and by extremely low or no sampling in some fisheries. In several instances, the modeled landed catch within [FRAM](#FRAM) did not match the catch reported to [RMIS](#RMIS)' Catch/Sample data, and mismatches persisted even after fisheries mapping issues were ruled out. Further investigation showed that official catch estimates were updated in [RMIS](#RMIS) after the [CoTC](#CoTC) completes the annual [ER](#ER) analyses using [BkFRAM](#BkFRAM). Therefore, we suggest annually re-querying and updating the most current 2-3 years to populate [BkFRAM](#BkFRAM). We also suggest similar procedures for catches reported to [RMIS](#RMIS).

Just as catches, stock specific escapements are also essential inputs to post-season [FRAM](#FRAM) to reconstruct pre-fishing cohort abundances. Missing and uncertain estimates of escapements related to lack of sampling and/or reporting hindered the number of ER comparisons that were possible. Another source of escapement error originates from inputs to [BkFRAM](#BkFRAM) runs. The [BkFRAM](#BkFRAM) option to apply [mark rates](#MR) (often from pre-season forecasts) to total (marked + unmarked) escapements has the potential to create inaccuracies in stock-specific and potentially whole fishery [mark rates](#MR), such as those assessed in this project. We recommend further review of [mark rate](#MR) estimates by regional experts, with special attention to instances where post-season escapement estimates are not provided by mark status. 

Accurate sampling data is critical for the integrity of the [CWT](#CWT) system, execution of mark-selective fisheries, production of [mark rate](#MR) and catch estimates, and collection of biological information. Hundreds of dedicated samplers annually collect these data, often under difficult circumstances. Several factors can lead to inaccurate results such as sampler error, equipment failure, sampling stratification that differs from fishery stratification, low sampling rates, etc. We encountered many instances where sampling was inadequate, relied on voluntary tag recoveries, or was conducted without electronic CWT detection equipment. Additionally, some fisheries, even if well sampled, exhibited highly improbable mark or tag rates.<br>

Since sampling is the backbone of collecting CWT information which flows into [exploitation rate](#ER) estimates, poor sampling can have a very detrimental effect on this type of analysis.<br>
Perhaps a feedback loop between RMIS and sampling programs could address some of these shortcomings by providing sampling supervisors with annual [RMIS](#RMIS) summary reports of sampling rates, CWT recoveries, and [mark rates](#MR). 

### [RMIS Data](#RMIS)

This project relied heavily on [RMIS](#RMIS) which houses a vast amount of information associated with [CWTs](#CWT). The complexity of these data creates considerable risk of misinterpretation, and this project underscored the importance of [RMIS](#RMIS) providing explanatory materials and training sessions to reduce this risk.

[RMIS](#RMIS) release data contains information about all salmon released from hatcheries as well as wild fish that were caught and marked or tagged prior to release. The [exploitation rate](#ER) analyses were conducted solely on [adipose-fin-clipped (marked)](#M) and [coded-wire-tagged](#CWT) fish. However, some of the salmon with tag codes from marked releases were assessed to be predominantly unmarked by field samplers. Either the tag codes were assigned the wrong mark status or the adipose clip was insufficient. Other problems encountered in the release data were optional fields left blank that would have allowed us to much more easily assign [CWTs](#CWT) to [FRAM](#FRAM) stocks. For example, the optional and often blank "Run" data field was needed to easily designate a Columbia River or WA Coast Coho stock as a "Summer" or "Fall" run or "Early" versus "Late".

Every salmon recovered with a [coded-wire tag](#CWT) (and subsequently decoded and reported) is a data record in the recovery database. This database contains many fields with fairly complex flagging options. There are a number of pitfalls that can lead an uninformed user to the wrong conclusions. For the database containing recovery data in particular, user training is highly recommended. In addition to the complexity of the database, deficiencies encountered in the recovery data were often caused by sampling issues. The principal concerns with this database are problematic [CWT](#CWT) expansions (see next bullet) and/or missing recoveries. Some entities fail to report CWT recoveries to [RMIS](#RMIS), sometimes resulting in entire fisheries and escapements missing from the database. In the case of missing escapements this can leave a large amount of uncertainty when conducting cohort reconstructions or [exploitation rate](#ER) analyses. We also found many instances of errors in recovery data fields, such as the `adclip_selective_fishery` field, where mark-selective fisheries (S) were incorrectly flagged as non-selective (N).

The Catch/Sample data is vital to the recovery data as [coded-wire-tag](#CWT) sampling expansions are calculated from the data. Inaccurately or non-reported catches or escapements can have a huge impact on [CWT](#CWT) analyses. Additionally, the stratification used in [RMIS](#RMIS) can be problematic. Generally, [RMIS](#RMIS) data is stratified by management or statistical weeks and months, while fisheries are commonly scheduled on a calendar month basis. This can result in miss-alignment of strata. Moreover, un-sampled strata can lead to missing CWT expansions (see [Handling of RMIS Catch/Sample data](#Handling)), and hence missing expanded recoveries.   

### [MSM](#MSM) User Guide and Code Documentation

The [Mixed Stock Model](#MSM) was rewritten in Visual Basic in 2007. Although the document ["Coho FRAM Base Development"](#Pac1) describes the project, the current [MSM](#MSM) lacks code documentation and a user guide, thus requiring a great deal of initiative by a new user to understand the code, the data requirements, and re-(run) analyses. A draft [user guide](#MSMGuide) of the [MSM](#MSM) was developed during the early stages of this project.

### Stock Mapping

Mapping hatchery releases to [FRAM](#FRAM) stocks was one of the most timing consuming tasks of this project. Thorough descriptions of the Coho stocks in [FRAM](#FRAM) would be extremely beneficial for future stock mapping exercises as well as understanding the stocks in the model. We therefore suggest compiling detailed stock profiles, describing the locations, hatcheries, calibration [CWT](#CWT) codes, expansions, associated wild components, modeling parameters, [exploitation rate](#ER) profiles, etc. 

### Fisheries Mapping

The tasks of mapping [RMIS](#RMIS) "recovery location codes" and "fisheries" to [FRAM](#FRAM) fisheries was similarly labor intensive as the stock mapping. We followed the mapping approach from the original [MSM](#MSM) to allow future [MSM](#MSM) tasks to build on the mapping efforts of this project. While mapping to [RMIS](#RMIS) fisheries, we discovered that **the Southeast Alaska sport fishery is entirely missing from FRAM**. We were unable to determine the reason for this omission and urge additional investigations. As with the stock mapping, developing detailed [FRAM](#FRAM) fishery profiles is highly recommended. 


### Handling of Fisheries in Response to [RMIS](#RMIS) Catch/Sample data {#Handling}

[RMIS](#RMIS) Catch/Sample data requires carefully examined before analyzing the [coded-wire-tags](#CWT) associated with a sample. 

The Catch/Sample data is stratified into time periods ranging from a few days to several months depending on the sample design and sample size. 

Coded-wire tags are expanded by the ratio of reported catch to sampled catch within a stratum, i.e., the `Estimated Number` in the [CWT](#CWT) recoveries tables reflects the number of CWTs that would have been recovered if 100% of the landed catch in a fishery were sampled. 

Fisheries sampling and catch reporting can be fraught with many problems and errors such as missing sampling, missing catch reporting, recording the wrong catch area, combining samples from more than one fishery, equipment failure, etc. 

The following sections describe how the most significant errors were addressed.

One of the biggest sources of error in tag expansions is introduced when a stratum with catch is not sampled.

[RMIS](#RMIS) does not address this error by expanding strata into sampled time periods and recomputing sampling expansions. 

This results in underestimating the [CWTs](#CWT) that should have been recovered had the catch been sampled and consequently in [ERs](#ER) that are too low for affected time periods.

Approximately 30% of the 2010-20 Coho catch occurred in unsampled strata. 

Depending on the fishery and the magnitude of the problem the following approaches were taken:

#### Exclude Entire Fisheries

Fisheries that were never or rarely sampled were excluded from the [exploitation rate](#ER) analysis. 

We also eliminated fisheries that were not modeled in [FRAM](#FRAM). 

Examples of unmodeled fisheries were some [South of Falcon](#SOF) fisheries that report a small amount of likely illegal catch. 

Additionally, all freshwater fisheries were automatically excluded, because they were frequently not or insufficiently sampled. As these fisheries are predominately catching the local stock(s), they were of minor informational value for this analysis. 

Excluded fisheries were removed from [FRAM](#FRAM) as well as [RMIS](#RMIS).

[Table 3.4]{#ExcludeFish} contains a list of excluded fisheries and the principal reason for exclusion (several fisheries were excluded for more than one reason). 

```{r t44, results='asis'}
d_exclude_fish|>
  knitr::kable() |>
  kableExtra::add_header_above(data.frame("Table 3.4. Fisheries Excluded from CWT Analysis and Main Reason for Exclusion", 3), TRUE) |>
  kableExtra::kable_styling() |>
  kableExtra::scroll_box( width="100%",height = "250px")
```

#### Recompute the Sampling Expansion

To address unsampled strata, a coefficient was calculated to adjust annual [CWT](#CWT) recoveries to the number of recoveries that should have occurred if all the strata had been sampled. 

This coefficient is referred to as an "Expansion" in Appendix table [8.3](#CWTAdj). The coefficient is computed by summing the reported [RMIS](#RMIS) catch over all the strata within a year divided by the sum of the [RMIS](#RMIS) catch associated with a sampling event.

Appendix table [8.3](#CWTAdj) contains a summary of [FRAM](#FRAM) catches, [RMIS](#RMIS) catches, sampled catches and fishery expansions by fishery and year. The expansions are computed as below. Fisheries or years that are eliminated contain a zero in the "Expansion" row in table 8.3.<br>

Equation 3.1:
$$AdjustedCWT = EstimatedCWT * \frac {\sum_{strata}NumberCaught}{\sum_{strata}(NumberCaught*b)}$$

Where,    

*   b = 1 if the stratum was sampled and 0 if the stratum was unsampled
*   EstimatedCWT = original expanded for sampling [CWT](#CWT) recoveries
*   NumberCaught = Number Caught in [RMIS](#RMIS) CatchSample database

#### Exclude Years

Years were eliminated for the following reasons:

* [FRAM](#FRAM) catch less than or equal to 10
* [RMIS](#RMIS) catch less than or equal to 10
* Catch associated with a sampling event less than or equal to 10
* Large discrepancy between [RMIS](#RMIS) and [FRAM](#FRAM) catch. Generally, years were excluded when [RMIS](#RMIS) and [FRAM](#FRAM) catch differed by more than 9.5% and 400 fish. 
* Adjustment to sampling expansion greater than 8.
Years were removed from [RMIS](#RMIS) as well as [FRAM](#FRAM).

#### Combine Fisheries

Some fisheries are conducted in the freshwater as well as the adjacent marine [terminal](#Term) area. The catch or the sampling data for these areas may be reported in one or the other area interchangeably. This was the case for Bellingham Bay and Samish River net (fisheries 97 & 98). Once combined, the [FRAM](#FRAM) catch closely matched the [RMIS](#RMIS) catch.

#### Special Cases

**Area 5 sport (fishery 91):** 
For Area 5 sport neither the [FRAM](#FRAM) modeled catch nor the [RMIS](#RMIS) catch matched the official Catch Record Card [(CRC)](#CRC) estimates. This was likely caused by CRC revisions. However, all strata of this fishery were sampled each year, eliminating the need for a sampling adjustment. For this fishery the [CWTs](#CWT) were simply multiplied by the ratio of "FRAM Catch" to "RMIS Catch".  

**Southeast Alaska net (fishery 198):**
The average annual [FRAM](#FRAM) catch for southeast Alaska of 648,000 was significantly smaller than the [RMIS](#RMIS) catch of 784,000. Fishery mapping issues were ruled out as a source of the discrepancy. Since the catch associated with sampling was similar to the [FRAM](#FRAM) catch each year (average 644,000), the [CWTs](#CWT) were not adjusted.

[**Columbia River Net and Sport (fisheries 24-32):**]{#CR}
Columbia River stock abundances (StockID 165-176) enter [FRAM](#FRAM) in units of "Return to the River". Pre-fishing abundances are calculated by adding marine fisheries mortalities to this abundance. In-river fisheries are not modeled in a [FRAM](#FRAM) forward run. Therefore, the resulting [FRAM](#FRAM) escapement is actually the abundance entering the Columbia River, prior to in-river sport and net fisheries. When comparing [FRAM](#FRAM) [exploitation rates](#ER) to those derived using [RMIS](#RMIS) [CWT](#CWT) recoveries, [RMIS](#RMIS) recoveries from in-river fisheries are added to escapement to reconstruct a "Return to the River" run size that is in the same units as [FRAM](#FRAM).

\newpage

# Sensitivity of [Exploitation Rates](#ER) of Key Stocks to Abundance Changes

## Introduction

Many West Coast Coho stocks have overlapping migration paths or natal streams in close proximity to one another. They can therefore be intercepted in the same mixed stock fisheries, where the abundance of one stock can have a pronounced effect on the possibility to encounter another stock.

The goal of this project is to evaluate the effects of abundance changes of influential Coho stocks on the [exploitation rate](#ER) of other model stocks to determine which forecasts of abundance have the greatest effect on key wild stocks. Results of this analysis can provide decision makers with valuable information about prioritizing the timing and availability of forecast information and allocating limited stock assessment resources.

## Methods

For this analysis the abundances of stock aggregates (test stocks) in the [FRAM](#FRAM) model was adjusted one aggregate at a time in order to evaluate the resulting mortalities of key wild model stocks (response stocks). For more information about the [FRAM](#FRAM) model please see [FRAM documentation](https://framverse.github.io/fram_doc/).

The Coho Technical Committee [(CoTC)](#CoTC) annually produces post-season Coho model runs. These runs start with the final pre-season model run and are subsequently updated with best estimates of stock-specific escapements and post-season assessments of landed and non-landed fisheries mortalities. The most current database is housed on the Pacific Salmon Commission's [(PSC)](#PSC) [Extranet website](https://extranet.psc.org/SitePages/Home.aspx). Model runs are currently available through 2020 and stored in Access database "PSC_CoTC_PostSeason_CohoFRAMDB_thru2020_021622.mdb". 

This evaluation was conducted using post-season Coho model runs [(BkFRAM runs)](#BkFRAM) by adjusting stock specific starting cohorts. Post-season runs model catches at fixed values regardless of the abundance of model stocks. Each stock aggregate was modeled at a 'low' and a 'high' abundance (see [Selection of Low and High Stock Abundances](#HighLow)) as well as in a 'low' and 'high' abundance run year (see [Selection of Fishing Years](#YrSelect)). Examining two run years, doubled the number of model runs, but added valuable information about the sensitivity of abundance effects to fishing conditions and run sizes during generally low (2015) versus high (2013) abundance years.

### Selection of Stocks

The Coho [FRAM](#FRAM) model contains all West Coast [stocks]($S) from Alaska to California. For the purpose of this project the influence of abundance changes of pre-selected stock aggregates (test stocks) on wild stocks (response stocks) managed under the [PST](#PST) was evaluated. 

[**Table 4.1. List of Response Stock Aggregates**]{#ResponseStk}
<br>
A production region (PR) is comprised of one to several stocks with similar migration paths and [exploitation rates](#ER).

**Management Unit**|**Short Name**|**Production Region**|**FRAM Stock Name**|**FRAM Stock ID**
:---|:---|:---:|:---|:---
Skagit|Skagit|2|Skagit, Baker|17,23
Stillaguamish/Snohomish|Stilly/Sno|3|Stillaguamish, Snohomish|29, 35
Hood Canal|HC|4|12/12B, 12C/D, Skokomish|45,55,59
Strait of Juan de Fuca|JDF|6|East JDF Misc., West JDF Misc.|115,117
Quillayute|Quill|8|Quillayute|131
Hoh|Hoh|9|Hoh|135
Queets|Queets|10|Queets|139
Grays Harbor|GH|12|Chehalis, Humptulips, Misc. Wild|149,153,157
Georgia Strait Mainland|GSML|17|Georgia Strait Mainland|207
Georgia Strait Vancouver Island|GSVI|18|Georgia Strait Vancouver Island|211
Lower Fraser|LwFraser|22|Lower Fraser|227
Upper Fraser|UpFraser|23|Upper Fraser|231

Test stocks were selected based on origin and the likelihood of affecting response stocks. 

[**Table 4.2. List of Test Stock Aggregates**]{#TestStk}

```{r t52, echo = FALSE, results='asis'}
# **Management Unit**|**Short Name**|**Production Region**|**FRAM Name**|**FRAM ID**
# :---|:---|:---:|:---|:---:
# Nooksack/Samish|NOSA|1|Nooksack Wild, Kendall Hat, Skookum Hat, Lummi Ponds, Bellingham Bay Net Pens, Samish Wild, 7-7A Independent Wild, Whatcom Hat|1-16
# Skagit|Skagit|2|Skagit Wild, Skagit Hat, Baker Wild, Baker Hat, Swinomish Hat, Oak Harbor Net Pens|17-28
# Stillaguamish/Snohomish|Stilly/Sno|3|Stilly Hat, Stilly Wild, Tulalip Hat, Sno Wild, Sno Hat, 8A Net Pens|29-40
# Hood Canal|HC|4|Port Gamble Net Pens, Port Gamble Wild, 12/12B Wild, Quilcene Hat, Quilcene Net Pens, 12A wild, Hoodsport Hat, 12C/12D Wild, George Adams Hat, Skokomish Wild|41-60
# South Puget Sound|SPS|5|13B Misc Wild, Deschutes, SPS Net Pens, Nisqually Hat, Nisqually Wild, Fox Is. Net Pens, Minter Hat, 13 Misc Wild, 13 Misc Wild, 13Misc Hat, Chambers Hat, 13A Misc Wild, Puyallup Hat, Puyallup Wild, 11 Hat, 11 Misc Wild, 10E Hat, 10E Misc Wild, Green Hat, Green Wild, Lk Washington Hat, Lk Washington Wild, 10 Hat, 10 Misc Wild|61-106
# Strait of Juan de Fuca|JDF|6|Elwah Wild, Dungeness Wild, East JDF Misc Wild, West JDF Misc. Wild, Port Angeles Net Pens, 9 Misc Wild|107-122
# Quillayute|Quill|8|Quillayute Summer Wild, Quillayute Summer Hat,Quillayute Fall Wild,Quillayute Fall Hat|127-134
# Hoh|Hoh|9|Hoh Wild, Hoh Hat|135-138
# Queets|Queets|10|Queets Fall Wild, Queets Fall Hat, Queets Suppl. Hat|139-144
# Grays Harbor|GH|12|Chehalis Wild, Chehalis Hat, Humptulips Wild, Humptulips Hat, GH Misc Wild, GH Net Pens|149-160
# Georgia Strait Mainland|GSML|17|Georgia Strait Mainland Hatchery and Wild|205-208
# Georgia Strait Vancouver Island|GSVI|18|Georgia Strait Vancouver Island Hatchery and Wild|209-212
# SW Vancouver Island|SWVI|20|SW Vancouver Island Hatchery and Wild|217-220
# Lower Fraser|LwFraser|22|LwFraser Hat, LwFraser Wild|225-228
# Upper Fraser|UpFraser|23|UpFraser Hat, UpFraser  Wild|229-232

table2_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_2_section1c.xlsx"))

reactable::reactable(table2_1c,
                     bordered = TRUE,
                     striped = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     highlight = TRUE,
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       align = "center", 
                                                       headerStyle = list(background = "#f7f7f8")),
                     columns = list(
                       `Management Unit` = reactable::colDef(minWidth = 120),
                       `Short Name` = reactable::colDef(minWidth = 80),
                       `Production Region` = reactable::colDef(minWidth = 95),
                       `FRAM Name` = reactable::colDef(minWidth = 350),
                       `FRAM ID` = reactable::colDef(minWidth = 80)
                     )
                    )
```
<br>

### [Selection of Fishing Years]{#YrSelect}

This project focused on recent 2010-2018 fishing years, because recent years are most pertinent to current environmental conditions and because this time series had received additional review. 2013 was selected to represent a 'high' Coho abundance year and 2015 to represent a 'low' Coho abundance year. The aggregate starting cohort abundances for all test stocks combined was 2,240,000 million and 760,000 respectively. The remaining years were solely used to calculate average starting abundances, but not to run the actual simulations (see description under [Model Runs](#ModelRuns)).

### [Selection of Low and High Stock Abundances]{#HighLow}

For each test stock aggregate (Table 4.2) the low abundance was set at 30% and the high abundance at 170% of the 2010-2018 average [starting cohort](#SC). Starting cohorts are abundances pre fishing and pre natural mortality. Lows and highs were chosen to be equidistant from the average. This limited the range, as abundances below 30% of average could result in negative escapements given the fisheries in the 2013 and 2015 model runs.

[**Table 4.3. 2010 to 2018 Starting Cohorts of Coho Test Stocks**]{#StartingCohorts}<br>
Table includes average, low (30% of average) and high (170% of average) starting cohorts for each test stock. In [FRAM](#FRAM) abundances are modeled using stock recruit scalers. These scalers linearly relate the abundance in a given year to the abundance in [FRAM](#FRAM) [base period](#BP) years. 

```{r t53, echo = FALSE, results = 'asis'}
table3_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_3_section1c.xlsx"))

reactable::reactable(table3_1c,
                     bordered = TRUE,
                     striped = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     highlight = TRUE,
                     columns = list(
                       PR = reactable::colDef(minWidth = 40, sticky = "left"),
                       `Stock ID` = reactable::colDef(minWidth = 60, sticky = "left"),
                       `2010` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2011` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2012` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2013` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2014` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2015` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2016` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2017` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       `2018` = reactable::colDef(minWidth = 70,
                                                  format = reactable::colFormat(separators = TRUE)),
                       Average = reactable::colDef(minWidth = 75,
                                                   style = list(background = "#ADE6D1"),
                                                   format = reactable::colFormat(separators = TRUE)),
                       `Low (30%)` = reactable::colDef(minWidth = 70,
                                                       style = list(background = "#E6BFAD"),
                                                       format = reactable::colFormat(separators = TRUE)),
                       `High (170%)` = reactable::colDef(minWidth = 70,
                                                         style = list(background = "#E6BFAD"),
                                                         format = reactable::colFormat(separators = TRUE)),
                       `Low Stock Recruit` = reactable::colDef(minWidth = 70,
                                                               style = list(background = "#ADD8E6"),
                                                               format = reactable::colFormat(digits = 3)),
                       `High Stock Recruit` = reactable::colDef(minWidth = 70,
                                                                style = list(background = "#ADD8E6"),
                                                                format = reactable::colFormat(digits = 3))
                     ),
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       align = "center", 
                                                       headerStyle = list(background = "#f7f7f8"))
                    )
```

### [Model Runs]{#ModelRuns}

Abundances of stock aggregates were adjusted one production region at a time. Each stock within a production region was modeled at a low (30% of average) and at a high (170% of average) abundance. This was done for the 2013 as well as the 2015 post-season runs, resulting in a total of 64 model runs (16 production regions * 2 years * 2 abundance levels). Appendix [table 8.8](#RunInfo) contains a list of model runs.

To avoid negative escapements, [extreme terminal](#ET) and freshwater fisheries were set to zero. These fisheries predominantly intercept local stocks with minimal interference from other test stocks; thus, providing little relevance to the sensitivity analysis.

Appendix [table 8.9](#FZero) contains a list of fisheries set to zero in the model runs.

Additionally, several fisheries were modified for individual run pairs (low/high pairs) to avoid negative escapements (see Appendix [table 8.13](#NegEsc)).

### Analysis of Results

After running the [FRAM](#FRAM) model, total fisheries related mortalities and escapements from response stocks were extracted from the Access database to compute [exploitation rates](#ER). For each model run, mortalities and escapements were summed over all response stocks within a production region (PR).

Equation 4.1:

$$ExploitationRate_{PR} = \frac{\sum_{stk \in PR}Mortalities}{\sum_{stk \in PR}Mortalities + \sum_{stk \in PR}Escapement}$$

This resulted in four exploitation rates (2013 Low, 2013 High, 2015 Low, 2015 High) per production region.

## Results

The magnitude of the [exploitation rate](#ER) difference between low (abundance at 30% of average) versus high (abundance at 170% of average) test stock abundances, henceforth referred to as the 'effect', is a function of two factors, (1) the test stock abundance and (2) overlapping fisheries distributions. A test stock aggregate with a high abundance will affect the response stock to a greater degree than one with a low abundance. If the test stock and result stock have a similar distribution in key fisheries, the effect will also be greater.

Equation 4.2: The "effect is calculated as:

$$ResultStock.ER_{TestStockLow} - ResultStock.ER_{TestStockHigh}$$

Results are displayed in three tables (Tables 4.5-4.7) and two graphs (Figures 4.1-4.2).

Table 4.5 displays the [exploitation rates](#ER) of key wild response stock resulting from varying the abundances of test stocks. Results are presented for every test/response stock combination in the 2013 and 2015 post-season runs. Exploitation rates from table 4.5 are graphed in Figure 4.1 with different color bars for 'Low' and 'High' test stock abundances for the 2013 and 2015 run years (4 bars per response stocks). Each test stock is displayed as a separate panel in the graph. Table 4.6 summarizes detailed results from table 4.5 by presenting only exploitation rate differences (as calculated in equation 4.2). Differences are averaged over 2013 and 2015 run years. Figure 4.2 displays exploitation rate differences from table 4.5 with different color bars for the 2013 and 2015 run years (2 bars per response stock). Each test stock is displayed as a separate panel in the graph. Table 4.7 lists the top five fisheries responsible for exploitation rate differences of notable test/response stock combinations, i.e., test stocks with a significant or otherwise notable impact on response stocks.

[**Table 4.5. [Exploitation Rates](#ER) of Key Wild Response Stock Resulting from Varying Abundance of Test Stocks**]{#ERKey}<br>
“Low” denotes test stock abundance at 30% of average and “High” at 170% of average. Run year 2013 represents a high abundance year, and 2015 a low abundance year. The columns “Diff” display the difference in response stock exploitation rates resulting from varying the abundance of test stocks in the 2013 and 2015 model run. The column 'Avg Diff' averages 2013 and 2015 differences.

```{r t56, echo = FALSE, results = 'asis'}
table6_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_6_section1c.xlsx"), range = "A2:K194")

reactable::reactable(table6_1c |> 
                     dplyr::rename(`Stock ID` = `Stk ID`),
                     bordered = TRUE,
                     striped = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     highlight = TRUE,
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       align = "center",
                                                       headerStyle = list(background = "#f7f7f8")),
                     columnGroups = list(
                       reactable::colGroup(name = "Test Stock", 
                                           columns = c("PR #", "PR Name")),
                       reactable::colGroup(name = "Response Stock", 
                                           columns = c("Stock ID", "Stock Name")),
                       reactable::colGroup(name = "Exploitation Rates", 
                                           columns = c("2013 Low", "2013 High",
                                                       "2015 Low", "2015 High",
                                                       "2013 Diff", "2015 Diff",
                                                       "Avg Diff"))
                     ),
                     columns = list(
                       `PR #` = reactable::colDef(minWidth = 70, filterable = TRUE, sticky = "left"),
                       `PR Name` = reactable::colDef(minWidth = 120, filterable = TRUE, sticky = "left"),
                       `Stock ID` = reactable::colDef(minWidth = 85, filterable = TRUE, sticky = "left"),
                       `Stock Name` = reactable::colDef(minWidth = 130, filterable = TRUE, sticky = "left"),
                       `2013 Low` = reactable::colDef(minWidth = 95, sortable = TRUE,
                                                      format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `2013 High` = reactable::colDef(minWidth = 97, sortable = TRUE,
                                                      format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `2015 Low` = reactable::colDef(minWidth = 95, sortable = TRUE,
                                                      format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `2015 High` = reactable::colDef(minWidth = 97, sortable = TRUE,
                                                       format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `2013 Diff` = reactable::colDef(minWidth = 95, sortable = TRUE,
                                                       style = list(background = "#D3D3D3"),
                                                       format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `2015 Diff` = reactable::colDef(minWidth = 95, sortable = TRUE,
                                                       style = list(background = "#D3D3D3"),
                                                       format = reactable::colFormat(percent = TRUE, digits = 2)),
                       `Avg Diff` = reactable::colDef(minWidth = 95, sortable = TRUE,
                                                      style = list(background = "#A9A9A9"),
                                                      format = reactable::colFormat(percent = TRUE, digits = 2))
                     ),
                     defaultSorted = list(`PR #` = "asc"),
                     rowStyle = reactable::JS("
                      function(rowInfo, state) {
                        if (!rowInfo) return
                        const firstSorted = state.sorted[0]
                        if (firstSorted && firstSorted.id === 'PR #') {
                          const nextRow = state.pageRows[rowInfo.viewIndex + 1]
                          if (nextRow && rowInfo.values['PR #'] !== nextRow['PR #']) {
                            return { boxShadow: 'inset 0 -2px 0 rgba(0, 0, 0, 0.1)' }
                          }
                        }
                      }
                    ")
)
```

**Table 4.6. Average [Exploitation Rate](#ER) Differences of Key Wild Response Stocks Resulting from Varying the Abundances of Test [Stocks]{#ERDiff}**<br>

The red bars display the relative magnitude of the effect on 'foreign' stocks, whereas blue bars show the relative magnitude of the effect on 'home' stocks belonging to the same production unit as the test stock. Green shading is used if test and response stock belong to the same country, yellow for US test stocks on Canadian response stocks or vice versa. The yellow bars on the last row index production region abundances to the largest (South Puget Sound (SPS)) stock, i.e. average 2010-18 abundance of production region divided by average SPS production region abundance.

```{r t57, echo = FALSE, message = FALSE, warning = FALSE, results = 'asis'}
table7_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_7_section1c.xlsx"))

table7_1c %>% 
    dplyr::mutate(
        nooksambarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          TRUE ~ "#FF0000A0"
          ),
        skagitbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("skagit", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        stillysnobarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("stilly", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        hcbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("hood canal", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        spsbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          TRUE ~ "#FF0000A0"
          ),
        jdfbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("jdf", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        quillbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("quill", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        hohbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("hoh", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        queetsbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("queets", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        quinbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          TRUE ~ "#FF0000F0"
          ),
        ghbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("gh", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        gsmlbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("gsml", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        gsvibarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("gsvi", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        swvibarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          TRUE ~ "#FF0000A0"
          ),
        lwfrasbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("lw fras", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
        upfrasbarcol = dplyr::case_when(
          `Response Stock` == "Abundance Index" ~ "#FFD700",
          grepl("up fras", tolower(`Response Stock`)) ~ "#0000FF80",
          TRUE ~ "#FF0000A0"
          ),
    ) %>% reactable::reactable(., 
                   defaultPageSize = 13,
                   highlight = TRUE,
                   bordered = TRUE,
                   defaultColDef = reactable::colDef(
                       align = "center",
                       headerStyle = list(background = "#f7f7f8")
                   ),
                   columns = list(`Response Stock` = reactable::colDef(minWidth = 145,
                                                                       align = "left",
                                                                       sticky = "left",
                                                                       style = function(value, index){
                                                                         if(value == "Abundance Index"){
                                                                           list(fontWeight = "bold")
                                                                         }
                                                                       }),
                                    `Nook Samish` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "nooksambarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    Skagit = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                    },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "skagitbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120), 
                                    `Stilly Sno` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "stillysnobarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    `Hood Canal` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "hcbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    SPS = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "spsbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    JDF = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "jdfbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    Quill = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "quillbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    Hoh = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "hohbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    Queets = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "queetsbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    Quinault = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "quinbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    `Grays Harbor` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#8FBC8F")  
                                        } else if(index %in% 9:12){
                                          list(background = "#FFFACD")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "ghbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    GSML = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#FFFACD")  
                                        } else if(index %in% 9:12){
                                          list(background = "#8FBC8F")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "gsmlbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    GSVI = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#FFFACD")  
                                        } else if(index %in% 9:12){
                                          list(background = "#8FBC8F")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "gsvibarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    SWVI = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#FFFACD")  
                                        } else if(index %in% 9:12){
                                          list(background = "#8FBC8F")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "swvibarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    `Lw Fraser` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#FFFACD")  
                                        } else if(index %in% 9:12){
                                          list(background = "#8FBC8F")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "lwfrasbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    `Up Fraser` = reactable::colDef(style = function(value, index){
                                        if (index < 9){
                                           list(background = "#FFFACD")  
                                        } else if(index %in% 9:12){
                                          list(background = "#8FBC8F")
                                        } else{
                                          list(background = "none")
                                        }
                                      },
                                      cell = reactablefmtr::data_bars(.,
                                                       text_position = "inside-base",
                                                       brighten_text = FALSE,
                                                       max_value = 1,
                                                       min_value = -0.01,
                                                       background = "#D3D3D380",
                                                       fill_color_ref = "upfrasbarcol",
                                                       number_fmt = scales::label_percent(accuracy = 0.1)),
                                      minWidth = 120),
                                    nooksambarcol = reactable::colDef(show = FALSE),
                                    skagitbarcol = reactable::colDef(show = FALSE),
                                    stillysnobarcol = reactable::colDef(show = FALSE),
                                    hcbarcol = reactable::colDef(show = FALSE),
                                    spsbarcol = reactable::colDef(show = FALSE),
                                    jdfbarcol = reactable::colDef(show = FALSE),
                                    quillbarcol = reactable::colDef(show = FALSE),
                                    hohbarcol = reactable::colDef(show = FALSE),
                                    queetsbarcol = reactable::colDef(show = FALSE),
                                    quinbarcol = reactable::colDef(show = FALSE),
                                    ghbarcol = reactable::colDef(show = FALSE),
                                    gsmlbarcol = reactable::colDef(show = FALSE),
                                    gsvibarcol = reactable::colDef(show = FALSE),
                                    swvibarcol = reactable::colDef(show = FALSE),
                                    lwfrasbarcol = reactable::colDef(show = FALSE),
                                    upfrasbarcol = reactable::colDef(show = FALSE)
                                  )
                                )
```

[**Figure 4.1. [Exploitation Rates](#ER) of Key Wild Stocks Resulting from Varying the Abundance of Test Stocks**]{#ERWild}<br>
Figure has separate panels for each test stock. Response stocks are graphed along the x-axis. The exploitation rate for each response stock is plotted for the four modeling scenarios (2013Low/High, 2015Low/High)

```{r f51}
SenA.Results=read.csv(paste0(table_dir,"SA_Results - Copy.csv"))
SenA.Results<-data.frame(SenA.Results)
SenA.Results1<-SenA.Results[,-c(9:11)]

SenA.Results1<-reshape::melt(SenA.Results1, id=c(1:2,7:8))
colnames(SenA.Results1)[5]<-"Run"
colnames(SenA.Results1)[6]<-"Expl.Rate"

p<-ggplot(SenA.Results1,aes(Result, Expl.Rate))
SenA.Graph1<-p+geom_bar(stat = "identity",aes(fill=Run), position="dodge")+
  xlab("Response Stock")+ylab("Expl.Rate")+ggtitle ("Exploitation Rate as Function of Abundance")+ 
  theme(axis.text.x = element_text(angle = 90),legend.position = "bottom")+ facet_wrap( ~ Adjusted.Abundance )+
  scale_fill_manual(values=c("bisque","burlywood1","cadetblue1","cadetblue3"),
  labels =c("2013 Low", "2013 High", "2015 Low", "2015 High"))
SenA.Graph1
#ggsave("SenA.Graph1.png",width=8, height=10)
```

[**Figure 4.2. [Exploitation Rate](#ER) Differences Resulting from Varying the Abundances of Test Stocks**]{#ERWildDiff} <br>
Figure has separate panels for each test stock. Response stocks are graphed along the x-axis. The [exploitation rate](#ER) differences for each response stock are plotted for 2013 and 2015 run years. Graph is excluding scenarios where test and response stocks are the same. Y-Axes are on the same scale for all panels. 

```{r f52}
SenA.Results2<-SenA.Results
SenA.Results2<-SenA.Results2[,-c(3:6,11)]
SenA.Results2<-reshape::melt(SenA.Results2, id=(1:4))
colnames(SenA.Results2)[6]<-"ER.Diff"
colnames(SenA.Results2)[5]<-"Run"
#Set rows to zero if response stock = adjusted.abundance stock to avoid unreadable y axis"
a<-index<-SenA.Results2$Stk1==SenA.Results2$Stk2
SenA.Results2$ER.Diff[index]<-0
SenA.Results2

q=ggplot(SenA.Results2,aes(Result, ER.Diff))
SenA.Graph2<-q+geom_bar(stat = "identity",aes(fill=Run), position="dodge")+
  xlab("Response Stock")+ylab("Expl.Rate Difference")+
  ggtitle ("Exploitation Rate Difference as Function of Abundance")+ 
  theme(axis.text.x = element_text(angle = 90), legend.position = "bottom")+ facet_wrap( ~ Adjusted.Abundance)+
  scale_fill_manual(values=c("burlywood1","cadetblue3"),
                    labels =c("2013 Diff", "2015 Diff"))
SenA.Graph2

#ggsave("SenA.Graph3.png",width=8, height=10)
```

[**Table 4.7. Top 5 Fisheries Responsible for Total Mortality Differences Caused by Varying the Abundances of Test Stocks**]{#Top5} <br> 
Table only displays difference for test/response stock combinations with notable differences. The difference is calculated as the average 2013 and 2015 difference in mortalities between low and high abundances of test stocks.

```{r t58, echo = FALSE, message = FALSE, warning = FALSE, results='asis'}
table8_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_8_section1c.xlsx"))

reactable::reactable(table8_1c,
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       align = "center",
                                                       minWidth = 150,
                                                       headerStyle = list(background = "#f7f7f8")),
                     highlight = TRUE,
                     outlined = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     columns = list(
                      `Test Stock` = reactable::colDef( 
                        style = reactable::JS("function(rowInfo, column, state) {
                          const firstSorted = state.sorted[0]
                          if (!firstSorted || firstSorted.id === 'Test Stock') {
                            const prevRow = state.pageRows[rowInfo.viewIndex - 1]
                            if (prevRow && 
                              rowInfo.values['Test Stock'] === prevRow['Test Stock'] && 
                              rowInfo.values['Response Stock'] === prevRow['Response Stock']) {
                                return { visibility: 'hidden' }
                            }
                          }
                        }")
                      ),
                      `Response Stock` = reactable::colDef( 
                        style = reactable::JS("function(rowInfo, column, state) {
                          const firstSorted = state.sorted[0]
                          if (!firstSorted || firstSorted.id === 'Response Stock') {
                            const prevRow = state.pageRows[rowInfo.viewIndex - 1]
                            if (prevRow && 
                              rowInfo.values['Test Stock'] === prevRow['Test Stock'] && 
                              rowInfo.values['Response Stock'] === prevRow['Response Stock']) {
                                return { visibility: 'hidden' }
                            }
                          }
                        }")
                      )
                     )
)
```

## Discussion

Two run years were selected to evaluate the effects of abundance changes of individual stock aggregates (test stocks) on key wild stocks (response stocks). The abundance of the test stock aggregate varied between 0.3 and 1.7 of the average abundance. The [exploitation rate](#ER) difference was on average 2.6 times greater during the low abundance year. The nature of exploitation rate calculations and post-season model run inputs results in unidirectional exploitation rate changes, i.e., exploitation rates of response stocks are always higher when the test stock is at low abundance. For almost all scenarios the exploitation rate difference was greater for the 2015 post-season run (low abundance year); i.e., the same response stock catch increase caused by a lower test stock abundance results in a higher exploitation rate when the denominator (abundance of the response stock) is lower. For the few test/response stock combinations where this didn't apply, having higher exploitation rate differences during high abundance years was usually the result of a few key fisheries that were much larger in 2013 than in 2015.

Key Points:
**Percentages in brackets below represent average [exploitation rate](#ER) differences from table 4.6.**<br>
- Nooksack/Samish: This stock aggregate has a considerable impact on Canadian stocks, especially Lower Fraser (9%). The Nooksack/Samish Coho stock is fairly large and has substantial overlap with other stocks in large Coho fisheries. Nooksack/Samish Coho abundance is also important for Skagit (6.6%), Hood Canal (4.1%), and Strait of Juan de Fuca (JDF) (2.5%) stock [exploitation rates](#ER).
- Skagit, Stillaguamish, Snohomish: These stocks have a relatively small effect on Canadian stocks (1.1%-1.8%) but can have a significant effect on Hood Canal Coho (4.8%, 7.7%) and each other; i.e., Stilly/Sno on Skagit (6%). 
- Hood Canal and South Puget Sound: These two stock aggregates are important contributors to several large fisheries. Their abundances can have a significant effect on Skagit, Stillaguamish/Snohomish, U.S. Strait of Juan de Fuca as well as Lower and Upper Fraser stocks. These two stocks also substantially affect each other. Although South Puget Sound Coho are not tracked as part of US/Canada agreements, due to its size this stock aggregate plays an important role in [exploitation rate](#ER) outcomes of wild stocks monitored under the Treaty. 
- U.S Strait of Juan de Fuca: This stock aggregate has a small effect on other U.S and Canadian stocks, mainly due to it low abundance.
- U.S. Coastal Stocks (Quillcene, Hoh, Queets, Quinault, Grays Harbor): Due to their distribution, these stocks have a small effect on Puget Sound and Canadian stocks but can have a large effect on each other.
- Georgia Strait (Mainland and Vancouver Island): Due to their relatively small abundances, these stocks have a fairly small effect on Puget Sound and Coastal stocks yet can have a pronounced effect on Lower Fraser stocks (3.1% Mainland and 2.3% Vancouver Island).
- Southwest Vancouver Island: This stock's abundances can have a significant effect (1.3%) on the Queets wild stock. Queets is a frequent limiter of U.S. ocean and coastal fisheries. This effect is due in large part to overlapping distribution in the West Coast Vancouver Island sport fishery. This stock aggregate also has a 1.2% [exploitation rate](#ER) effect on Skagit and a 1.4% effect on Hood Canal stocks.
- Lower Fraser: Of the Canadian stocks tested, the Lower Fraser stock aggregate has the greatest effect on U.S. stocks, in particular Skagit, Hood Canal, and Strait of Juan de Fuca stocks. The abundance of the Lower Fraser stock aggregate also significantly impacts Upper Fraser (4.2%) and Georgia Strait stocks (7.6%).
- Upper Fraser: This stock aggregate has a relatively small effect on Puget Sound and Coastal stocks (0.2%-0.8%) as well as Canadian stocks (1.4%-1.5%).
- The following 7 fisheries are listed over 70% in the top 5 categories (Table 4.7):
  - Area 5 Sport
  - WCVI Sport
  - Area 6 Sport
  - Area 4/4B Treaty Troll
  - Area 9 Sport
  - 6/7/7A Net Tr
  - BC JDF Sport

Except for Lower Fraser, the Canadian Coho stocks included in the [Pacific Salmon Treaty](#PST) have a relatively minor effect on the [exploitation rates](#ER) of U.S. Puget Sound and Coastal stocks.

Several of the Puget Sound stocks have a significant effect on Canadian stocks, as well as other Puget Sound stocks, in particular Nooksack/Samish, Hood Canal, and South Puget Sound stock aggregates.

With the effect on response stocks in mind, investing in improved forecast assessments to increase the accuracy of [FRAM](#FRAM) [exploitation rate](#ER) estimates would be most beneficial for Nooksack/Samish, Hood Canal, South Puget Sound, and Lower Fraser Coho stocks. 
\newpage

# Comparison of Fishery Adipose [Mark Rates](#MR) from [FRAM](#FRAM) [Validation Runs](#VR) to Mark Rates from Sampling

## Introduction

In the nineties, the Washington Department of Fish and Wildlife (WDFW)
began to remove the [adipose fin](#M) of most hatchery Coho and Chinook
salmon, a practice referred to as mass-marking. A missing adipose fin
became a visual indicator for a hatchery salmon, enabling anglers to
distinguish a wild from a hatchery fish. Subsequently, the Department
introduced mark-selective fisheries ([MSF](#MSF)) as a fisheries management
tool. In MSFs, salmon with a missing adipose fin (marked) can be
retained, while salmon with an intact fin (unmarked) must be released,
thus reducing the impact on wild salmon.

A fishery [mark rate](#MR) is calculated as the number of marked salmon divided
by all the salmon (marked + unmarked) in a fishery.

[Mark rates](#MR) are an indirect measure of the Fisheries Regulation and
Assessment Model's ([FRAM](#FRAM)) ability to accurately predict fishery stock
composition and stock abundances. If modeled mark rates deviate
significantly from observed mark rates, assessments of mortalities on
hatchery and wild stock components will be incorrect.

A 2011 study conducted by Robert Cope [(Cope 2011)](#Cope) found that pre-season
[FRAM](#FRAM) consistently overestimates [mark rates](#MR) in Washington ocean fisheries, resulting
in too few unmarked/wild mortalities. The study speculates on
forecasting error as the cause.

This project uses [mark rates](#MR) from post-season model runs. As these runs are populated with best estimates of actual returns rather than pre-season forecasts, the abundance error is significantly reduced. Mark rates from post-season [FRAM](#FRAM) runs, using reconstructed Coho abundances, are either compared to mark rates reported to the Regional Mark Information System ([RMIS](#RMIS)) for non-selective fisheries or to mark rates of encounters (retained plus released Coho) from sampling programs, for mark-selective fisheries.

## Methods

For this study, [mark rates](#MR) from 2010- 2018 [FRAM](#FRAM) were compared to mark
rates from fishery sampling programs. Most Coho fisheries along the
American West Coast are sampled by fishery technicians (samplers). For
many sport fisheries, samplers inspect the landed (retained) catch for
the presence of an [adipose mark](#M). They also interview anglers about the
mark status of released fish, allowing the calculation of a [mark rate](#MR) of
Coho encounters (landed plus released). This is important, as the mark
rate of the landed catch in mark-selective fisheries is almost 100% by
design. For most commercial fisheries, the sampling programs only report
the mark status of the landed catch, because the mark-status of released
Coho is rarely enumerated on commercial vessels. Therefore, this study
only evaluates the landed [mark rate](#MR) in non-selective commercial
fisheries.

Fishery sampling programs along the entire West Coast report [adipose
mark](#M) information to the Regional Mark Information System [(RMIS)](#RMIS). [RMIS](#RMIS)
houses fishery catch, sampling, and coded-wire-tag [(CWT)](#CWT) information.
[RMIS](#RMIS)' "CatchSample" data provides a one-stop-shop for coast-wide mark
rate information. For this analysis adipose mark rate information was
downloaded from
[RMIS](https://www.rmis.org/cgi-bin/queryfrm.mpl?Table=catch_sample&Version=4.1).

[RMIS](#RMIS) [adipose mark](#M) rates are only reported for landed catches. [Mark rates](#MR)
for [mark-selective fisheries](#MSF) were provided by data stewards from
Washington State's ocean and Puget Sound sampling programs. For
mark-selective sport fisheries, these programs collect additional
information on the number and mark-status of released salmon that allows
the computation of the [mark rate](#MR) of encounters.

[FRAM](#FRAM) [mark rates](#MR) are model projections using
historical stock-/fishery-specific [exploitation rates](#ER) in combination
with post-season estimates of fishery mortalities and stock-specific
abundances (see [FRAM
documentation](https://framverse.github.io/fram_doc/)).

[FRAM](#FRAM) fishery [adipose mark](#M) rates were either compared to [mark rates](#MR) that
came directly from sampling programs or mark rates from [RMIS](#RMIS).
Additionally, a comparison between adipose mark rates from sampling
programs, the primary data source, with mark rates from [RMIS](#RMIS), the
secondary data source, was also conducted to assess whether the data
manipulations performed to import sampling data into [RMIS](#RMIS) can result in
differences[.]{#MRSource}

```{r t61, include = TRUE, results='asis'}
tibble(
  `Source 1` = c("FRAM", "FRAM", "RMIS"),
  `Source 2` = c("RMIS", "Sampling", "Sampling"),
  Fisheries = c("Non-selective sport and commercial fisheries",
                "Mark-selective sport fisheries",
                "Non-selective sport and commercial fisheries")
) |> 
  gt::gt() |> 
  gt::cols_label(.list = list(`Source 1` = gt::md("**Source 1**"), 
                              `Source 2` = gt::md("**Source 2**"),
                              Fisheries = gt::md("**Fisheries**"))) |> 
  gt::tab_header(title = "Table 5.1. Sources for Mark Rate Comparisons by Fishery Type")
```

### General Data Manipulations

For compatibility purposes, all three data sources underwent the
following data manipulations:

**Select Run Years 2010 to 2018** <br>
The data analysis was limited to run years 2010 to 2018. Post-season [FRAM](#FRAM) data after 2009 has received additional vetting. The limitation to 2018 was due to incomplete data reporting post 2018 at the beginning of this project.

**Eliminate sampling estimates with fewer than 20 observations per [time step](#TS)** <br>
[RMIS](#RMIS) and sampling program estimates were removed from the analysis if fewer than 20 Coho were sampled per [time step](#TS). 

**Sum Over [Time Steps](#TS)** <br>
Landed catch and encounter data was summed over [time steps](#TS) within
a year. Coho [FRAM](#FRAM) is stratified into 5 time steps ranging from 1 month
to 6 months in duration (see Appendix [table 8.11](#TStep)). Almost all the analyses are
conducted on annual data to avoid slight differences in [time step](#TS)
definitions between [RMIS](#RMIS) and [FRAM](#FRAM) and to increase the power of the
evaluations.

**Combine Treaty/Non-treaty Net Fishery Pairs** <br>
In [FRAM](#FRAM), most net fisheries exists in treaty, non-treaty pairs; i.e., Area 10 treaty and Area 10 non-treaty net. Differentiating a treaty from a non-treaty net fishery in [RMIS](#RMIS) is often not possible. Additionally, since these
fisheries occur in the same geographic area with similar gear, they
should experience similar [mark rates](#MR). Therefore, catch information from
treaty/non-treaty fishing pairs was summed for [FRAM](#FRAM) and sampling
analyses. See Appendix [table 8.10](#Pairs) for a list of fishery pairs.

**Wilcoxon signed-rank test** <br>
The Wilcoxon signed-rank test is a
non-parametric hypothesis test. Similar to the Student's t-test, it is a
paired difference test. It can be a good alternative to the t-tests when
results are not normally distributed. This test was used to compare mark
rates from [FRAM](#FRAM) and sampling. P-values smaller than 0.05 were selected to be statistically significant, i.e., FRAM [mark rates](#MR) are significantly different from sampling mark rates. The test was applied to annual mark rates for fisheries with at least 5 years of data. It was run as a "paired" test.

### [Mark Rates](#MR) from Sampling

Mark rates from Washington ocean areas 1-4 and Puget Sound fisheries
were provided by the data stewards of the Ocean and Puget Sound Sampling
Programs in January of 2022. The number of marked landed, unmarked
landed, marked releases, and unmarked released were pre-summarized by
fishery, year, and month.

Calculation of [mark rates](#MR) from landed catch: 
Equation 5.1: 
$$MarkRate_{Landed} = \frac{Marked_{Landed}}{Marked_{Landed}+Unmarked_{Landed}}$$

Calculation of [mark rates](#MR) from encounters: 
Equation 5.2: 
$$MarkRate_{Encounters} = \frac{Marked_{Landed}+Marked_{Released}}{Marked_{Landed}+Unmarked_{Landed}+Marked_{Released}+Unmarked_{Released}}$$ 
For each fishery and year combination numbers of Marked_Landed,
Unmarked_Landed, Marked_Released and Unmarked_Released were first summed
over months withing a [time step](#TS) and subsequently summed over [time steps](#TS)
within a year. Release information was only available for
sport fisheries. 

Mark-selective ocean troll fisheries were removed from the sampling
data, because troll fisheries are only sampled for landed catch,
resulting in non-informative mark-rates that were at or near 100% [adipose marked](#M).

Processing Steps:

* Read in data

* Eliminate records that can't be assigned to a [FRAM](#FRAM) fishery

* Combine Treaty/Non-treaty fishing pairs

* Group and summarize data by run year, fishery ID, and [time step](#TS)

* Calculate [mark rates](#MR) by year, fishery, and [time step](#TS)

  + Calculate encounter mark rates for mark-selective fisheries   

  + Calculate landed mark-rates for non-selective fisheries

* Exclude

  + Records with less than 20 sampled per year, fishery, and [time step](#TS)

  + Records for fisheries that did not have a mark selective or non-selective fishery flag (8, 2) in [FRAM](#FRAM)

* Summarize data by run year and fishery ID by grouping over [time steps](#TS)


### [Mark Rates](#MR) from [RMIS](#RMIS)

On May 17, 2022, [RMIS](#RMIS) "CatchSample" data was queried for Coho data
between 2010 and 2020.

These data were mapped to [FRAM](#FRAM) fisheries based on information in the
`catch_location_code`, `fishery`, and `gear` fields. See [fishery mapping description](#fishmap).

[RMIS](#RMIS) fields `period_type` and `period` were utilized to map data to [FRAM](#FRAM)
[time steps](#TS) (see Appendix table [8.11](#TStep)).

During quality control procedures records where `number_caught` or
`number_sampled` was either "na" or zero were removed. Additionally,
records with missing "mark_rates" were eliminated. Designated "M" mixed
mark-selective/non-selective fisheries also resulted in discarding the
record. Records were also removed if the fishery was listed as non-selective, but the [time step](#TS) [mark rate](#MR) was greater than 98%, strongly indicating a mark-selective fishery. 

For each record the number [marked](#M) were calculated as `number_caught`
times [`mark_rate`](#MR).

For each fishery the number caught and the number marked were first
summed over months within a [time step](#TS) and then over [time steps](#TS) within a
year.

The [mark rate](#MR) was calculated as: <br>
Equation 5.3: 
$$MarkRate_{Landed} = \frac{Marked_{Landed}}{Marked_{Landed}+Unmarked_{Landed}}$$

### [Mark Rates](#MR) from [FRAM](#FRAM)

The [Coho Technical Committee](#CoTC) annually produces post-season Coho model
runs. These runs start with the final pre-season model runs and are
subsequently updated with best estimates of stock-specific escapements
and post-season assessments of landed and non-landed fisheries
mortalities. The database is housed on the Pacific Salmon
Commission's (PSC) [Extranet
website](https://extranet.psc.org/SitePages/Home.aspx). Model runs are
updated through 2020 and contained in an Access database called
"PSC_CoTC_PostSeason_CohoFRAMDB_thru2020_021622.mdb".

2010 to 2018 data was extracted from the "Mortality" table using
[framr](https://github.com/FRAMverse/framr), an r package that processes [FRAM](#FRAM) output data.

`LandedCatch`, `MSFLandedCatch`, `MSFNonRetention` data were summarized
by year, fishery, [time step](#TS), and [mark](#M) status. In [FRAM](#FRAM), the `LandedCatch`
column contains the landed catch in a non-selective fishery, the
`MSFLandedCatch` column contains the landed catch in a mark-selective
fishery, and the `MSFNonRetention` column contains the mortality of Coho
released due to mark-selective fishing regulations. `MSFNonRetention`
mortalities were divided by the release mortality rate to convert
mortalities into released encounters.

Equation 5.4: 
$$Encounters_{MSF.NonRetention} = \frac{Mortality_{MSF.NonRetention}}{ReleaseMortalityRate}$$

Mixed mark-selective/non-selective fisheries were excluded from the
data.

For all fisheries, regardless of whether they were non-selective or
mark-selective, the `MSFLandedCatch` in combination with `MSFNonRetention`
were used to compute encounter [mark rates](#MR). In Coho [FRAM](#FRAM) the landed catch
[mark rate](#MR) is identical to the encounter [mark rate](#MR).

Equation 5.5: 
$$Encounters_{MSF.Fishery} = MSFLandedCatch + Encounters_{MSF.NonRetention}$$

[Mark rates](#MR) of encounters were computed as in equation 5.2.

### Combining [FRAM](#FRAM) [Mark Rates](#MR) with Sampling Mark Rates

*As RMIS mark rates are derived from data provided by the sampling programs, this document often refers to RMIS as well as sampling program data as "sampling", unless the two data sources are explicitly considered separately.*

[Mark rates](#MR) provided by Washington sampling programs were filtered to mark selective sport fisheries and then combined via a [full join](#fulljoin) to [FRAM](#FRAM) mark selective fisheries and [time steps](#TS).

Combining [FRAM](#FRAM) [mark rates](#MR) with those provided by [RMIS](#RMIS) required additional steps for annual comparisons. [RMIS](#RMIS) only reports landed catch information. Therefore, only non-selective fisheries were compared using a [full join](#fulljoin) of both data sources. In [FRAM](#FRAM) [terminal fisheries](#Term) are often modeled in the final [time step](#TS) (step 5) even if they occur in time 4. Additionally, [RMIS](#RMIS) fisheries are often reported by management instead of calendar months. This can result in the first or last few days of a month being reported in a different management month and time step. The designation of a fishery as mark-selective or non-selective came from FRAM, because this field is frequently incorrect in [RMIS](#RMIS). If a fishery is only modeled in time 5 then [FRAM](#FRAM) will not report the status ([MSF](#MSF) or non-selective) of a fishery in non-modeled [time steps](#TS). Therefore, when conducting a [full join](#fulljoin) of [RMIS](#RMIS) and [FRAM](#FRAM) data based on [FRAM](#FRAM) fisheries designations, significant amounts of valuable [RMIS](#RMIS) data may be excluded. To prevent omitting RMIS data from non-selective fisheries reported in a different this step than FRAM, [RMIS](#RMIS) fisheries in [time steps](#TS) meeting the following conditions were included in the annual [mark rate](#MR) estimates:<br>
* Mark rate < 0.96% & <br>
* Missing [FRAM](#FRAM) fishery designation (FisheryFlag = "NA") &<br>
* Adjoining [time step](#TS) (time step +/- 1) is non-selective (FisheryFlag = 2).


## Results

### Compare Sampling to [FRAM](#FRAM)

For the first comparison, annual [mark rates](#MR) from sampling or [RMIS](#RMIS)
depending on whether the fishery was mark selective sport or
non-selective, respectively, were compared to [FRAM](#FRAM) [mark rates](#MR). [Mark rates](#MR)
from [RMIS](#RMIS) or sampling programs will be henceforth referred to as either sampled
or observed [mark rates](#MR). Incorporating salmon encounter information from
mark-selective sampling programs added valuable information about mark-selective sport fisheries to the analysis.

Tabular results of [mark rates](#MR) for all available fisheries and [time steps](#TS)
can be found in Appendix table [8.12](#markrates).

Initial [adipose mark](#M) rate comparisons were grouped into 7 geographical
regions. Additionally, Puget Sound marine was separated into commercial
and sport fisheries.

Figures 5.1 and 5.2 were created by summing total landed catch or
encounters as well as total marked catch or encounters over all time
steps withing a year and fishery before computing [mark rates](#MR). Box plot
5.1 displays the difference between the sampling [mark rate](#MR) and the [FRAM](#FRAM)
[mark rate](#MR). The vertical line in the box represents the median
difference, the boundaries of the box represent the first quartile
(0.25) and the "whiskers" (outer lines) the third quartile (0.75).
Outliers are plotted as individual points beyond the whiskers.

The Wilcoxon test was used to compare [mark rates](#MR) from [FRAM](#FRAM) and sampling ([RMIS](#RMIS), sampling programs) for fisheries with at least 5 years of data. P-values smaller than 0.05 were selected to be statistically significant, i.e., FRAM [mark rates](#MR) are significantly different from sampling [mark rates](#MR). The test was applied to annual [mark rates](#MR). 
Of the 51 fisheries evaluated for [mark rate](#MR) differences between [FRAM](#FRAM) and sampling, 28 fisheries had a p-value of less than 0.05. Figures 5.3 and Table 5.2 display results from the Wilcoxon test. 

<br> **Alaska** <br> Sampling and [FRAM](#FRAM) mark rates are low and within a
few percentage points of each other. [FRAM](#FRAM) [mark rate](#MR) are higher than
sampling [mark rates](#MR) for all fisheries in this region. Differences are
statistically significant for 4 of the 5 Alaskan fisheries.

**Canada** <br> [Mark rates](#MR) in Canadian fisheries can vary widely with
very low [mark rates](#MR) in Northern and Central British Columbia and mark
rates around 50% for the West Coast of Vancouver Island. [FRAM](#FRAM) [mark rates](#MR)
can be higher as well as lower for individual fisheries within this
area. [FRAM](#FRAM) and sampling [mark rates](#MR) are within a similar order of
magnitude, except for Johnstone Strait net. [Mark rate](#MR) differences for
this fishery are large between [FRAM](#FRAM) (32.6%) and sampling (6.8%), but due
to low sample size statistical significance cannot be evaluated. An
evaluation of significance is also confounded by small sample sizes for
four other British Columbia fisheries. Only North British Columbia
troll and net, and North-Central troll had enough data points
to evaluate significance. Differences for these fisheries were
significant for North British Columbia net.

**Washington Coastal** <br> [FRAM](#FRAM) [mark rates](#MR) can be higher as well as
lower for individual fisheries within this area. Six of the eight fisheries had higher [FRAM](#FRAM) [mark rates](#MR). 
Four of the Washington Coast fisheries had
enough data points to test statistical significance. Of these,
Quillayute and Grays Harbor net exhibited significant [mark rate](#MR)
differences.

**Washington Ocean** <br> For all six ocean fisheries [FRAM](#FRAM) [mark rates](#MR) were
higher than observed [mark rates](#MR), with five of these fisheries evaluating as
significantly different, and 1 fishery being excluded from the analysis
of significance. 
[Mark rate](#MR) differences are larger in northern areas (ocean areas 3 and 4)
than in southern areas (ocean areas 1 and 2). Four of the fisheries with
significant mark-rate differences were well-sampled mark-selective sport
fisheries. These findings confirm earlier analysis by Robert Cope [(Cope
2011)](#Cope) and warrant further evaluations. See the [case study](#casestudy) in this document for possible causes of the [mark rate](#MR) discrepancies for the Area 4/4B Treaty troll fishery.

**Puget Sound Commercial** <br> Puget Sound fisheries generally have
high [mark rates](#MR), because many Puget Sound stocks are dominated by
mass-marked hatchery fish. Notable exceptions are Stillaguamish, Skagit,
and Snohomish Coho. [FRAM](#FRAM) [mark rates](#MR) can be higher as well as lower for
individual fisheries within this area. Of the 20 fisheries in this
aggregate, 16 contained a sufficient number of data points to be
evaluated using the Wilcoxon test. Of these, six fisheries had [mark rates](#MR)
that differed significantly between [FRAM](#FRAM) and sampling. Differences were
particularly pronounced for Areas 12/12B net, 10E net, and 6/7/7A net.

**Puget Sound Freshwater** <br> [FRAM](#FRAM) [mark rates](#MR) can be higher as well as
lower for individual fisheries within this area. Of the 13 fisheries in
this aggregate, only five contained a sufficient number of data points to
be evaluated using the Wilcoxon test. Of these, three fisheries, Skokomish
River net, Nisqually River net, and Puyallup River net had [mark rates](#MR) that differed
significantly between [FRAM](#FRAM) and sampling. Differences were particularly
pronounced for Skokomish River net. This fishery also had a large
variance of the [mark rate](#MR) estimates. Mark rate differences in freshwater fisheries are a notable result as these fisheries are typically well understood and predominantly comprised of the local stock, eliminating or reducing errors resulting from stock composition differences. This suggests errors in the accounting of hatchery and natural stock components.

**Puget Sound Sport** <br> [FRAM](#FRAM) [mark rates](#MR) can be higher as well as
lower for individual fisheries within this area. All ten fisheries in
this aggregate had sufficient number of data points to be evaluated
using the Wilcoxon test. Of these, six fisheries had [mark rates](#MR) that
differed significantly between [FRAM](#FRAM) and sampling.

**South of Falcon [(SoF)](#SOF)** <br> South of Falcon encompasses ocean areas
between Cape Falcon (Oregon) and the California/Mexico border. All five SoF
fisheries had higher [FRAM](#FRAM) [mark rates](#MR). Three fisheries in this aggregate had
sufficient number of data points to be evaluated using the Wilcoxon
test. Of these, Newport sport had a [mark rate](#MR) that differed
significantly between [FRAM](#FRAM) and sampling[.]{#MRDiff}
<br>

```{r f91}
#produce box plots of absolute % differences

#find duplicates
# dplyr::group_by(RunYear, FisheryID, FisheryName, TimeStep, type) %>%
# dplyr::summarise(n = dplyr::n(), .groups = "drop") %>%
# dplyr::filter(n > 1L) 
  
box1<-mr_yr_fram_samp_rmis|> 
 left_join(d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))

box1out<-ggplot(box1 |> mutate(Label = paste(FisheryArea, "-", FisheryName)),aes(x = FisheryName, y = AbsPctDiff)) +
  coord_flip() +
  geom_boxplot(color = "orange", fill = "orange", alpha = 0.5, outlier.shape = NA) +
  geom_hline(yintercept = 0) +
  scale_x_discrete("") +
  scale_y_continuous("Sampling mark rate - FRAM mark rate") +
  labs(title =str_wrap("Figure 5.1. Differences between Post-season Coho FRAM Modeled Mark Rates and Sampling Observations, 2010-18",50)) +
  theme(axis.text=element_text(size=8), 
        axis.title= element_text(size=10)) +
  facet_wrap(~FisheryArea, scales = "free_y", nrow = 2)

box1out  

ggsave("box1.png", plot=box1out, height=7, width=10)

```
*Differences were calculated as [RMIS](#RMIS) minus [FRAM](#FRAM) mark rates. Boxes to the left of the vertical black/zero line result when [FRAM](#FRAM) mark rates are greater than [RMIS](#RMIS) mark rates.*

Figure 5.2 displays mean absolute [mark rates](#MR) from [FRAM](#FRAM) (blue bars) and sampling (yellow bars) by fishery with a separate panel for each region[.]{#AvgMR} 

```{r f92}
#Bar plots of average annual catch by Fishery and Ocean Area



bar1<-mr_yr_fram_samp_rmis|>
  select(RunYear,FisheryID,FisheryName,FRAM_MrkRate = MR,Samp_MrkRate = MR_sampl)|>
  left_join(d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))|>
  pivot_longer(names_to = "type", values_to = "val", cols =  c("FRAM_MrkRate","Samp_MrkRate"))|>
  group_by(FisheryID,FisheryName,FisheryArea,type)|>
  summarise(across(val, mean),.groups="drop")

ggplot(bar1, aes(FisheryName, val, color = type, fill = type)) +
  coord_flip() +
  geom_col(position = position_dodge(), show.legend = T) +
  scale_color_manual(values = pal, aesthetics = c("color", "fill")) +
  facet_wrap(~FisheryArea, labeller = label_wrap_gen(multi_line = F), scales = "free_y", nrow = 2) +
  #labs(subtitle = d$FisheryName[1]) +
  theme(legend.position = "top",
        legend.text = element_text(size=8),
        legend.key.size = unit(0.4, 'cm'),
        legend.title=element_blank(),
        axis.text=element_text(size=6),
        axis.text.x = element_text(angle = 75),
        strip.text.x = element_text(size = 10))+
  ggtitle("Figure 5.2. 2010-18 Average Fishery Mark Rates from FRAM and Sampling Observations")+
  xlab("Fishery")+
  ylab("Mark Rate")

#calculate grand total mark rates summed over all the stocks excluding Alaska
mr_yr<-mr_yr_fram_samp_rmis|>
  left_join(d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))

mr_yr_no_Ak<-filter(mr_yr,FisheryArea!="Alaska")|>
  select(RunYear,Marked,Total,Marked_samp,Total_samp)|>
  group_by(RunYear)|>
  summarise(across(c(Marked,Total,Marked_samp,Total_samp),sum),.groups="drop")|>
  mutate(MR_rmis=Marked_samp/Total_samp,
         MR_fram=Marked/Total)

mr_yr_no_Ak_Ca<-filter(mr_yr,FisheryArea!=c("Alaska","Canada"))|>
  select(RunYear,Marked,Total,Marked_samp,Total_samp)|>
  group_by(RunYear)|>
  summarise(across(c(Marked,Total,Marked_samp,Total_samp),sum),.groups="drop")|>
  mutate(MR_rmis=Marked_samp/Total_samp,
         MR_fram=Marked/Total)

```

<br>
Figure 5.3 displays Wilcoxon p-values by fishery sorted by magnitude[.]{#Wilcox1} 
```{r f93_t92, results='asis'}
# test whether mark rate differences are statistically significant for each fishery; use annual mark rates for comparisons
FisheryWilcoxTest<-mr_yr_fram_samp_rmis|>
  select(RunYear, FisheryID, FisheryName, FRAM_MrkRate =MR ,Samp_MrkRate = MR_sampl) |>
  # pivot_longer(names_to = "type", values_to = "val", cols =  c("FRAM_MrkRate","Samp_MrkRate"))|>
  group_by(FisheryID)|>
  filter(FisheryID %in% (mr_yr_fram_samp_rmis|>
                           select(RunYear, FisheryID) |> group_by(FisheryID) %>% summarize(n = n()) |> filter(n >= 5))$FisheryID) %>%
  do(tidy(wilcox.test(.$FRAM_MrkRate, .$Samp_MrkRate, alt = "two.sided", paired = TRUE)))
x<-nrow(FisheryWilcoxTest)
  
# graph fishery in order of wilcox significance
WilcoxGraph<-FisheryWilcoxTest |> 
  left_join(mr_yr_fram_samp_rmis |> select(FisheryID, FisheryName), "FisheryID") |>
  mutate(Significant = if_else(p.value > 0.05, "Non-significant", "Significant")) |>
  distinct()

ggplot(WilcoxGraph, aes(x = reorder(FisheryName, -p.value), y = p.value, fill = Significant, color = Significant)) +
  geom_col() +
  geom_hline(yintercept = 0.05, linetype = "dotted") +
  theme_bw() +
  labs(title ="Figure 5.3. Wilcoxon p-values by Fishery", subtitle= "Fisheries with significant mark rate differences (p<0.05) in turquoise") +
  xlab("Fisheries")+
  ylab("p-value")+
  coord_flip()


```
[Table 5.2]{#Wilcox2} displays the number of fisheries evaluated with Wilcox's test and the number of fisheries with significant mark rate differences (p<0.05) within each regional fishery aggregate. Note that all ocean fishery comparisons resulted in significant FRAM versus RMIS mark rate differences.
```{r results ='asis'}
# summarize Wilcox significance by region
WilcoxTable1<-left_join(WilcoxGraph,
                        d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))|>
  select(FisheryID,FisheryArea,Significant)|>
  group_by(FisheryArea)|>
  summarise(n=n(),.groups= "drop")

WilcoxTable2<-left_join(WilcoxGraph,
                        d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))|>
  select(FisheryID,FisheryArea,Significant)|>
  group_by(FisheryArea,Significant)|>
  summarise(n=n(),.groups= "drop")|>
  group_by(FisheryArea)|>
  mutate(Fisheries = sum(n))|>
  mutate(`#Significant`=n)|>
  ungroup()|>
  filter(Significant =="Significant")|>
  dplyr::rename(`#Fisheries` = Fisheries)

WilcoxTable2<-WilcoxTable2[,c(1,4,5)] |>
  gt::gt() |> 
  gt::tab_header(title = "Table 5.2. Number of Fisheries Evaluated with Wilcoxon and Number of Fisheries with Siginificant Mark Rate Differences in Each Area")
  
WilcoxTable2
```

<br>

Figure 5.3 displays annual [mark rates](#MR) for fisheries that had at least 3 annual values. Generally, the variability is higher for observed [mark rates](#MR) from sampling than modeled [mark rate](#MR) from [FRAM](#FRAM)[.]{#Wilcox1}

```{r f93_annualMR}
# Bar plots of annual mark rates by fishery FRAM/Sampling
 
#exclude_fish=c(18,20,21,44,54,63,74,77,84,98,108,113,117,128,131,138,162,166,172,174,179,181,182,187,194,195,196,197)
# eliminate fisheries with less than 3 years of data
x<-mr_yr_fram_samp_rmis|>
  count(FisheryName)|>
  filter(n>2)
mr_yr_fram_samp_rmis1<-mr_yr_fram_samp_rmis|>
  filter(FisheryName %in% x$FisheryName)|>
  dplyr::rename(FRAM_MrkRate = MR, Samp_MrkRate = MR_sampl)


bar2<-mr_yr_fram_samp_rmis1|> 
 left_join(d_uniquefishery2 %>% select(FisheryID,FisheryArea), by = c("FisheryID"))|>
  pivot_longer(names_to = "type", values_to = "val", cols =  c("FRAM_MrkRate","Samp_MrkRate"))


ggplot(bar2, aes(RunYear, val, color = type, fill = type)) +
  geom_col(position = position_dodge(), show.legend = T) +
  #scale_x_continuous("") +
  #scale_y_continuous("Mark rate", labels = scales::percent) +
  scale_color_manual(values = pal, aesthetics = c("color", "fill")) +
  facet_wrap(~FisheryID+FisheryName, labeller = label_wrap_gen(multi_line = F), nrow = 7) +
  #labs(subtitle = d$FisheryName[1]) +
  theme(legend.position = "top",
        legend.text = element_text(size=8),
        legend.key.size = unit(0.4, 'cm'),
        legend.title=element_blank(),
        axis.text=element_text(size=7),
        axis.text.x = element_text(angle = 90),
        strip.text.x = element_text(size = 7))+
  ggtitle("Figure 5.3. 2010-18 Fishery Mark Rates from FRAM and Sampling by Year")+
  xlab("Year")+
  ylab("Mark Rate")

    
#Overview side-by-side bars
#sum over time steps and years. Group into fishery aggregates from same #region with similar FRAM/Sampling relationships
selectfish = c(17,33,37,40,41,39,43,44,47,50,71,81,83,88,91,92,93,97,118,120,122,129,136,138,142,152,154,171,172,174,178,179,181,194,195,196,197,198)    

mr_fram_samp_aggs<-mr_fram_samp_long|>
  group_by(FisheryID,FisheryName,type)|>
    summarize(val=mean(val))|>
  filter(FisheryID %in% selectfish)|>
  mutate(
     Area = case_when(FisheryID==17~ "SoF",
                      FisheryID==33|FisheryID==37 ~ "A1,2Spt",
                      FisheryID==40|FisheryID==41 ~ "A3,4Spt",
                      FisheryID==39|FisheryID==43|FisheryID==44 ~ "T_OcnTrl",
                      FisheryID==47|FisheryID==50|FisheryID==71 ~ "CoastNet",
                      FisheryID==81 ~ "JDF_Net",
                      FisheryID==88 ~ "7,7A_Net",
                      FisheryID==91|FisheryID==92|FisheryID==93 ~ "JDF&SJI_Spt",
                      FisheryID==97 ~ "B_ham_Net",
                      FisheryID==118|FisheryID==120|FisheryID==122|FisheryID==129 ~ "CntrPS_Net",
                      FisheryID==136|FisheryID==138|FisheryID==142 ~ "SPS_Net",
                      FisheryID==152 ~ "HC_Spt",
                      FisheryID==154 ~ "12,12B_Net",
                      FisheryID==171|FisheryID==172|FisheryID==178|FisheryID==179 ~ "NCBC_Net",
                      FisheryID==174|FisheryID==181 ~ "WCVI_Net",
                      FisheryID==194|FisheryID==195|FisheryID==196|FisheryID==197|FisheryID==198 ~ "Alaska"
                      )
     )

mr_fram_samp_aggs<-mr_fram_samp_aggs|>
  group_by(Area,type)|>
  summarize(val=mean(val)) |>
  filter(!is.na(Area))

# Create vector of ordered fisheries by FRAM value
FRAMorder <- mr_fram_samp_aggs |> 
  ungroup() |>
  filter(type == "MR_fram") |>
  arrange(val) |>
  mutate(Area = factor(Area))

# ggplot(mr_fram_samp_aggs |> ungroup() |> mutate(Area = factor(Area, levels = FRAMorder$Area)), aes(val,Area, color = type, fill = type)) +
#   geom_col(position = position_dodge(), show.legend = T) +
#   scale_color_manual(values = pal, aesthetics = c("color", "fill")) +
#   theme(legend.position = "top",legend.text = element_text(size=8),
#         legend.key.size = unit(0.4, 'cm'),legend.title=element_blank())+
#   ggtitle("Figure 6.4. 2010-18 Aggregate Fishery Mark Rates from FRAM and Sampling")+
#   xlab("Mark Rate")+
#   ylab("Fishery")
```

Table 5.3 places fisheries in three categories. Agreement is listed as 'good' when [mark rates](#MR) have an absolute difference of less than 10%, 'moderate' below 30%, and poor above 30%[.]{#GMP}

```{r t93, results='asis'}
# Calculates the relative difference by fishery and ranks fisheries

Rank_RelDiff_Fish <-bar1|>
pivot_wider(names_from = type, values_from = val)|>
mutate(MR_Diff = FRAM_MrkRate-Samp_MrkRate)|>
filter(FRAM_MrkRate>0.1)|>
mutate(AbsRelDiff = abs((FRAM_MrkRate-Samp_MrkRate)/FRAM_MrkRate) )|>
mutate(Agreement = if_else(AbsRelDiff<0.10,"good",if_else(AbsRelDiff<0.30,"moderate","poor")))
  Rank_RelDiff_Fish<-Rank_RelDiff_Fish[order(Rank_RelDiff_Fish$AbsRelDiff),]
  #ggplot(Rank_RelDiff_Fish,aes(x=RelDiff))+geom_histogram(binwidth=0.1) #+
   # scale_x_continuous(limits = c(0, 1))
  # Format as a table
 Rank_RelDiff_Fish$FRAM_MrkRate<-round(Rank_RelDiff_Fish$FRAM_MrkRate,3)
 Rank_RelDiff_Fish$Samp_MrkRate<-round(Rank_RelDiff_Fish$Samp_MrkRate,3)
 Rank_RelDiff_Fish$AbsRelDiff<-round(Rank_RelDiff_Fish$AbsRelDiff,3)

Rank_RelDiff_Fish<- Rank_RelDiff_Fish[,-6]|>
  gt::gt() |>
  gt::tab_header(title = "Table 5.3. 'Good', 'Moderate', 'Poor' Agreement between FRAM and Sampling Mark Rates")

Rank_RelDiff_Fish
```

<br>

Table 5.4 displays the mean [mark rate](#MR) difference ([FRAM](#FRAM) - Sampling) by fishery area. 
[Mark rates](#MR) were first calculated for each fishery and year, then averaged over years, and finally averaged over fishing areas. When averaging [mark rates](#MR) over broad geographic areas such as Southeast Alaska, Canada, or Puget Sound [mark rate](#MR) differences are quite small (<5%) with the exception of Washington Coastal and Ocean fisheries (11.7% and 14.8% respectively)[.]{#MeanMRDiff}

```{r t94, results='asis'}
# Calculates the mean mark rate difference by fishery aggregate

MR_RMSE<-bar1|>
pivot_wider(names_from = type, values_from = val)|>
mutate(MR_Diff = FRAM_MrkRate-Samp_MrkRate)

#mutate(Abs_Diff = abs(MR_Diff))
#write.csv(MR_RMSE,"MR_RMSE.csv")
MR_RMSE<-MR_RMSE|>
group_by(FisheryArea)|>
summarise(across(MR_Diff,mean),.groups="drop")
#mutate(Root_SQ_Diff = sqrt(SQ_Diff))
#rounds to 2 decimals
#numeric_columns<-sapply(MR_RMSE,mode) =='numeric'
MR_RMSE$MR_Diff<-round(MR_RMSE$MR_Diff,3)

MR_RMSE<-MR_RMSE|>
  gt::gt() |>
  gt::tab_header(title = "Table 5.4. Mean Mark Rate Difference by Fishery Area")

MR_RMSE

```

### Compare Sampling to [RMIS](#RMIS)

Adipose [mark rates](#MR) from sampling programs along the West Coast are
reported to [RMIS](#RMIS). [Mark rates](#MR) in [RMIS](#RMIS) are stratified into several time
intervals, with weekly and monthly strata being the most common,
depending on the sample size and the sample design. [RMIS](#RMIS) [mark rates](#MR) are
computed by dividing the number of [adipose marked](#M) Coho by the number of
total Coho sampled for a given fishery and time stratum. For the mark
rate analysis, we used [RMIS](#RMIS) and sampling data as the [mark rate](#MR) source.
As [RMIS](#RMIS) data is derived from sampling data, verification that these data
match became part of the data validation of this project.

[Figure 5.4 Displays 2010-2018 Adipose [Mark Rates](#MR) by Fishery in
Side-by-Side Bar Graphs.]{#MRbyYr}

```{r f94}
#compare sampling to RMIS

# create annual sums of landed and marked and then calculate mark rate
# eliminate fishery 33 (Area 1 sport) because Oregon data is not included in ocean sampling report
  mr_rmis_samp<-filter(mr_rmis_samp, FisheryID !=33)
mr_rmis_samp_annual<-mr_rmis_samp|>
  group_by(RunYear,FisheryID,FisheryName)|>
    summarize_at(c("Total_rmis","Marked_rmis","Total_samp","Marked_samp"), sum)|>
  mutate(MR_RMIS = Marked_rmis/Total_rmis, MR_samp = Marked_samp/Total_samp)|>
select(RunYear, FisheryID, FisheryName, MR_samp ,MR_RMIS) |>
  pivot_longer(names_to = "Source", values_to = "Mark_Rate", cols =  c("MR_samp","MR_RMIS"))
 
ggplot(mr_rmis_samp_annual, aes(RunYear, Mark_Rate, color = Source, fill = Source)) +
    geom_col(position = position_dodge(), show.legend = T) +
    scale_x_continuous(breaks = seq(2010, 2018, by = 1)) +
    scale_y_continuous("Mark Rate", labels = scales::percent) +
    scale_color_manual(values = pal, aesthetics = c("color", "fill")) +
    facet_wrap(~FisheryID+FisheryName, labeller = label_wrap_gen(multi_line = F), nrow = 5) +
    #labs(subtitle = d$FisheryName[1]) +
    theme(legend.position = "top",legend.text = element_text(size=8),
          legend.key.size = unit(0.4, 'cm'),
          legend.title=element_blank(),axis.text=element_text(size=8),
          strip.text.x = element_text(size = 10))+
    ggtitle("Figure 5.4. 2010-18 Fishery Mark Rates from RMIS and Sampling")+
    xlab("Year")

```

As expected, [RMIS](#RMIS) and sampling [mark rates](#MR) closely match. Upon
examination, the small differences in [mark rates](#MR) were found to be caused
by slight differences in [time step](#TS) definitions; i.e., sampling data was
summarized on a calendar month basis and [RMIS](#RMIS) data on a statistical
month basis. Since [time steps](#TS) with small sample sizes were excluded from
the analysis, this resulted in some [time step](#TS) exceeding the sample size
threshold for one data source and being retained while falling short in
the other data source and being eliminated.

## Case Study: Area 4/4B Treaty Troll {#casestudy}
The Washington North Coast Treat troll fishery in ocean areas 4/4B is a well sampled, non-selective fishery. Average Coho catch exceeds 20,000 (excluding 2015, 2016) in the past decade. 
Due to these characteristics, the fishery was chosen as a test case to investigate [mark rate](#MR) differences between [FRAM](#FRAM) and [RMIS](#RMIS). 

From 2010 to 2020 the [FRAM](#FRAM) [mark rate](#MR) averaged 49% and the [RMIS](#RMIS) [mark rate](#MR) 33%. This study examines whether the lower observed sampling [mark rates](#MR) are a function of underestimating unmarked wild fish, overestimating marked hatchery fish, or substantial differences in the stock composition of marked hatchery Coho.


### Steps for Conducting the Analysis
<br>
1. Summarize the marked [FRAM](#FRAM) abundances by stock and year by summing  catch from fisheries that were part of the [exploitation rate](#ER) analysis ([see chapter 6.1](#ERAnalysis)) and adding these to escapement. These abundances are later used to compute stock specific RMIS catches by multiplying them times RMIS exploitation rates as calculated in chapter [6.1](#ERAnalysis). For a valid comparison, it is important that both quantities are in the same units. 

Equation 5.6:
$$MarkedAbundance_{Stk,Yr} = {\sum_{fish \in Study}Mortalities_{Stk,Yr} + Escapement_{Stk,Yr}}$$
<br>
2. Calculate associated unmarked abundances. 
Unmarked abundances could not be calculated as mortalities plus catch, because of the bias introduced by releasing unmarked Coho in mark-selective fisheries. 
Instead, the stock- and year-specific [mark rate](#MR) was summarized from post-season [FRAM](#FRAM) [validation runs](#VR) using the time 4 'working cohort' from [FRAM](#FRAM)'s "Cohort" table. This [mark rate](#MR) was then applied to the marked abundance from step 1 to calculate the associated unmarked abundance.

Equation 5.7:
$$UnmarkedAbundance_{Stk,Yr} = \frac{MarkedAbundance_{Stk,Yr}}{MarkRate_{Stk,Yr}}-MarkedAbundance_{Stk,Yr}$$
<br>
3. Compute [FRAM](#FRAM) and [RMIS](#RMIS) catches as abundance times [pseudo-exploitation rate](#Pseudo). 
For the purpose of computing [pseudo-exploitation rates](#Pseudo) the same subset of fisheries was utilized in the denominator as for the [exploitation rate](#ER) analysis (see paragraph 6.1).

Equation 5.8:
$$MarkedFRAMCatch_{Stk,Yr}=MarkedAbundance_{Stk,Yr} * ER.pseudo.FRAM_{Stk,Yr}$$
Equation 5.9:
$$UnMarkedFRAMCatch_{Stk,Yr}=UnMarkedAbundance_{Stk,Yr} * ER.pseudo.FRAM_{Stk,Yr}$$
Equation 5.10:
$$MarkedRMISCatch_{Stk,Yr}=MarkedAbundance_{Stk,Yr} * ER.pseudo.RMIS_{Stk,Yr}$$
Equation 5.11:
$$UnMarkedRMISCatch_{Stk,Yr}=UnMarkedAbundance_{Stk,Yr} * ER.pseudo.RMIS_{Stk,Yr}$$
<br>
4. Sum catches over stocks by year for Area 4/4B treaty troll for marked and unmarked [FRAM](#FRAM) and [RMIS](#RMIS) catches.

Equation 5.12:
$$Catch_{Yr}= {\sum_{Stk}Catch_{Stk,Yr} }$$
<br>
5. Summarize [FRAM](#FRAM) catches of marked and unmarked hatchery fish and marked and unmarked wild fish in the treaty troll fishery. The catches can be retrieved from standard [FRAM](#FRAM) output and were summarized by year.
The total annual [FRAM](#FRAM) catch was summed over hatchery as well as wild stocks landed in the fishery. Values summarized from post-season [FRAM](#FRAM) runs matched the values calculated in step 4 above. This provides a validation of the methods uses to calculate [FRAM](#FRAM) [pseudo-exploitation rates](#Pseudo) and reconstruct associated abundances. Furthermore, the [FRAM](#FRAM) database also provided the estimated wild catch and total catch of the fishery. The total annual landed catch modeled in FRAM for the treaty troll fishery also matches the catch reported by RMIS. Thus, FRAMTotCatch can be substituted with RMISTotCatch in below equations.

Equation 5.13:
$$FRAMTotCatch_{Yr}=\sum{(Catch_{HatMrkd}+Catch_{HatUM}+Catch_{WildMrkd}+Catch_{WildUM})}$$
<br>
6. Calculate [FRAM](#FRAM) and [RMIS](#RMIS) [mark rates](#MR) by dividing the marked catch by the total catch. 

Equation 5.14:
$$MR_{FRAM}=\frac{FRAMMrkdCatch_{Yr}}{FRAMTotCatch_{Yr}}$$
Equation 5.15:
$$MR_{RMIS}=\frac{RMISMrkdCatch_{Yr}}{FRAMTotCatch_{Yr}}$$
Where, <br>
  Stk       = Stock<br>
  Yr        = Year<br>
  Mrkd      = Marked<br>
  UM        = Unmarked

### Results

Table 5.5 displays [FRAM](#FRAM) and [RMIS](#RMIS) landed catches by mark status, observed [mark rates](#MR) from sampling programs, and [FRAM](#FRAM) and [RMIS](#RMIS) [mark rates](#MR) as calculated in steps 1-6 above[.]{#MarkStatus} 
```{r results ='asis'}
d_ERtrollsum<-read.csv(file.path(table_dir,"ER_Ttrollsum.csv"))
CaseStudyTable<-d_ERtrollsum|>
   # DT::datatable(filter = 'top', 
   #                   extensions = 'Buttons', 
   #                   options = list(
   #                     dom = 'Bfrtip',
   #                buttons = c('copy', 'csv', 'excel', 'pdf', 'print')))
  
  #DA help. Can't get table to render
  
  gt::gt() |> 
  gt::cols_label(.list = list(`X1.Year` = gt::md("**(1)Year**"),
                              `X2.FRAM.Mrkd.Hatchery.Catch` = gt::md("**(2) FRAM Mrkd Hatchery Catch**"), 
                              `X3.FRAM.UM.Hatchery.Catch` = gt::md("**(3) FRAM UM Hatchery Catch**"),
                              `X4.FRAM.Hatchery.Mark.Rate`=gt::md("**(4) FRAM Hatchery Mark Rate**"),
                              `X5.Reconstructed.RMIS.Mrkd.Catch` = gt::md("**(5) Reconstr. RMIS Mrkd Catch**"), 
                              `X6.RMIS.UM.Hatchery.Catch` = gt::md("**(6) RMIS UM Hatchery Catch**"),
                              `X7.RMIS.Hatchery.Mark.Rate`=gt::md("**(7) RMIS Hatchery Mark Rate**"),
                              `X8.All.Stocks.FRAM.Catch` = gt::md("**(8) All Stocks FRAM Catch**"), 
                              `X9.Observed.Mark.Rate` = gt::md("**(9) Observed Mark Rate**"),
                              `X10.FRAM.Mark.Rate`=gt::md("**(10) FRAM Mark Rate**"),
                              `X11.RMIS.Marked.Catch`=gt::md("**(11) RMIS Marked Catch**")
                              )) |> 
  gt::tab_header(title = "Table 5.5.Comparison of FRAM and RMIS Catches by Mark Status ") 
CaseStudyTable
```
<br>
Although the hatchery stock composition displayed some differences between [FRAM](#FRAM) and [RMIS](#RMIS) (see table 5.6), the [mark rate](#MR) of the hatchery stocks was very similar and averaged 89% for both data sources[.]{#MarkedCatch} 

```{r results ='asis'}
#DA help. can't get table to render
tbltop10<-d_top10|>
  gt::gt()|>
   gt::tab_header(title = "Table 5.6. Comparison of FRAM and RMIS Marked Catches in the Area 4/4B Treaty Troll Fishery ")
tbltop10
```
<BR>
The mean percent error between RMIS "observed" marked catches (Table 5.5 Col 11) and [FRAM](#FRAM) (Table 5.5 Col 2) [marked](#M) catches is 55% and between RMIS "observed" marked catches and [RMIS](#RMIS) (Table 5.5 Col 5) reconstructed catches is -2%. To clarify, the difference between RMIS "observed" and RMIS reconstructed catches is the following: (1) the former is calculated as reported catch times observed (sampled) mark rate, while (2) the later is computed independently as the FRAM abundance of marked hatchery stocks caught in the fishery times RMIS exploitation rates as calculated in chapter [6.1.](#C61). The average [RMIS](#RMIS) marked catch deviates only 100 Coho from the observed catch of 6700. There is also good agreement in "observed" and reconstructed marked [RMIS](#RMIS) catches on an annual basis (R Square=0.89). It is noteworthy that RMIS ERs are likely too high for some stocks, because of incomplete escapement accounting. However, this bias results in higher RMIS marked catches and mark rates. Therefore, if corrected, would exacerbate FRAM/RMIS mark rate disparities.
Even though [RMIS](#RMIS) and [FRAM](#FRAM) stock compositions differ, there is good agreement between hatchery [mark rates](#MR) of both sources. This is perhaps not too surprising as many stocks caught in the treaty troll fishery originate from hatchery programs with mass marking, where almost all released hatchery Coho are [adipose marked](#M). Additionally, [RMIS](#RMIS) [marked](#M) catches as computed using equation 5.10, while much lower than FRAM estimates of marked catches, provide good accounting for observed marked catches (Table 5.5 Col 11). Therefore, if the stock composition of hatchery Coho does not result in signification mark rate differences and RMIS reconstructed marked catches provide good accounting of "observed" catches, unmarked, wild fish are likely the missing component. This points towards missing/underestimated wild abundance as the source of the discrepancy. The average composition of the catch based on sampling and hatchery [mark rates](#MR) is approximately 32% hatchery marked, 2% hatchery unmarked, and 64% wild unmarked. Given these proportions and assumptions, [FRAM](#FRAM) misses approximately 30% of the wild abundance. Since the [mark rate](#MR) discrepancies decrease for fisheries further South, Washington South Coast and Columbia River stocks are less likely to be the cause of the [mark rate](#MR) differences.  
We recommend additional evaluations and if possible Genetic Stock Identification studies to investigate further.

## Conclusions

For many fisheries, agreement between [FRAM](#FRAM) [mark rates](#MR) and sampling [mark rates](#MR) is moderate or poor. Although when aggregating over fisheries and years within large geographic areas, [mark rate](#MR) differences can be quite small (see table 9.4).

For some fisheries [FRAM](#FRAM) [mark rates](#MR) are consistently higher than sampling [mark rates](#MR). This is especially apparent for Washington coastal and ocean fisheries and warrants [further investigation](#casestudy).

In non-selective fisheries, [FRAM](#FRAM) [mark rates](#MR) are largely a function of modeled marked and unmarked stock abundances and [base period](#BP) [exploitation rates](#ER). Since [base period](#BP) exploitation rates are static, [FRAM](#FRAM) [mark rates](#MR) are generally less variable across years than sampling [mark rates](#MR). 

Some of the greatest [mark rate](#MR) differences can be observed in freshwater fisheries, i.e., Skokomish net, Skagit sport. Due to the [terminal](#Term) nature of these fisheries, [mark rates](#MR) are usually well known and driven by the [mark rate](#MR) of the local stock. Therefore, differences may be due to sampling error or fishing in locations where the hatchery (predominantly [adipose marked](#M)) and natural (predominantly unmarked) stock components segregate.    

[Mark rate](#MR) comparisons are just one parameter useful in assessing the accuracy of the [FRAM](#FRAM) model. As mentioned earlier, [FRAM](#FRAM) [mark rates](#MR) are a function of base period [exploitation rates](#ER) and stock specific abundances. It is often difficult to assess whether [mark rate](#MR) differences are caused by a change in stock distributions, forecast errors (see [case study](#casestudy)), or other data errors. Both data sources, [FRAM](#FRAM) and [RMIS](#RMIS) have many potential sources of error. We found mistakes in [FRAM](#FRAM) abundance data, an annual input, which resulted in significant [mark rate](#MR) errors for some stocks and years. Conversely, the sampling data had many instances of problematic sampling sizes or highly improbable [mark rate](#MR) observations. Well sampled fisheries are crucial to mark rate and exploitation rate assessments.  

\newpage

# Comparison of Post-season [FRAM](#FRAM) Fishery Impacts with [RMIS](#RMIS) Recovery Patterns {#ERAnalysis}

This analysis examines post-season ("validation") estimates of stock and fishery specific mortality generated by [Coho FRAM](https://framverse.github.io/fram_doc/calcs_data_coho.html), comparing them to patterns of coded-wire-tag [(CWT)](#CWT) recoveries housed in [RMIS](https://www.rmpc.org/).

Coho [FRAM](#FRAM) predicts impacts to modeled stocks according to a "[base period](#BP)" parameterization from [CWT](#CWT) recoveries in catch years 1986-1992 (see https://framverse.github.io/fram_doc/index.html for additional description). The continued applicability of this base period is of interest to fisheries managers concerned with meeting conservation and harvest allocation objectives.

Comparing a dataset of recent recovery patterns with results based on older CWT recoveries does not address the fundamental assumption that inferences about the harvest of non-hatchery stocks are adequately informed by [coded-wire tagging](#CWT) that is overwhelmingly focused on hatchery production. Nor does this analysis attempt to determine whether recent year Coho [FRAM](#FRAM) post-season exploitation rates [(ER)](#ER) are quantitatively accurate relative to those that might be reconstructed outside of [FRAM](#FRAM): arriving at a "true" stock ER requires data from robust sampling throughout all phases of harvest and post-harvest escapement, and these data are not readily available for most of the Coho [FRAM](#FRAM) stock units.

Instead, this investigation focuses on (dis)agreements in the relative assignment of catch mortality in primarily [pre-terminal](#PTerm), highly mixed stock fisheries, seeking to understand how well the [FRAM](#FRAM) representation agrees with current observations in [RMIS](#RMIS). That is, we ask, "For a given stock, how well do the magnitude and relative ranks of annual per-fishery exploitation align between datasets?" Attempting to answer this seemingly straightforward question also highlights discrepancies that warrant additional auditing to ensure that data are sufficiently high quality to support management decisions and risk considerations.

## [Methods]{#C61}

Coho [FRAM](#FRAM) includes net, troll and sport fisheries that occur in [pre-terminal](#PTerm), [terminal](#Term) and freshwater areas. Freshwater sport, most freshwater net and some near-terminal net fisheries were excluded from the analysis due to the practical constraint of limited and inconsistent sampling and the conceptual presumption that the stock composition of these fisheries is subject to much less uncertainty, notwithstanding some amount of straying. In addition, [RMIS](#RMIS) Catch & Sample data were inspected and compared with fishery catch values used in post-season [FRAM](#FRAM) runs. Fisheries were excluded if the proportion of catch that was sampled was chronically low. Where smaller annual differences occurred more sporadically, an expansion coefficient was generated. These coefficients were then applied to the sampling-expansion based "estimated_number" values reported for individual recoveries, thereby better aligning the fishery totals between [FRAM](#FRAM) and [RMIS](#RMIS). See section [Handling of Fisheries in Response to RMIS Catch/Sample data](#Handling) in this report for detail descriptions of how fisheries were selected and/or adjusted for this project. 

```{r included_fisheries, eval=FALSE}
# #excluding:
# - freshwater sport (little/no sampling)
# - most terminal (stock comp presumed less varied; very uneven sampling)
# - others with little/no recent year catch and/or sampling
# included:
fishery_ids <- c(
  15, 17, 18, 19, 20, 21, 22, #CA/OR
  23:33, 34, 35, 37, 38, 39, 40, 41, 42, 43, #B10/WA Ocn
  44, 45, 47, 48, 50, 81, 83, 88, 91, 92, 93, 
  97, 98, 102, 106, 107, 110, 112, 115,
  118, 120, 122, 124, 125, 129, 136, 142, 144,
  152, 154, 156, 158, 160, #161, 162, #161 Skok R 162 Quilcene R
  171, 172, 174, 175, 178, 179, 182, 189, #BC
  194, 195, 196, 197, 198, #AK
  #NT members of Tr/NT pairs missing from fram_samp_fulla
  80, 82, 87, 96, 101, 109, 111, 119, 121, 123, 
  130, 132, 137, 139, 141, 143, 145, 153, 155, 157, 159
  )

#per AHB request to combine NT & Tr net fisheries
#dropping several NT of now dropped pairs: 130, 132, 137, 139, 146 
fishery_ids_paired <- c(
  80, 82, 87, 96, 101, 109, 111, 119, 121, 123, 
  #130, 132, 137, 139,
  141, 143,
  #145,
  153, 155, 157, 159
  )

#per-fishery-year coefficients based on the proportion sampled of RMIS-reported catch
#these are applied to adjust the "estimated_number" values of the stock-specific recoveries in that fishery-year
#if a fishery-year reporting a catch of 1000 had sampling for 500, the coefficient value is 2
#attempting to "expand" the sampled portion of the total catch in RMIS
#the ideal situation would have the FRAM landed catch input match the RMIS total and sampled values
samp_catch_coef <- read_csv(fp$samp_catch_coef) |> select(FisheryID:Action)

```

### [FRAM](#FRAM)

Estimates of landed catch and post-fishing escapement (from the `Mortality` and `Escapement` database tables respectively) were queried from the Pacific Salmon Commission [Coho Technical Committee](#CoTC) post-season [FRAM](#FRAM) database housed at PSC's [Extranet website](https://extranet.psc.org/SitePages/Home.aspx). Model runs are currently available through 2020 and stored in Access database "PSC_CoTC_PostSeason_CohoFRAMDB_thru2020_021622.mdb". **The 2015 "crash" and 2016 "post-crash" years were included in the dataset for some analyses but were generally considered uninformative of "typical" patterns.**

```{r get_starting_fram_dataset, eval=FALSE}
fram <- list()

# ## initial focus on relatively data-rich units, particularly with PST relevance
# #Skagit, Sno, Quilcene, George Adams, Nisqually, Green, Quillayute, Queets, Quinault, Chehalis, Willapa, Columbia Early, Columbia Late, Fraser Lower, Fraser Upper
# stock_ids <- c(20, 38, 48, 58, 68, 96, 134, 142, 148, 152, 164, 166, 176, 226, 230)

# #all hatchery stocks, M and UM, including Northern hat/wild
# tbl(db_con, "Stock") |> filter(Species=="COHO") |> collect() |> filter(str_detect(StockLongName, "atch|ens")) |> pluck("StockID") |> paste(collapse = ", ")
stock_ids <- c(3, 4, 5, 6, 7, 8, 9, 10, 15, 16, 19, 20, 21, 22, 25, 26, 27, 28, 
  31, 32, 33, 34, 37, 38, 39, 40, 41, 42, 47, 48, 49, 50, 53, 54, 
  57, 58, 65, 66, 67, 68, 71, 72, 73, 74, 77, 78, 79, 80, 83, 84, 
  87, 88, 91, 92, 95, 96, 99, 100, 103, 104, 109, 110, 113, 114, 
  119, 120, 125, 126, 129, 130, 133, 134, 137, 138, 141, 142, 143, 144,
  147, 148, 151, 152, 155, 156, 159, 160, 163, 164, 165, 166, 167, 168,
  175, 176, 177, 178, 181, 182, 185, 186, 189, 190, 191, 192, 193, 194, 
  197, 198, 201, 202, 205, 206, 209, 210, 213, 214, 217, 218, 221, 222, 
  225, 226, 229, 230, 233, 234, 235, 236, 237, 238)

run_ids <- 34:44 #2010:2020

db_con <- DBI::dbConnect(
  drv = odbc::odbc(),
  .connection_string = paste0(
    "Driver={Microsoft Access Driver (*.mdb, *.accdb)};DBQ=",
    fp$mdb_pst, ";")
  )

fram$fishery <- tbl(db_con, "Fishery") |> 
  filter(Species == "COHO") |> 
  select(FisheryID:FisheryTitle) |> 
  collect()

DBI::dbDisconnect(db_con)


fram$escp <- framr::read_coho_escp(fp$mdb_pst, runs = run_ids, stocks = stock_ids)

fram$mort <- framr::read_coho_mort(fp$mdb_pst, runs = run_ids, stocks = stock_ids)

# fram$mort |> distinct(StockID, StockLongName) |> arrange(StockID) |> gt()
# fram$fishery |> filter(FisheryID %in% fishery_ids) |> gt()

```

### [RMIS](#RMIS)

All recoveries of [marked](#M) Coho between the years 1986 and 2020 were queried from [RMIS](#RMIS), and these were joined with the release data associated with the distinct tag codes in this set. This core dataset was augmented with recent year Canadian stock escapement recovery records provided by DFO. An age field was calculated as the difference between `run_year` and `brood_year`, allowing later limitation to the 3-year-old fish tracked in [FRAM](#FRAM), and records with a `sample_type` 5 were excluded as uninformative. Finally, the individual records were mapped to [FRAM](#FRAM) stocks according to a `tag_code` look-up table (see section [Assigning Coded-Wire-Tagged Groups to FRAM Model Stocks](#SM) for more information) and to [FRAM](#FRAM) fisheries via reference to the `recovery_location_code` and `fishery` fields (see section [Assigning Coded-Wire-Tagged Recoveries to FRAM Fisheries](#FM)). Appendix table [8.7](#TagRec) contains a list of [CWT](#CWT) recoveries in fisheries and escapement from 2010-2019 by [FRAM](#FRAM) stock . <br>


We then examined the recovery patterns for each stock to select a group of focal stocks for this analysis. These stocks were selected based on their importance to the [Pacific Salmon Treaty](#PST) and the number of [CWT](#CWT) recoveries in mixed stock fisheries and escapement.

<!-- The first three chunks preserve one-time operations to generate starting datasets from [RMIS](#RMIS), and the fourth performs the stock and fishery mapping. -->

```{r read_in_rmis_recoveries, eval=FALSE}
dir_rmis_recs <- "T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/ERAnalysis/RMISRecoveries/"

focal_cols_recs <- c(
  "recovery_id", "tag_code", "run_year", 
  "recovery_location_code", "recovery_location_name", 
  "fishery", "gear", "estimated_number", "recorded_mark",
  "catch_sample_id", "sample_type")

rmis_recs <- purrr::map_dfr(
  list.files(dir_rmis_recs, pattern = ".TXT", full.names = T)
  ,
  ~readr::read_csv(.x, col_select = all_of(focal_cols_recs)) |> 
    dplyr::mutate(recorded_mark = as.character(recorded_mark))
  )

#saveRDS(rmis_recs, "T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/rmis_coho_recoveries_86_2020.rds")

```

```{r get_releases_by_tag_code, eval=FALSE}
rmis_recs|>
 distinct(recovery_location_code)
rmis_recs |>
  distinct(tag_code) |> rownames_to_column(var = "id") |>
  mutate(
    id = as.numeric(id),
    block = cut(id, breaks = c(seq(0,20000,by=2000))),
    block = letters[as.integer(block)]
    ) |> #count(block)
  split(~block) |>
  writexl::write_xlsx("T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/ERAnalysis/tag_codes_by_block.xlsx")

f <- "T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/ERAnalysis/tag_codes_by_block.xlsx"

readxl::read_excel(f, sheet = "a") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "b") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "c") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "d") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "e") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "f") |> pluck("tag_code") |> writeLines("clipboard")
readxl::read_excel(f, sheet = "g") |> pluck("tag_code") |> writeLines("clipboard")

focal_cols_rel <- c(
  "tag_code_or_release_id", "brood_year",
  "release_location_code", "hatchery_location_code", "stock_location_code",
  "release_location_name", "hatchery_location_name", "stock_location_name",
  "release_location_state", "release_location_rmis_region", "release_location_rmis_basin",
  "cwt_1st_mark", "cwt_1st_mark_count", "cwt_2nd_mark", "cwt_2nd_mark_count"
  )

rmis_rel <- map_df(
  c(
  "https://www.rmis.org/reports/CSV13497.txt",
  "https://www.rmis.org/reports/CSV13656.txt",
  "https://www.rmis.org/reports/CSV13738.txt",
  "https://www.rmis.org/reports/CSV13814.txt",
  "https://www.rmis.org/reports/CSV13854.txt",
  "https://www.rmis.org/reports/CSV14030.txt",
  "https://www.rmis.org/reports/CSV14116.txt"
  ),
  ~read_csv(.x, col_select = all_of(focal_cols_rel)) |> 
    mutate(
      cwt_1st_mark = as.character(cwt_1st_mark),
      cwt_2nd_mark = as.character(cwt_2nd_mark)
      )
  )

#saveRDS(rmis_rel, "T:/DFW-Salmon Mgmt Modeling Team - General/Southern Fund Projects/Coho FRAM MSM CWT analysis/rmis_coho_releases_for_rec_86_2020.rds")

```

```{r join_rel_to_rec, eval=FALSE}
left_join(
  readRDS(fp$rmis_recs) |> 
    mutate(
      rlsp = str_sub(recovery_location_code, 1, 1),
      rlwtr = str_sub(recovery_location_code, 2, 2),
      rlsctr = str_sub(recovery_location_code, 3, 3),
      rlrgn = str_sub(recovery_location_code, 4, 5),
      rlarea = str_sub(recovery_location_code, 6, 9),
      rllocn = str_sub(recovery_location_code, 10, 16),
      rlsubl = str_sub(recovery_location_code, 17, 19),
      f_type = case_when(
        between(fishery, 10, 19) ~ "troll",
        between(fishery, 20, 39) ~ "net",
        between(fishery, 40, 49) ~ "sport",
        between(fishery, 50, 59) ~ "escp"
      )
    )
  ,
  readRDS(fp$rmis_rel) |> select(-contains("n_code")),
  by = c("tag_code" = "tag_code_or_release_id")
  ) |> #saveRDS(file.path(dir_proj, "rmis_coho_rec_rel_86_2020.rds"))

```

```{r rmis_rr, eval=FALSE}
#primary dataset, constructed and preserved in chunks above
rmis_rr <- readRDS(fp$rmis_rec_rel) |> 
  filter(
    sample_type != 5,
    !is.na(release_location_state)
    ) |>  #drop 3 recent
  mutate(age = run_year - brood_year) |> 
  rowid_to_column("uid")

#adding missing Canadian escapements
rmis_rr <- bind_rows(
  rmis_rr,
  readxl::read_xlsx(
    file.path(dir_proj, "ERAnalysis/2022-02-23 StAD Escapement Report_Coho_EPAD.xlsx"),
    sheet = "Obs-Est-Exp"
    ) |>
    rowid_to_column(var = "uid") |> 
    mutate(
      uid = uid + last(rmis_rr$uid),
      age = `Recovery Year` - `Brood Year`,
      fishery = 50,
      f_type = "escp"
    ) |> 
    select(
      uid,
      tag_code = Tagcode, 
      recovery_location_code = `PSC Location Code`,
      recovery_location_name = `Recovery Site Name`,
      run_year = `Recovery Year`,
      brood_year = `Brood Year`,
      age,
      fishery, f_type,
      estimated_number = `Final Expanded`
    ) |> 
    filter(
      tag_code %in% unique(rmis_rr$tag_code),
      age == 3
    ) 
  )


#add CCT mapped StockIDs
rmis_rr <- rmis_rr |> left_join(lu$cwt_stk, by = c("tag_code"))

#add AHB mapped FisheryIDs
#cannot directly join on 'PSC_Code'...have to use the Num_Chars field
#but that means using many possible lengths of recovery_location_code
#brute force iterator over Num_chars + uid to allow recombination

#first declare temp object of only key fields
rmis_rr_ahb_map <- rmis_rr |> 
  select(uid, recovery_location_code, fishery) |> 
  mutate(fishery_single = floor(fishery/10)) 
#now overwrite as split-then-bind pattern on Num_Chars
#per AHB request, added explicit right string padding to specified string length
#even when actual strings in RMIS and LU are/use fewer char
rmis_rr_ahb_map <- map_df(
  sort(unique(lu$ahb_fishmap$Num_Chars)),
  ~inner_join(
    rmis_rr_ahb_map |> 
      mutate(PSC_Code = str_sub(recovery_location_code, 1, .x) |> str_pad(width = .x, side = "right"))
    ,
    lu$ahb_fishmap |> filter(Num_Chars == .x) |> 
      mutate(PSC_Code = str_pad(PSC_Code, width = .x, side = "right"))
    ,
    by = c("PSC_Code", "fishery_single")
  )
)

#some cases of same fishery from different n-chars being considered on same recovery_location_code
#some cases initially from simple row duplication in look-up, now fixed on read-in
overmapped <- rmis_rr_ahb_map |> count(uid) |> filter(n>1) |> select(uid)

#now rejoin to full dataset, after making FisheryID distinct per UID
rmis_rr <- left_join(
  rmis_rr, 
  #split out multiple-mapped, make distinct per UID and rebind
  bind_rows(
    anti_join(rmis_rr_ahb_map, overmapped, by = "uid") |> 
      select(uid, FisheryID)
    ,
    semi_join(rmis_rr_ahb_map, overmapped, by = "uid") |> 
      group_by(uid) |> 
      slice_max(Num_Chars, n = 1, with_ties = F) |> 
      ungroup() |> 
      select(uid, FisheryID) #|> count(uid) |> filter(n>1)
  )
  ,
  by = "uid") |> 
  left_join(fram$fishery, by = "FisheryID")

rm(rmis_rr_ahb_map, overmapped)

#split Ocean troll into treaty/non-treaty
rmis_rr <- rmis_rr |> 
  mutate(
    FisheryID = if_else(FisheryID %in% c(35, 38, 42) & fishery==15, FisheryID+1, FisheryID),
    FisheryName = if_else(FisheryID==36, "Area2TrlTR", FisheryName),
    FisheryName = if_else(FisheryID==39, "Area3TrlTR", FisheryName),
    FisheryName = if_else(FisheryID==43, "A4/4BTrlTR", FisheryName)
  ) 

# rmis_rr |> filter(FisheryID %in% c(35, 36, 38, 39, 42, 43)) |> distinct(FisheryID, FisheryName, fishery, gear)
```

#### Recovery Patterns Since 1986 {#RP}

Recent years have generally seen fewer Coho [CWT](#CWT) recoveries than during the Coho [FRAM](#FRAM) [base period](#BP) years, within both sampled fisheries and escapement [(to hatcheries and in-river)]{#RecCount}. 

```{r rmis_rr_rgn, eval=FALSE}
rmis_rr_rgn <- rmis_rr |> 
  filter(between(run_year, 1986, 2019), between(as.integer(cwt_1st_mark), 5000, 5500)) |> 
  count(run_year, f_type, release_location_state, release_location_rmis_region, release_location_rmis_basin) |> 
  mutate(
    rgn = case_when(
      release_location_state == "AK" ~ "AK",
      release_location_state == "BC" | (is.na(release_location_state) & f_type == "escp") ~ "BC",
      
      release_location_rmis_region == "CECR" | release_location_rmis_region == "CRGN" |
        release_location_rmis_region == "LOCR" | release_location_rmis_region == "UPCR" ~ "Columbia R.",
      
      release_location_rmis_region == "HOOD" | release_location_rmis_region == "JUAN" |
        release_location_rmis_region == "MPS" | release_location_rmis_region == "NOWA" |
        release_location_rmis_region == "NPS" | release_location_rmis_region == "SKAG" |
        release_location_rmis_region == "SPS" ~ "WA Puget Sound",
      
      release_location_rmis_region == "GRAY" | release_location_rmis_region  == "NWC" |
        release_location_rmis_region =="WAGN" | release_location_rmis_region == "WILP" ~ "WA Coastal",
      
      release_location_rmis_region == "SOOR" | release_location_rmis_region == "KLTR" |
        release_location_rmis_region == "NOCA"  ~ "SONCC",
      
      release_location_rmis_region == "NOOR" | release_location_rmis_region == "ORGN" |
        release_location_rmis_region == "CECA" | release_location_rmis_region == "SAFA" |
        release_location_rmis_region == "SNAK" ~ "OR,CA,SNAK"
    ),
    yrs = case_when(
      between(run_year, 1986, 1992) ~ "1986-1992",
      between(run_year, 2010, 2019) ~ "2010-2019"
    )
  )

```

```{r recs_by_yrs, fig.height=6, fig.width=7}
rmis_rr_rgn |> 
  group_by(yrs, f_type) |> 
  summarise(n = sum(n), .groups = "drop") |> 
  drop_na(yrs, f_type) |> 
  ggplot(aes(yrs, n, fill = f_type)) +
  geom_col(position = position_dodge()) +
  scale_fill_manual(values = pal_ftype, name = "") +
  scale_x_discrete("") +
  scale_y_continuous("CWT count", labels = scales::comma) +
  labs(title = "Figure 6.1. Recovery Counts by Time Period and Type")
```

These decreases can be differentiated relative to the locations of tagged releases.<br>
[**Figure 6.2. Annual CWT Recoveries by Release Location in Fisheries and Escapement**]{#AnnualRecCount}
```{r recs_by_rgn_series_cols, fig.height=11, fig.width=9}
(
rmis_rr_rgn |>  
  filter(f_type != "escp") |> 
  ggplot(aes(run_year, n, fill = f_type)) + 
  geom_col() +
  scale_x_continuous("", n.breaks = 8, guide = guide_axis(n.dodge = 2)) +
  scale_y_continuous("CWT count", labels = scales::comma) +
  scale_fill_manual(values = pal_ftype, name = "") + #wacolors::scale_fill_wa_d(palette = "larch", name = "") +
  facet_wrap(~rgn, ncol = 1) +
  labs(title = "Fishery recoveries")
) + (
rmis_rr_rgn |>  
  filter(f_type == "escp") |> 
  ggplot(aes(run_year, n, fill = f_type)) + 
  geom_col(show.legend = F) +
  scale_x_continuous("", n.breaks = 8, guide = guide_axis(n.dodge = 2)) +
  scale_y_continuous("", labels = scales::comma) +
  scale_fill_manual(values = pal_ftype) + #wacolors::scale_fill_wa_d(palette = "uw") +
  facet_wrap(~rgn, ncol = 1) +
  labs(title = "Escapement recoveries") #+ theme(legend = element_blank())
) +
  plot_layout(nrow = 1, guides = 'collect')
```

#### Recent Year Recovery Counts for Assessed Stocks

Relevance to the [PST](#PST) and the number of [CWT](#CWT) recoveries during the years 2010-2019 were primary determinants of the focal stocks assessed and presented [here]{#FocalCWT}. 

```{r results='asis'}
# rmis_n_recs_per_stock <- rmis_rr |> 
#   drop_na(estimated_number, StockID) |> 
#   filter(between(run_year, 2010, 2019)) |> 
#   count(StockID, f_type) |> 
#   drop_na(f_type) |> 
#   pivot_wider(names_from = f_type, values_from = n) |> 
#   left_join(fram$mort |> distinct(StockID, StockLongName), by = "StockID") |> 
#   drop_na(StockLongName) |> 
#   select(StockID, StockLongName, everything())

rmis_n_recs_per_stock |> 
  filter(StockID %in% c(20, 34, 48, 58, 68, 84, 96, 134, 142, 148, 152, 164, 166, 176, 210, 226, 230)) |> 
  mutate(total = escp + net + sport + troll) |> 
  arrange(desc(total)) |> 
  gt::gt() |> 
  gt::sub_missing() |> 
  gt::fmt_number(columns = escp:total, decimals = 0)|>
  gt::tab_header(title = "Table 6.1. Number of CWT Recoveries for Assessed Stocks, 2010-1019")
```


### Exploitation Rate Calculations {#Pseudo}

Per stock, per year, a ["pseudo-ER"](#Pseudo) was calculated by dividing the annual per-fishery landed catch ([FRAM](#FRAM), non-selective + mark-selective, time steps 1 to 5) or `estimated_number` (RMIS, summed over individual recoveries) by the sum of the total across fisheries and escapement. 

Equation 6.1 (for the FRAM case):
$pseudoER_{s,f,y} = \frac{LC_{s,f,y}}{Escp_{s,y} + \sum_{f}LC_{s,f,y}}$

In addition, the "total pseudo-ER" was calculated similarly by dividing the across-fishery total by the sum of the total across fisheries and escapement. 

Equation 6.2 (for the RMIS case): 
$pseudoER_{s,y} = \frac{\sum_{f}EstNum_{s,f,y}}{Escp_{s,y} + \sum_{f}EstNum_{s,f,y}}$

Using identical stocks, years, and fisheries for a direct comparison between RMIS and [FRAM](#FRAM) data sources, these values were then joined by year, stock and fishery.

Where,<br>
  Escp    = Escapement<br>
  EstNum  = Estimated Number of CWT Recoveries from RMIS<br>
  f       = Fishery<br>
  LC      = FRAM Landed Catch<br>
  s       = Stock<br>
  y       = Year

```{r er_objects_included_fisheries, eval=FALSE}
#limiting expansion to 8x, drops several Canadian fishery-years
samp_catch_coef_vals <- samp_catch_coef |>
  filter(Source == "Xpansion", (Action != "Remove  fishery" | is.na(Action))) |>
  select(FisheryID, `2010`:`2019`) |> 
  pivot_longer(cols = -FisheryID, names_to = "run_year", values_to = "coef") |> 
  mutate(
    run_year = as.numeric(run_year),
    coef = if_else(coef > 8, 0, coef),
    coef = if_else(FisheryID %in% c(39, 44, 142, 179, 182, 189), 0, coef)
    )


er_fram <- full_join(
  fram$mort |>
    mutate(
      run_year = as.integer(RunYear), RunYear = NULL,
      FisheryID = if_else(FisheryID %in% fishery_ids_paired, FisheryID+1, as.double(FisheryID)),
      FisheryName = if_else(FisheryID %in% c(fishery_ids_paired, fishery_ids_paired+1),
                            str_sub(FisheryName, 1,8), FisheryName),
      FisheryID = if_else(FisheryID == 98, 97, as.double(FisheryID)),
      FisheryName = if_else(FisheryID == 97, "A7BCDNet", FisheryName)
    ) |> 
    filter(FisheryID %in% fishery_ids) |> 
    group_by(run_year, StockID, StockLongName, FisheryID, FisheryName) |> 
    summarise(
      across(c(LandedCatch, MSFLandedCatch), sum), #sum over timesteps
      fram_lnd_yr_fsh = LandedCatch + MSFLandedCatch,
      .groups = "drop") |> 
    left_join(
      samp_catch_coef_vals |> filter(coef == 0),
      by = c("FisheryID", "run_year")
    ) |> 
    mutate(fram_lnd_yr_fsh = if_else(!is.na(coef), coef * fram_lnd_yr_fsh, fram_lnd_yr_fsh)) |> 
    group_by(run_year, StockID, StockLongName) |> 
    mutate(
      LandedCatch = NULL, MSFLandedCatch = NULL,
      fram_lnd_yr_tot = sum(fram_lnd_yr_fsh)
    ) |> 
    ungroup()
  ,
  fram$escp |> 
    rename(fram_escp = escp) |> 
    mutate(run_year = as.integer(RunYear), RunYear = NULL) |> 
    group_by(run_year, StockID) |> 
    summarise(across(fram_escp, sum), .groups = "drop") 
  ,
  by = c("run_year", "StockID")
) |>
  mutate(
    fram_abnd = fram_lnd_yr_tot + fram_escp,
    fram_er_fsh = fram_lnd_yr_fsh / fram_abnd,
    fram_er_tot = fram_lnd_yr_tot / fram_abnd
  ) |> 
  select(run_year, everything()) |> 
  filter(between(run_year, 2010, 2019))


# #large majority of non-escp NA fisheries are 91 from AK
# rmis_rr |> 
#   drop_na(estimated_number, StockID) |>
#   filter(is.na(f_type)) |> count(fishery, rlsp)

er_rmis <- full_join(
  #analogous to FRAM mortality table
  rmis_rr |> 
    drop_na(estimated_number, StockID, FisheryID) |>
    mutate(
      FisheryID = if_else(FisheryID %in% fishery_ids_paired, FisheryID+1, as.double(FisheryID)),
      FisheryName = if_else(FisheryID %in% c(fishery_ids_paired, fishery_ids_paired+1),
                            str_sub(FisheryName, 1,8), FisheryName),
      FisheryID = if_else(FisheryID == 98, 97, as.double(FisheryID)),
      FisheryName = if_else(FisheryID == 97, "A7BCDNet", FisheryName)
    ) |> 
    filter(
      age == 3, !between(FisheryID, 24, 32),
      between(run_year, min(er_fram$run_year), max(er_fram$run_year)),
      StockID %in% stock_ids,
      FisheryID %in% fishery_ids
    ) |> #count(sample_type)
    group_by(run_year, StockID, f_type, FisheryID, FisheryName) |> 
    summarise(rmis_est_num_yr_fsh = sum(estimated_number), .groups = "drop") |> 
    left_join(
      samp_catch_coef_vals, 
      by = c("FisheryID", "run_year")
    ) |>
    mutate(rmis_est_num_yr_fsh = if_else(!is.na(coef), rmis_est_num_yr_fsh * coef, rmis_est_num_yr_fsh)) |> 
    group_by(run_year, StockID) |> 
    mutate(rmis_est_num_yr_tot = sum(rmis_est_num_yr_fsh)) |> 
    ungroup()
  ,
  #analogous to FRAM escapement table
  rmis_rr |> 
    drop_na(estimated_number, StockID) |> 
    left_join(
      samp_catch_coef_vals |> filter(between(FisheryID, 24, 32)), 
      by = c("FisheryID", "run_year")
    ) |>
    mutate(estimated_number = if_else(!is.na(coef), estimated_number * coef, estimated_number)) |> 
    filter(
      age == 3, (f_type == "escp" | between(FisheryID, 24, 32)),
      between(run_year, min(er_fram$run_year), max(er_fram$run_year)),
      StockID %in% stock_ids
    ) |>
    group_by(run_year, StockID) |> 
    summarise(rmis_escp = sum(estimated_number), .groups = "drop")
  ,
  by = c("run_year", "StockID")
  ) |> 
  mutate(
    run_year = as.integer(run_year),
    rmis_abnd = rmis_est_num_yr_tot + rmis_escp,
    rmis_er_fsh = rmis_est_num_yr_fsh / rmis_abnd,
    rmis_er_tot = rmis_est_num_yr_tot / rmis_abnd
  )


er_fram_rmis <- full_join(
  er_fram |> select(-StockLongName, -FisheryName, -coef),
  er_rmis |> select(-FisheryName, -coef),
  by = c("run_year", "StockID", "FisheryID")
) |> 
  left_join(fram$mort |> distinct(StockID, StockLongName), by = "StockID") |> 
  left_join(fram$mort |> distinct(FisheryID, FisheryName), by = "FisheryID") |> 
  mutate(
    fidnm = paste(FisheryID, FisheryName, sep = ":") |> fct_reorder(.x = FisheryID, .fun = min),
    sidnm = paste(StockID, StockLongName, sep = ":") |> fct_reorder(.x = StockID, .fun = min),
    d_yr_fsh = fram_lnd_yr_fsh - rmis_est_num_yr_fsh,
    d_yr_tot = fram_lnd_yr_tot - rmis_est_num_yr_tot,
    d_escp_fram_rmis = fram_escp - rmis_escp,
    d_er_fsh = fram_er_fsh - rmis_er_fsh,
    d_er_tot = fram_er_tot - rmis_er_tot
    ) |> 
  filter(between(run_year, 2010, 2019))

# #a trace smattering of NA fisheries when mapped stock-year escapement recoveries without fishery recoveries 
# #prevents comparison even though FRAM has both mort and escp
# er_fram_rmis |> filter(is.na(FisheryID))
# er_rmis |> filter(run_year==2014, StockID==194) #only escapement
# rmis_rr |> filter(run_year==2014, StockID==194) |> count(f_type) #confirming only escapement
# fram$mort |> filter(RunYear==2014, StockID==194)
# fram$escp |> filter(RunYear==2014, StockID==194)
# er_fram_rmis |> filter(run_year==2014, StockID==194)

```


## Results

The [FRAM](#FRAM) and [RMIS](#RMIS) [pseudo-ERs](#Pseudo) can be compared from "stock-centric" and "fishery-centric" perspectives. 

Differences across individual fisheries for single stocks are presented first, followed by differences organized by fishery, across stocks.

The followings sections contain several plots per stock. First, bar charts of "total" [pseudo-ER](#Pseudo) per year illustrate the values summed across fishery for each data source ([FRAM](#FRAM) and [RMIS](#RMIS)). The number of years in which the absolute difference was smaller than each of several thresholds is noted. The overall pattern of differences offers an indication of whether the stock might be affected by sampling issues in escapement and/or fisheries catches. 

These plots are supplemented by an illustration of the per-year fishery-specific [pseudo-ER](#Pseudo) differences ([FRAM](#FRAM) minus [RMIS](#RMIS)), which can help to clarify whether "total" differences are driven primarily by one or two key fisheries or are shared across many fisheries (with the latter perhaps indicating issues related to escapement sampling or input estimates). In these stock-year-fishery-specific plots, the points show both the relative magnitude of differences (as point size and color saturation) and whether [FRAM](#FRAM) (orange, upward pointing) or [RMIS](#RMIS) (purple, downward pointing) produced the larger value. In addition, values are shown as black dots for combinations in which a [FRAM](#FRAM)-generated impact estimate lacks a corresponding record in [RMIS](#RMIS). Arithmetic means across fisheries per year (top of columns) and across years per fishery (right-most column) also suggest characteristic tendencies in each of these dimensions (noting the caveat that the means of opposing signs can approach 0 despite non-zero positive and negative values). 

Following these [pseudo-ER](#Pseudo) plots, two panels display the ranked proportional fishery impacts from each data source. These are based on the same subset of values, but emphasize the relative magnitudes among fisheries, as bars of greater and lesser lengths, or the differences in rank as "ladder rungs". The bars represent the proportion of the total landed catch that can be attributed to a fishery. This value is calculated as the landed catch in a fishery divided by the total catch summed over all fisheries for each method ([FRAM](#FRAM),[RMIS](#RMIS)). For the "ladder panel", equivalence in rank is shown as a flat or level line between the [FRAM](#FRAM) and [RMIS](#RMIS) columns, with more steeply angled lines indicating greater divergence in ranks. For example, line sloping sharply down and right could highlight a fishery in which [FRAM](#FRAM) assigned a high percentage of cumulative catch mortality but which lacked a corresponding degree of recoveries in [RMIS](#RMIS). In these plots, fisheries with assigned [FRAM](#FRAM) impacts but without any recoveries throughout 2010-2019 appear as light gray lines, typically sloping up from [FRAM](#FRAM) ranks that do not exist for [RMIS](#RMIS). Because there are usually several RMIS fisheries with zero recoveries where FRAM assigns an impact, the final RMIS fishery label (at 0.0%) becomes difficult to read. The size of the dots in the "ladder panel" are proportional to the magnitude of the fishery. 


```{r gg_stk_function, eval=FALSE}

gg_stk <- function(sid = 20, min_fram_landed = 10){
    tot <- filter(er_fram_rmis, StockID == sid) |> 
      drop_na(d_er_tot) |>
      mutate(run_year = as.character(run_year)) |> 
      distinct(sidnm, run_year, FRAM = fram_er_tot, RMIS = rmis_er_tot) |>
      pivot_longer(-c(sidnm, run_year), names_to = "var", values_to = "val")
    
    tot_d_yr <- filter(er_fram_rmis, StockID == sid) |> 
      drop_na(d_er_tot) |>
      distinct(sidnm, run_year, d_er_tot) |> 
      pluck("d_er_tot") |> abs()
      
    tot_d_yr = map_chr(
          c(0.01, 0.05, 0.1),
          ~paste0("d_yr<",.x,": ", sum(tot_d_yr < .x))
          ) |> paste(collapse = ";  ")
 
    fsh <- bind_rows(
      filter(er_fram_rmis, StockID == sid, fram_lnd_yr_fsh >= min_fram_landed) |> 
        mutate(run_year = as.character(run_year)),
      #across years, end of fishery row
      filter(er_fram_rmis, StockID == sid, fram_lnd_yr_fsh >= min_fram_landed) |> 
        mutate(
          d_er_fsh = if_else(is.na(rmis_er_fsh), fram_er_fsh, d_er_fsh),
          d_er_fsh = if_else(is.na(fram_er_fsh), -1*rmis_er_fsh, d_er_fsh)
          ) |> 
        group_by(sidnm, fidnm) |>
        summarise(
          run_year = "mean",
          d_er_fsh = mean(d_er_fsh, na.rm=T), .groups = "drop"),
      #across fisheries, bottom of year col
      filter(er_fram_rmis, StockID == sid, fram_lnd_yr_fsh >= min_fram_landed) |> 
        mutate(
          d_er_fsh = if_else(is.na(rmis_er_fsh), fram_er_fsh, d_er_fsh),
          d_er_fsh = if_else(is.na(fram_er_fsh), -1*rmis_er_fsh, d_er_fsh)
        ) |> 
        mutate(run_year = as.character(run_year)) |> 
        group_by(sidnm, run_year) |>
        summarise(
          fidnm = "0: annual mean",
          d_er_fsh = mean(d_er_fsh, na.rm=T), .groups = "drop")
      ) 
    
    fsh <- bind_rows(
      fsh,
      fsh |> filter(fidnm == "0: annual mean") |> 
        group_by(sidnm, fidnm) |> 
        summarise(run_year = "mean", d_er_fsh = mean(d_er_fsh, na.rm=T), .groups = "drop")
    )
    
    #fishery ranks as all-years summed percentages: sum per fishery / sum all fisheries
    #meant to stabilize year-to-year fluctuations in sampling etc.
    rnk <- filter(er_fram_rmis, StockID == sid, fram_lnd_yr_fsh > 0) |> 
      group_by(sidnm, fidnm) |> 
      summarise(
        fram_fsh = sum(fram_lnd_yr_fsh, na.rm = T),
        rmis_fsh = sum(rmis_est_num_yr_fsh, na.rm=T),
        .groups = "drop") |> 
      group_by(sidnm) |> 
      mutate(
        fram_tot = sum(fram_fsh),
        rmis_tot = sum(rmis_fsh),
        fram_fsh_tot = fram_fsh / fram_tot,
        rmis_fsh_tot = rmis_fsh / rmis_tot,
        across(ends_with("_fsh_tot"), list(rnk = ~min_rank(-.))),
        d_rnk = fram_fsh_tot_rnk - rmis_fsh_tot_rnk,
        d_rnk_color = case_when(
          fram_fsh_tot_rnk == rmis_fsh_tot_rnk ~ "darkgreen", #d_rnk == 0
          fram_fsh_tot_rnk < rmis_fsh_tot_rnk  ~ "orange", #d_rnk < 0
          fram_fsh_tot_rnk > rmis_fsh_tot_rnk  ~ "purple" #d_rnk > 0
        ),
        d_rnk_color = if_else(round(fram_fsh_tot, 3) == 0 & rmis_fsh_tot == 0, "grey50", d_rnk_color)
      ) |> 
      ungroup() |> 
      arrange(sidnm, fram_fsh_tot_rnk) |> 
      mutate(across(c(fram_fsh_tot, rmis_fsh_tot), ~round(.,3))) |> 
      select(sidnm, fidnm, contains("fsh_tot"), contains("d_rnk")) |> 
      unite("FRAM", starts_with("fram_")) |> 
      unite("RMIS", starts_with("rmis_")) |> 
      pivot_longer(c(FRAM, RMIS), names_to = "var", values_to = "val") |> 
      separate(val, into = c("pct", "rnk"), sep = "_") |> 
      mutate(
        pct = as.numeric(pct),
        rnk = factor(as.numeric(rnk))
      )
    
    snm <- tot$sidnm[1]
      
    gg_tot_bars <- bind_rows(
        tot,
        tot |> group_by(sidnm, var) |> summarise(run_year = "mean", across(val, mean), .groups = "drop") 
        ) |> 
      ggplot(aes(run_year, val, fill=var)) +
      geom_col(position = position_dodge()) +
      geom_text(aes(label = if_else(is.na(val), "", paste0(round(100*val, digits = 1), "%"))),
                position = position_dodge2(width = 1, preserve = "single"),
                size = 2, vjust = -0.5) +
      scale_fill_grey() +
      scale_x_discrete("", guide = guide_axis(n.dodge = 2)) +
      scale_y_continuous("", labels = scales::percent) +
      labs(
        title = "Annual total pseudo-ER, FRAM & RMIS",
        subtitle = tot_d_yr) +
      theme(legend.position = "right")
    
    gg_fsh <- fsh |>
      mutate(
        fidnm = factor(fidnm, levels = c(levels(er_fram_rmis$fidnm), "0: annual mean")),
        shp = case_when(
          d_er_fsh > 0 ~ "FRAM > RMIS",
          d_er_fsh < 0 ~ "FRAM < RMIS",
          d_er_fsh == 0 ~ "FRAM = RMIS",
          is.na(rmis_er_fsh) ~ "FRAM-only",
          is.na(fram_er_fsh) ~ "RMIS-only"
        ),
        #sze = if_else(!is.na(d_er_fsh), abs(d_er_fsh), 0)
        sze = if_else(!is.na(d_er_fsh), abs(d_er_fsh), fram_er_fsh)
      ) |>
      ggplot(aes(run_year, fidnm)) +
      geom_point(aes(shape = shp), size = 1) +
      geom_point(aes(shape = shp, size = sze, fill = d_er_fsh)) +
      geom_text(aes(label = if_else(is.na(d_er_fsh), "", paste0(round(100*d_er_fsh, digits = 0), "%"))), hjust = -0.5, size = 3) +
      scale_shape_manual(name = "", values = c("FRAM > RMIS" = 24, "FRAM < RMIS" = 25, "FRAM = RMIS" = 22, "FRAM-only" = 20, "RMIS-only" = 3)) +
      scale_size("ER diff. abs.", range = c(1.5, 5)) +
      scale_fill_fermenter(type = "div", palette = "PuOr", breaks = c(-0.1, -0.05, 0, 0.05, 0.1), limits = c(-.2, .2)) +
      geom_hline(yintercept = length(levels(factor(fsh$fidnm))) - 0.5) + 
      guides(
        fill = guide_colorbar(title="ER difference"),
        x = guide_axis("Run Year", n.dodge = 2),
        y = guide_axis("")) +
      facet_wrap(~sidnm, ncol = 2, scales = "free_y") +
      labs(
        title = "Per-year per-fishery pseudo-ER difference",
        subtitle = paste("FRAM - RMIS, min. annual FRAM landed >=", min_fram_landed, "coho")
      )

    gg_rnk_col <- rnk |>
      ggplot(aes(fct_reorder(fidnm, pct, max), pct, group = var, fill = var)) +
      geom_col(position = position_dodge(width = 0.3), alpha = 0.5) +
      coord_flip() +
      scale_fill_manual(values = c(FRAM = "orange", RMIS = "purple")) + #scale_fill_grey() +
      theme(legend.position = "bottom")

    gg_rnk_ladder <- rnk |> 
      ggplot(aes(var, rnk, color = d_rnk_color, alpha = pct, group = fidnm)) + 
      geom_point(aes(size = pct), show.legend = F) + 
      geom_line(show.legend = F) +
      geom_text(
        data = filter(rnk, var == "FRAM"),
        aes(label = paste(fidnm, scales::percent(pct, accuracy = 0.1), sep = "  ")),
        size= 2.2, hjust = 1, vjust = -0.2, show.legend = F) +
      geom_text(
        data = filter(rnk, var == "RMIS"),
        aes(label = paste(fidnm, scales::percent(pct, accuracy = 0.1), sep = "  ")),
        size= 2.2, hjust = -0.1, vjust = -0.2, show.legend = F) +
      scale_x_discrete("") +
      scale_y_discrete("Rank", limits = rev) +
      scale_size(range = c(0.2, 4)) +
      scale_color_manual(values = c(darkgreen = "darkgreen", orange = "orange", purple = "purple")) +
      scale_alpha_continuous(range = c(0.4, 1)) #+labs(title = "Proportional ranks, 2010-2019")

    return(list(
      pseudo_er = gg_tot_bars + gg_fsh +
        plot_layout(ncol = 1, heights = c(1/4, 3/4)) +
        plot_annotation(title = snm),
      rel_rank = gg_rnk_col + gg_rnk_ladder + 
        plot_layout(ncol = 2) +
        plot_annotation(title = paste(snm, "Proportional ranking, 2010-2019"))
    ))
    
}

```

### Puget Sound

#### Skagit

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 2, 5 and 7 years respectively
  - FRAM value larger in 2013 and 2014, but otherwise RMIS value larger
  - RMIS escapement may be low in some years; in 2019 RMIS escp < fisheries, while FRAM escp ~2x fisheries
  
- **Fishery-specific**
  - annually varied +/- differences for most fisheries
  - contrasting differences did not fully offset in PS sport: lower FRAM for A10, A5, A9 and A8-2 not compensated by A8-1 and A7
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest impacts to A5, and FRAM assigns third-most impacts to A9, the second ranked RMIS fishery
  - RMIS indicates a considerably higher proportion of impact in A10, A8-2 and A2 than in FRAM
  - RMIS indicates a lower proportion of impact in 4/4B Treaty troll, A6, A8-1, and 8A net than in FRAM
  - Counter misalignment in 8-1 & 8-2 may be due to errors in datasets and/or mapping<br>
  
[**Figures 6.3.1.a-d. Skagit Pseudo ERs and Proportional Fishery Rankings**]{#Skagit1}
```{r}
g <- gg_stk(20) 
g$pseudo_er#; ggsave(file.path("results","f_20_Skagit.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_20_Skagit_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==20, between(run_year, 2017, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)
```


#### Tulalip

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.1, 0.05 and 0.01 in 1, 5 and 7 years respectively
  - FRAM value larger in 2016 only
  - RMIS escapement may be low in some years; 2019 escp. missing entirely preventing [ER](#ER) calculations (note no 2019 total column)
  
- **Fishery-specific**
  - most fisheries show negative differences where comparison possible. This is likely caused by missing escapement accounting and sampling at the hatchery and in the [extreme terminal](#ET) area.
  - FRAM appears to overestimate Tulalip impacts in Tulalip Bay (area 8D), while RMIS assigns more impacts to fisheries outside of the bay. 
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest impacts to 8D net (first) and 8A net (second)
  - both have A5 and A9 sport as third or fourth largest, and 8-2 as fifth; other fisheries are relatively minor in both <br>
  
[**Figures 6.3.2.a-d. Tulalip Pseudo ERs and Proportional Fishery Rankings**]{#Tulalip1}
```{r}
g <- gg_stk(34)
g$pseudo_er#; ggsave(file.path("results","f_34_Tulalip.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_34_Tulalip_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==34, between(run_year, 2017, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)
```


#### Quilcene

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.1, 0.05 and 0.01 in 0, 2 and 4 years respectively
  - FRAM value consistently larger; RMIS larger in 2010 only
  - Fishery recoveries may be missing in some years; for example
    - in 2019 RMIS escp is ~3x fisheries while FRAM is ~1:1
    - in 2012 FRAM fisheries ~3x escp while RMIS fisheries < escp

- **Fishery-specific**
  - possibly mapping or other dataset errors affecting Hood Canal net fisheries 12A, 9A and 12/12B (see rank comment)
  - possible mis-assignments:
    - sport: FRAM impacts higher in A5, A10, and A6, while A9, A12, A4 and A2 lower than RMIS
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest impacts to 12A net
  - 9/9A net second for RMIS but 5th for FRAM, while 12/12B is third for FRAM but minor for RMIS <br>
  
[**Figures 6.3.3.a-d. Quilcene Pseudo ERs and Proportional Fishery Rankings**]{#Quilcene1}
```{r}
g <- gg_stk(48)
g$pseudo_er#; ggsave(file.path("results","f_48_Quilcene.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_48_Quilcene_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==48, between(run_year, 2010, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)
```


#### George Adams

- **Total across fisheries** 
  - no years close; RMIS likely missing fishery recoveries and/or mapping error

- **Fishery-specific**
  - many positive differences may be related to missing fishery recoveries and/or possibly mapping errors (e.g., among Hood Canal net fisheries)
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest impacts to 12CD net and large portions to A5, A9 and A10 sport
  - A12 sport missing from RMIS in most years hinders comparison
  - Area 12/12B net (FRAM higher) and Area 9/9A Net (RMIS higher) countering differences suggests mapping or reporting error<br>
  
[**Figures 6.3.4.a-d. George Adams Pseudo ERs and Proportional Fishery Rankings**]{#GeorgeAdams1}
```{r}
g <- gg_stk(58)
g$pseudo_er#; ggsave(file.path("results","f_58_GeorgeAdams.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_58_GeorgeAdams_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==58, between(run_year, 2017, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)
```

#### Nisqually

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 4 and 7  years respectively
  - FRAM value larger in 2012-4 and 2017-8; RMIS larger otherwise; no dominant pattern

- **Fishery-specific**
  - largest differences in 2015 Puget Sound sport, subject to atypical year caveats

- **Rank agreement**
  - both FRAM and current RMIS associate the greatest impacts to A5, A10 and A9 sport (recall in-river net not included)
  - contrasting ranks for A13 (FRAM higher) and A6 (RMIS higher) may be related to sampling intensity<br>
  
[**Figures 6.3.5.a-d. Nisqually Pseudo ERs and Proportional Fishery Rankings**]{#Nisqually1}
```{r}
g <- gg_stk(68)
g$pseudo_er#; ggsave(file.path("results","f_68_Nisqually.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_68_Nisqually_rel_rank.png"), width = 9, height = 10)
```

#### Puyallup

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 3 and 5 years respectively
  - FRAM value larger in 2011 and 2012 only
  - Larger RMIS values may be related to missing escapement recoveries in some years
  
- **Fishery-specific**
  - year to year fishery differences primarily in context of larger RMIS [pseudo-ERs](#Pseudo)
  - contrasting differences among Puget Sound sport fisheries in some years suggests FRAM assignments may not match particular year sampled patterns (e.g, A10 vs A6 or A9)
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest [(pre-terminal)](#PTerm) impacts to A10, A5 and A9 sport (first, second and third respectively), with A6 and A11 also appreciable <br>
  
[**Figures 6.3.6.a-d. Puyallup Pseudo ERs and Proportional Fishery Rankings**]{#Puyallup1}
```{r}
g <- gg_stk(84)
g$pseudo_er#; ggsave(file.path("results","f_84_Puyallup.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_84_Puyallup_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==84, between(run_year, 2010, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)

#escp |> filter(StockID == 84, year == 2019) #2019 5869+1583=7452
```

#### Green 

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 4 and 8 years respectively
  - FRAM had larger value in all years but 2015
  - Fishery recoveries may be missing in some years; for example
    - in 2017 RMIS fisheries ~1/2 escp when FRAM ~1:1, vs 2018 and 2019 both similar
  
- **Fishery-specific**
  - As for total, larger differences may be related to missing RMIS fishery recoveries
  - Larger FRAM [pseudo-ER](#Pseudo) values in Area 10A net in 2017 align with possibility of fewer recoveries (or none, as in 2019)
  
- **Rank agreement**
  - Despite [pseudo-ER](#Pseudo) differences, both FRAM and current RMIS associate the greatest proportion of impacts to the same three Puget Sound sport fisheries (A10, A5, A9)
  - Area 10A net and Area 4/4B troll account for a relatively lower proportion of RMIS recoveries, possibly related to differences in sampling intensity.<br>
  
[**Figures 6.3.7.a-d. Green Pseudo ERs and Proportional Fishery Rankings**]{#Green1}
```{r}
g <- gg_stk(96)
g$pseudo_er#; ggsave(file.path("results","f_96_Green.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_96_Green_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==96, between(run_year, 2017, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)

er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |>
  filter(StockID==96, between(run_year, 2010, 2019), FisheryID %in% c(122)) |>
  select(run_year, sidnm, fidnm, fram_lnd_yr_fsh:fram_er_tot, rmis_est_num_yr_fsh:rmis_er_tot) |> 
  mutate(
    fram_fsh_tot = fram_lnd_yr_fsh / fram_lnd_yr_tot,
    rmis_fsh_tot = rmis_est_num_yr_fsh / rmis_est_num_yr_tot
    )

```

### WA Coast

#### Quillayute Fall

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 1 and 4 years respectively
  - RMIS value larger in all years; differences may be related to missing escapement recoveries and/or mapping errors
  
- **Fishery-specific**
  - largest differences concentrated in A2 sport, the closest largest [pre-terminal](#PTerm) fishery; potentially accounted for by differences in escapement recoveries and/or [BkFRAM](#BkFRAM) targets and flagging
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to the same three fisheries: A2 and A1 sport and Area 4/4B troll<br>
  
[**Figures 6.3.8.a-d. Quillayute ERs and Proportional Fishery Rankings**]{#Quillayute1}
```{r}
g <- gg_stk(134)
g$pseudo_er#; ggsave(file.path("results","f_134_QuillayuteFall.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_134_QuillayuteFall_rel_rank.png"), width = 9, height = 10)
```

#### Queets

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 1, 2 and 4 years respectively, comparison possible in 6 years
  - larger RMIS value in 2014 and 2015, also largest differences in these years
    - RMIS fishery recoveries may be missing in some years: notably higher total estimated number in 2014, driven by A2 sport
  - RMIS escapement entirely missing in several years (2012, 2013, 2016, 2019) and possibly missing in others
  - releases may be a contributing factor? varied mark and/or tagging rates?
  
- **Fishery-specific**
  - contrasting differences may be related to errors in dataset or mapping
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to the same three fisheries: A2 and A1 sport and Area 4/4B troll
  - no recent recoveries in SW VI troll despite FRAM assigning ~7% of impacts over 2010-2019 <br>
  
[**Figures 6.3.9.a-d. Queets ERs and Proportional Fishery Rankings**]{#Queets1}
```{r}
g <- gg_stk(142)
g$pseudo_er#; ggsave(file.path("results","f_142_Queets.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_142_Queets_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==142, between(run_year, 2010, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)

# #SW VI Trl only 2 recoveries, both in 2005, despite appreciable proportion in FRAM
# rmis_rr |> filter(StockID==142, FisheryID==175)
#at least in 2014, this fishery contributes appreciably to total catch and ER
er_fram_rmis |>
  filter(StockID==142, between(run_year, 2014, 2015), FisheryID==175)

er_fram_rmis |>
  filter(StockID==142, between(run_year, 2014, 2015)) |> 
  arrange(desc(rmis_est_num_yr_fsh))

```

#### Quinault

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 1, 5 and 7 years respectively
  - RMIS value larger in all years; differences may be related to missing escapement recoveries and/or BKFRAM errors

- **Fishery-specific**
  - largest differences concentrated in A2 sport, the closest largest [pre-terminal](#PTerm) fishery; potentially accounted for by differences in escapement recoveries and/or [BkFRAM](#BkFRAM) targets and flagging
  - consistently higher FRAM [pseudo-ER](#Pseudo) values in A3 and A5 sport despite lower FRAM total [pseudo-ER](#Pseudo); neither a major driver of total but may indicate shift in impact since [base period](#BP)
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to the same three fisheries: A2 and A1 sport and Area 4/4B troll
  - Area 3 sport, Newport sport and A5 sport show considerably lower proportion of recent RMIS recoveries relative to FRAM catch assignments
  - Area 2 and 3 NT troll show higher proportion of recent RMIS recoveries relative to FRAM catch assignments<br>
  
[**Figures 6.3.10.a-d. Quinault ERs and Proportional Fishery Rankings**]{#Quinault1}
```{r}
g <- gg_stk(148)
g$pseudo_er#; ggsave(file.path("results","f_148_Quinault.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_148_Quinault_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==148, between(run_year, 2010, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp)

#A3 sport ER over-represented?
er_fram_rmis |>
  filter(StockID==148, between(run_year, 2017, 2019), FisheryID==40)
#while A2 NT troll under-represented?
er_fram_rmis |>
  filter(StockID==148, between(run_year, 2017, 2019), FisheryID==35)
```


#### Chehalis

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 1, 6 and 10 years respectively
  - varied year to year in source of larger [pseudo-ER](#Pseudo) value; differences small in several years

- **Fishery-specific**
  - Area 4/4B troll consistently assigned higher FRAM [pseudo-ER](#Pseudo) than present in RMIS; not a major driver of total, but difference in relative rank may indicate shift in impact since [base period](#BP)
  - Grays Harbor net also a consistently assigned higher FRAM [pseudo-ER](#Pseudo) than present in RMIS but differences possibly more related to escapement errors given small rank difference
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to Area 2 sport, and similar high relative ranks to Grays Harbor net, A7 sport and Grays sport
  - Willapa Bay net exhibits higher proportion of recent RMIS recoveries relative to FRAM catch assignments, but year-fishery [pseudo-ERs](#Pseudo) are not uniformly higher in RMIS<br>
  
[**Figures 6.3.11.a-d. Chehalis ERs and Proportional Fishery Rankings**]{#Chehalis1}
```{r}
g <- gg_stk(152)
g$pseudo_er#; ggsave(file.path("results","f_152_Chehalis.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_152_Chehalis_rel_rank.png"), width = 9, height = 10)
```

```{r}
#4/4B troll
#not heavily driven by samp coef: samp_catch_coef |> filter(FisheryID==43)
er_fram_rmis |>
  filter(StockID==152, between(run_year, 2010, 2019), FisheryID==43) |> 
  select(sidnm, fidnm, run_year, contains("er_fsh"))

#also not heavily driven by samp coef:samp_catch_coef |> filter(FisheryID==50)
er_fram_rmis |>
  filter(StockID==152, between(run_year, 2010, 2019), FisheryID==50) |> 
  select(sidnm, fidnm, run_year, contains("er_fsh"))

#Willapa impacts under-represented?
#not subject to major coef: samp_catch_coef |> filter(FisheryID==47)
#relatively minor overall element of catch, and not uniformly higher in RMIS, but also not driven by 2015
er_fram_rmis |>
  filter(StockID==152, between(run_year, 2010, 2019), FisheryID==47) |> 
  select(sidnm, fidnm, run_year, contains("er_fsh"))

```

#### Willapa

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 3 and 5 years respectively
  - FRAM value larger in 2014 and 2019
  - differences possibly related to missing recoveries and/or BKFRAM errors
    - 2017 fewest RMIS escapement recoveries and lowest ratio to FRAM escapement
    - 2019 RMIS escapement ~3.4x FRAM escapement and RMIS Willapa net recovery ~1.6x FRAM catch

- **Fishery-specific**
  - differences dominated by Willapa net, with RMIS pseudo-ER substantially larger in earlier years (2010-12), but FRAM value larger after 2014, excepting 2017
    - note FRAM catch atypically low in 2019 in Willapa and across fisheries, while RMIS `Estimated Number` within prior range but higher than FRAM; suggests possible BKFRAM or catch assignment issue
    - contrast with 2014 in which escapement recovery ratio typical, but very large catches may have affected fishery sampling (second lowest ratio of RMIS fishery est_mun to FRAM catch)
    - possibly also missing fishery recoveries in 2018, with lowest ratio for fishery (est_mun/catch)

- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to the same three fisheries: Willapa net, A2 and A1 sport; Willapa sport also appreciable in both<br>
  
[**Figures 6.3.12.a-d. Willapa ERs and Proportional Fishery Rankings**]{#Willapa1}
```{r}
g <- gg_stk(164)
g$pseudo_er#; ggsave(file.path("results","f_164_Willapa.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_164_Willapa_rel_rank.png"), width = 9, height = 10)
```

```{r}
#not subject to major coef: samp_catch_coef |> filter(FisheryID==47)
er_fram_rmis |>
  filter(StockID==164, between(run_year, 2010, 2019), FisheryID==47) |> 
  select(sidnm, fidnm, run_year, contains("_yr_fsh"), contains("_yr_tot"), contains("escp"), contains("er_fsh")) |> 
  mutate(
    ratio_fsh = rmis_est_num_yr_fsh / fram_lnd_yr_fsh,
    ratio_escp = rmis_escp / fram_escp
    )

#note 2019 fram_escp 885 reflects TargetFlag==2 on BackwardsFRAM.TargetEscAge3 of 23995
```


### Columbia

#### Columbia River Earlies

[Note](#CR): We accounted for modeling river run size in [BkFRAM](#BkFRAM) rather than escapement by adding RMIS freshwater fishery recoveries to escapement.

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 2, 7 and 8 years respectively
  - differences larger since 2016, with RMIS values more consistently larger 

- **Fishery-specific**
  - largest difference in Buoy10 sport, but do not account for post-2016 larger RMIS total [pseudo-ER](#Pseudo)
  - Newport sport shifts from mostly larger FRAM values to larger RMIS after 2016
    - 2018 RMIS fishery [pseudo-ER](#Pseudo) ~10.8% considerably larger than previous and FRAM
    - 2019 RMIS fishery [pseudo-ER](#Pseudo) ~10.8% persists and FRAM nearly identical (with smaller difference in 2019 total)
  
- **Rank agreement**
  - both FRAM and current RMIS associate the greatest proportion of impacts to the same 6 fisheries
    - very little indication in recent RMIS recoveries of different distribution of impacts relative to FRAM<br>

[**Figures 6.3.13.a-d. Columbia Earlies ERs and Proportional Fishery Rankings**]{#ColumbiaEarlies1}
```{r}
g <- gg_stk(166)
g$pseudo_er#; ggsave(file.path("results","f_166_ColEarly.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_166_ColEarly_rel_rank.png"), width = 9, height = 10)
```


```{r col_early_eda, eval=FALSE}
#not subject to major coef: samp_catch_coef |> filter(FisheryID==23)
er_fram_rmis |>
  filter(StockID==166, between(run_year, 2010, 2019), FisheryID==23) |> 
  select(sidnm, fidnm, run_year, contains("_yr_fsh"), contains("_yr_tot"), contains("escp"), contains("er_fsh")) |> 
  mutate(
    ratio_fsh = rmis_est_num_yr_fsh / fram_lnd_yr_fsh,
    ratio_escp = rmis_escp / fram_escp
    )

#not subject to major coef: samp_catch_coef |> filter(FisheryID==17)
#but 2018 and 2019 jump out
er_fram_rmis |>
  filter(StockID==166, between(run_year, 2010, 2019), FisheryID==17) |> 
  select(sidnm, fidnm, run_year, contains("_yr_fsh"), contains("_yr_tot"), contains("escp"), contains("er_fsh")) |> 
  mutate(
    ratio_fsh = rmis_est_num_yr_fsh / fram_lnd_yr_fsh,
    ratio_escp = rmis_escp / fram_escp
    )


lu$cwt_stk |> filter(StockID==166) |> count(Hatchery_name, Stock_name) |> arrange(desc(n))

#careful summing cwt_1st_mark_count over recoveries - should distinct and join (similar for n_rel) 
rmis_rr |> drop_na(estimated_number) |> filter(StockID==166, age==3) |> #count(hatchery_location_name) |> arrange(desc(n))
  group_by(StockID,
           #tag_code,
           hatchery_location_name, 
           #Stock_name,
           run_year, f_type) |> 
  summarise(
    across(c(estimated_number), sum), .groups = "drop"
  ) |> 
  #filter(str_detect(hatchery_location_name, "DWOR"))
  drop_na(f_type) |> 
  filter(estimated_number > 10) |> 
  filter(f_type != "escp", between(run_year, 2010, 2019)) |> 
  ggplot(aes(run_year, estimated_number, fill = f_type)) + geom_col() + facet_wrap(~hatchery_location_name)

```

#### Columbia River Lates

[Note](#CR): We accounted for modeling river run size in [BkFRAM](#BkFRAM) rather than escapement by adding RMIS freshwater fishery recoveries to escapement.

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 0, 0 and 2 years respectively
  - RMIS value considerably larger in all years; differences likely related to missing escapement recoveries and/or stock mapping errors
    - ratio of RMIS:FRAM for fisheries larger than RMIS:FRAM for escapement in all years
    - largest deviation in these ratios occurred in 2017 and 2019, when total [pseudo-ER](#Pseudo) difference also largest 

- **Fishery-specific**
  - consistency across fisheries points to handling of stock in FRAM and/or RMIS

- **Rank agreement**
  - despite differences in [pseudo-ER](#Pseudo) magnitude, both FRAM and current RMIS associate the greatest proportion of impacts to similar fisheries<br>
  
[**Figures 6.3.14.a-d. Columbia Lates ERs and Proportional Fishery Rankings**]{#ColumbiaLates1}
```{r}
g <- gg_stk(176)
g$pseudo_er#; ggsave(file.path("results","f_176_ColLate.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_176_ColLate_rel_rank.png"), width = 9, height = 10)
```

```{r}
er_fram_rmis |>
  drop_na(rmis_est_num_yr_tot, fram_lnd_yr_tot) |> 
  filter(StockID==176, between(run_year, 2010, 2019)) |> 
  distinct(run_year, sidnm, fram_lnd_yr_tot, fram_escp, rmis_est_num_yr_tot, rmis_escp) |> 
  mutate(
    ratio_fsh = rmis_est_num_yr_tot / fram_lnd_yr_tot,
    ratio_escp = rmis_escp / fram_escp
    )
```

### British Columbia

#### Strait of Georgia Van. Isl.

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 1, 7 and 8 years respectively
  - overall low [pseudo-ER](#Pseudo) values result in small absolute differences but large proportional differences, particularly before 2017
  - missing fisheries recoveries appear to be partially addressed in 2017-19, but FRAM values still larger, suggesting mapping errors and/or missing fishery recoveries

- **Fishery-specific**
  - several Canadian fisheries excluded entirely due to limited sampling
  - compared to GS Mainland, fewer fishery-years with FRAM catch assigned in the absence of any RMIS
  
- **Rank agreement**
  - Of the fisheries included and with recoveries, both FRAM and current RMIS associate substantial proportions to N BC troll and Area 5 sport
  - Area 6-7 Net is assigned a considerably higher proportion of impacts in FRAM than is seen in recent RMIS recoveries
  - Conversely, NW AK troll is assigned a trace proportion of impacts in FRAM but has the third greatest proportion of recent RMIS recoveries<br>

[**Figures 6.3.15.a-d. Strait of Georgia Vancouver Island ERs and Proportional Fishery Rankings**]{#SGVI1}
```{r}
g <- gg_stk(210)
g$pseudo_er#; ggsave(file.path("results","f_210_GeorgiaVI.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_210_GeorgiaVI_rel_rank.png"), width = 9, height = 10)
```

#### Fraser Lower

- **Total across fisheries** 
  - [pseudo-ER](#Pseudo) differences less than 0.01, 0.05 and 0.1 in 2, 5 and 6 years respectively
  - larger differences before 2016 may be related to missing fishery recoveries

- **Fishery-specific**
  - several Canadian fisheries excluded entirely due to limited sampling
  - several fisheries missing recoveries in 2018-2019 prevent comparison
  - overall small magnitude of [pseudo-ER](#Pseudo) post-2016 results in small absolute differences
  
- **Rank agreement**
  - Of the fisheries included and with recoveries, both FRAM and current RMIS associate the most impact to the same three fisheries
  - Area 7BCD Net is assigned a considerably higher proportion of impacts in FRAM than is seen in recent RMIS recoveries
  - Conversely, Area 4, 6 and 2 sport exhibit higher proportion in RMIS than in FRAM<br>
  
[**Figures 6.3.16.a-d. Fraser Lower ERs and Proportional Fishery Rankings**]{#FraserLower1}
```{r}
g <- gg_stk(226)
g$pseudo_er#; ggsave(file.path("results","f_226_FraserLower.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_226_FraserLower_rel_rank.png"), width = 9, height = 10)
```

#### Fraser Upper

- **Total across fisheries** 
  - pseudo-ER differences less than 0.01, 0.05 and 0.1 in 2, 5 and 8 years respectively
  - since 2013, FRAM value larger; varied differences before 2016
  - Limited or missing RMIS fishery and/or escapement recoveries may affect differences

- **Fishery-specific**
  - several Canadian fisheries excluded entirely due to limited sampling, and only a few fishery-years allow comparison for those included
  - overall small magnitude of [pseudo-ER](#Pseudo) post-2016 results in small absolute differences
  
- **Rank agreement**
  - Of the fisheries included and with recoveries, both FRAM and current RMIS associate the most impact to Area 5 sport (but this fishery has the greatest number of years for comparison)
  - However, the relative ranks of many other included fisheries deviate appreciably
    - Area 4/4B troll and Area 7BCD net account for a considerably lower proportion of impact in RMIS than in FRAM
    - Area 2 and 4 sport account for a considerably higher proportion of impact in RMIS than in FRAM, but are driven by only a few years of comparison<br>
    
[**Figures 6.3.17.a-d. Fraser Upper ERs and Proportional Fishery Rankings**]{#FraserUpper1} 
```{r}
g <- gg_stk(230)
g$pseudo_er#; ggsave(file.path("results","f_230_FraserUpper.png"), width = 9, height = 11)
g$rel_rank#; ggsave(file.path("results","f_230_FraserUpper_rel_rank.png"), width = 9, height = 10)
```

### Across Stocks

Viewed across years, the distributions of differences in total [pseudo-ERs](#Pseudo) fell within +/- 5% for many of the focal assessed stocks (dashed vertical lines). As noted above, for the stocks with characteristic differences outside of +/- 10% and larger year to year variation (i.e., wider IQR), sampling issues in either significant fisheries (e.g., Hood Canal stocks) or escapement (e.g., Quillayute) likely contributed to the consistent biases (shown as dark colors further from the 0 line). <br>

**Figure 6.4.** 
```{r gg_total_er_boxes_focal}
er_fram_rmis |> 
  drop_na(d_er_tot) |>
  filter(
    run_year != 2015,
    StockID %in% c(20, 34, 48, 58, 68, 84, 96, 134, 142, 148, 152, 164, 166, 176, 210, 226, 230)
    ) |> 
  group_by(sidnm) |> 
  distinct(sidnm, run_year, fram_er_tot, rmis_er_tot, d_er_tot) |> 
  mutate(
    across(contains("er_tot"), list(med = median)),
    n = n()) |> 
  #filter(n > 9) |> 
  ggplot(aes(fct_rev(sidnm), d_er_tot, fill = d_er_tot_med)) +
  geom_boxplot(aes(weight = n), outlier.shape = NA, varwidth = T, show.legend = F) +
  geom_jitter(aes(color = d_er_tot_med), alpha = 0.6, show.legend = F) +
  coord_flip() +
  scale_y_continuous(labels = scales::percent) +
  geom_hline(yintercept = 0, size = 1.3) + 
  geom_hline(yintercept = c(-0.05, 0.05), linetype = 2) +
  geom_hline(yintercept = c(-0.1, 0.1), linetype = 3) + 
  scale_fill_fermenter(type = "div", palette = "PuOr", breaks = c(-0.1, -0.05, 0, 0.05, 0.1), limits = c(-.2, .2), aesthetics = c("fill", "color")) +
  guides(
    fill = guide_colorbar(title="ER difference"),
    x = guide_axis(""),
    y = guide_axis("")) +
  labs(
    title = "Differences in Annual Total Pre-terminal Pseudo-ER", subtitle = "Excluding 2015; FRAM - RMIS, positive value if FRAM larger"
  )

```

```{r total_er_annual_gt, eval=FALSE}
er_fram_rmis |> 
  drop_na(d_er_tot) |>
  group_by(sidnm) |> 
  distinct(sidnm, run_year, fram_er_tot, rmis_er_tot, d_er_tot) |>
  ungroup() |>
  arrange(sidnm) |> #print(n=50)
  gt(groupname_col = "sidnm", rowname_col = "run_year") |> 
  summary_rows(groups = T, fns = "mean", formatter = fmt_percent, decimals = 1) |> 
  fmt_percent(contains("_er"), decimals = 1) |> 
  tab_style(
    style = list(cell_fill(color = "purple", alpha = 0.5)),
    locations = cells_body(rows = d_er_tot < -0.05)) |> 
  tab_style(
    style = list(cell_fill(color = "orange", alpha = 0.5)),
    locations = cells_body(rows = d_er_tot > 0.05)) |> 
  tab_style(
    style = list(cell_text(weight = "bold")),
    locations = cells_summary()
  )
```

```{r total_er_annual_gg, eval=FALSE}
er_fram_rmis |> 
  drop_na(d_er_tot) |>
  group_by(sidnm) |> 
  distinct(sidnm, run_year, fram_er_tot, rmis_er_tot, d_er_tot) |> 
  ungroup() |>
  mutate(
    run_year = as.character(run_year),
    shp = case_when(
      d_er_tot > 0 ~ "FRAM > RMIS",
      d_er_tot < 0 ~ "FRAM < RMIS",
      d_er_tot == 0 ~ "FRAM = RMIS"
    ),
    sze = if_else(!is.na(d_er_tot), abs(d_er_tot), 0)
  ) |> 
  ggplot(aes(run_year, fct_rev(sidnm))) +
  geom_point(aes(shape = shp), size = 2) +
  geom_point(aes(shape = shp, size = sze, fill = d_er_tot)) +
  geom_text(aes(label = paste0(round(100*d_er_tot, digits = 0), "%")), nudge_y = 0.3) +
  scale_shape_manual(name = "", values = c("FRAM > RMIS" = 24, "FRAM < RMIS" = 25, "FRAM = RMIS" = 22, "FRAM-only" = 20, "RMIS-only" = 3)) +
  scale_size("ER diff. abs.", range = c(1, 10)) +
  scale_fill_fermenter(type = "div", palette = "PuOr", breaks = c(-0.1, -0.05, 0, 0.05, 0.1), limits = c(-.2, .2)) +
  guides(
    #size = guide_bins(range = c(2.5, 10)),
    fill = guide_colorbar(title="ER difference"),
    x = guide_axis("Run Year"),
    y = guide_axis("")) +
  labs(
    title = "Differences in Annual Total Pre-terminal Pseudo-ER", subtitle = "FRAM - RMIS"
  )

```

### Differences by Fishery

Most fisheries exhibit little or no consistent bias across stocks and/or years, with the distribution of per-year-per-stock differences with +/-0.05 [pseudo-ER](#Pseudo) units and medians near 0. 

Big differences are frequently observed in (near-)terminal fisheries where exploitation rates and consequently ER differences of the local stock(s) can be very high. As such, figure 6.5 is less illustrative of the accuracy of FRAM versus RMIS ER estimates, but more useful as a tool to identify potential estimation problems. Points to the far right of the zero line (FRAM>RMIS) can result when CWTs are missed during fishery sampling, while those to the far left of the line (RMIS>FRAM) are often caused by missing escapement of the dominant stock(s). In other cases, mapping or reporting errors may account for patterns in fisheries with similar stock compositions but roughly offset differences (e.g., 12A net vs 12B net).

**Figure 6.5.**
```{r d_er_fsh_boxes}
er_fram_rmis |> 
  filter(
    fram_lnd_yr_fsh >= 10,
    abs(d_er_fsh) > 0.001,
    abs(d_er_fsh) < 0.5
    ) |>
  mutate(run_year = as.character(run_year)) |> 
  ggplot(aes(d_er_fsh, fidnm, fill = d_er_fsh)) +
  geom_boxplot(outlier.shape = NA) +
  geom_jitter(aes(color = d_er_fsh), alpha = 0.4, position = position_jitter(height = 0.2), ) +
  scale_fill_fermenter(type = "div", palette = "PuOr", breaks = c(-0.1, -0.05, 0, 0.05, 0.1), limits = c(-.2, .2), aesthetics = c("fill", "color")) +
  geom_vline(xintercept = 0, size = 1.3) + 
  geom_vline(xintercept = c(-0.05, 0.05), linetype = 2) +
  geom_vline(xintercept = c(-0.1, 0.1), linetype = 3) +
  labs(
    title = "Fishery-specific Pseudo-ERs (FRAM-RMIS), Distributions Across Stocks and Years",
    subtitle = "Min. 10 landed Coho in FRAM, Selected Data Range: 0.001 < abs. diff < 0.5"
  )
```

The pattern of most absolute stock-year differences less than 0.1 or 0.05 [pseudo-ER](#Pseudo) units holds when looking more closely at several fisheries presumed to have highly mixed stock composition. The plot below shows the distributions of differences across stocks when stratified by year. No strong trends in the central tendency of bias were evident across years, nor did the tails of the distributions (reflecting larger differences for one or two stocks) display a clear, uniform increase or decrease between fisheries. 

<br>

**Figure 6.6. Histogram of Differences between FRAM and RMIS Pseudo-ERs**
```{r ggridges_4_fishery}
er_fram_rmis |> 
  filter(
    FisheryID %in% c(118, 91, 43, 37),
    fram_lnd_yr_fsh >= 10
    ) |>
  mutate(run_year = as.character(run_year)) |> 
  ggplot(aes(x = d_er_fsh, y = run_year)) +
  ggridges::geom_density_ridges_gradient(aes(fill = stat(x)), show.legend = F) +
  geom_rug(aes(color = d_er_fsh), alpha = 0.5) +
  geom_vline(xintercept = 0, size = 1.3) + 
  geom_vline(xintercept = c(-0.05, 0.05), linetype = 2) +
  geom_vline(xintercept = c(-0.1, 0.1), linetype = 3) + 
  scale_fill_fermenter(name = "", type = "div", palette = "PuOr", breaks = c(-0.1, -0.05, 0, 0.05, 0.1), limits = c(-.2, .2), aesthetics = c("fill", "color")) +
  scale_x_continuous("Difference in Pseudo-ER, FRAM - RMIS") +
  scale_y_discrete("") +
  facet_wrap(~fidnm, ncol = 1)

```

Even for fisheries with highly mixed stock compositions, larger differences at the whole-fishery or total scale can follow from a single stock, particularly if one or several years exert leverage on the overall distribution. Such single-year, single-stock deviations are perhaps more likely to reflect errors in FRAM (input values, input flagging) and/or RMIS (sampling, reporting) than to indicate a consistent issue in how the FRAM [base period](#BP) represents stock composition. Further, the fundamental approach in FRAM - the application of [base period](#BP) [exploitation rates](#ER) derived as averages over several fishing years - is inherently less likely to capture single-year fluctuations that may occur within the RMIS dataset. 
<br>

**Figure 6.7 Area 10 Sport Differences between FRAM and RMIS Pseudo-ERs by Stock and Year**
```{r a10_focal}
er_fram_rmis |> 
  filter(
    FisheryID %in% c(118),
    fram_lnd_yr_fsh >= 10
    ) |>
  mutate(run_year = as.character(run_year)) |> 
  drop_na(d_er_fsh) |> 
  ggplot(aes(x = sidnm, y = d_er_fsh)) +
  geom_col(aes(fill = run_year), position = position_dodge(width = 0.7, preserve = "single")) +
  coord_flip() + 
  geom_hline(yintercept = 0, size = 1.3) + 
  geom_hline(yintercept = c(-0.05, 0.05), linetype = 2) +
  geom_hline(yintercept = c(-0.1, 0.1), linetype = 3) + 
  scale_fill_brewer(palette = "Greens") +
  scale_x_discrete("") +
  scale_y_continuous("Difference in Pseudo-ER, FRAM - RMIS") +
  labs(
    title = "Fishery Specific Pseudo-ER Differences",
    subtitle = "Area 10 Sport"
  )

```

## Conclusions

We hypothesized that results of Coho [FRAM](#FRAM) [validation runs](#VR) for the fishing years 2010-2019 would correspond poorly with those calculated directly from RMIS recoveries during the same time span. Contrary to this expectation, for the stocks and fisheries with adequate data included in this analysis, we found that the 1986-1992 [FRAM](#FRAM) [base period](#BP) parameterization yielded output with mostly modest differences in [pseudo-ERs](#Pseudo) and with mostly good agreement in terms of the relative ranking of impacts among [pre-terminal](#PTerm) fisheries. In addition, most larger differences between the two datasets were plausibly explained by sampling limitations in fisheries, escapement or both. Nonetheless, determining with certainty the more accurate of the two values was difficult or impossible in almost all cases. 

This analysis may suggest a few specific stock-fishery interactions that warrant further scrutiny for possible adjustments in pre- and post-season modeling (e.g., how impacts are assigned among several Puget Sound sport fisheries). However, it clearly underscores the need for additional sampling, data reporting, and auditing of both the post-seaon Coho [FRAM](#FRAM) and RMIS datasets. At minimum, additional careful review is warranted for the [Backwards FRAM](#BkFRAM) escapement targets and catch values, the other primary input to post-season Coho [FRAM](#FRAM) runs.

The vast amount of data summarized and presented for this analysis can be the outset for further stock- and fishery- specific scrutiny and analysis by regional experts. Applying expert knowledge to some of the observed discrepancies has the potential to yield a substantial amount of data resolution.

\newpage

# Investigation of Potential Pathways to Updating the Existing Base-period 

## Introduction

The [Mixed Stock Model (MSM)](#MSM) computer program was developed to produce the [base period](#BP) [exploitation rates](#ER) that form the foundation of Coho [FRAM](#FRAM) [(Packer J., 2007)](#Pac1). Evaluations of Coho salmon [exploitation rates](#ER) in fisheries along the Canadian and U.S. West Coast rely on [FRAM](#FRAM). [CWT](#CWT) recoveries from 1986-1992, also referred to as the base-period, were used to develop the stock-fishery-time period specific exploitation rates that act as core model parameters driving [FRAM](#FRAM) calculations of the impacts associated with proposed and observed fisheries and returns. For this project, we investigated steps necessary for combining contemporary [coded-wire-tags](#CWT) with existing [base period](#BP) [coded-wire-tags](#CWT) to update [MSM](#MSM) [base period](#BP) exploitation rates for stocks with sufficient [CWT](#CWT) recoveries. <br>

For methods, we looked to the Chinook calibration that fairly recently updated the old [base period](#BP) using [CWTs](#CWT) from fishing years 1978-1984 to CWTS from 2008-2013.
For Chinook, the [FRAM](#FRAM) [base period](#BP) calibration procedures utilize CWTs outside of the [base period](#BP) years for stocks with missing CWT recoveries during the [base period](#BP) years.  However, the Chinook out-of-base procedure does not combine tags from outside of the [base period](#BP) years with tags from within the [base period](#BP) years, rather it combines stocks with recoveries from [base period](#BP) years with stocks represented by recoveries from a set of out-of-base years. To adapt the Chinook out-of-base procedure to several stocks at once for Coho and potentially combining data with current [base period](#BP) data, we identified several challenges and considerations detailed in the next paragraph.  

Although we didn't update the existing [base period](#BP) due to the difficulties and time frame associated with this project, to facilitate future development and updates, the following products were developed:<br>
 - Description of challenges and considerations when combining [base period](#BP) [CWTs](#CWT) with contemporary recoveries<br>
 - Development of [Mixed Stock Model](#MSM) documentation<br>
 - Updated instructions to map RMIS recoveries to [FRAM](#FRAM) stocks and fisheries

## Combining Contemporary [CWTS](#CWT) with the Existing [Base Period](#BP)

**1.)	Issue: Poor/No [CWT](#CWT) sampling in some contemporary fisheries**<br>
Description: To calculate [base period](#BP) [exploitation rates](#ER) for contemporary tags, modelers must develop estimates of cohort size (denominator in the [base period](#BP) exploitation rate calculation). Because cohort size is a combination of escapement and mortality (both fishing and natural mortality from later [time steps](#TS)), it is not possible to exclude fisheries from out-of-base development that are poorly sampled and still achieve correct cohort sizes.  Incorrect cohort sizes influence [base period](#BP) exploitation rate calculations for all fishery-[time steps](#TS) in the model, and therefore, even fisheries that are sampled well will require correct cohort sizes to estimate base exploitation rates.<br>
There are several potential causes of poor sampling:<br>
a.) Sampling did not occur and no [CWTs](#CWT) were recovered despite fishing. 
<br>
b.) Sampling did occur but there were sampling strata which included catch with 0 [CWT](#CWT) recoveries.  These strata would not contribute to CWT expansions to represent the full catch of the fishery. 
<br>
c.) Only a small portion of the fishery catch was sampled, resulting in a small number of [CWTs](#CWT) returned which increases uncertainty in stock compositions. 
<br>
d.) The total catch of a fishery is unknown.  Catch did occur and CWTs were recovered but were expanded to represent a [CWT](#CWT) estimated number by using a surrogate area or time period. 
<br>
Potential Solutions: The Chinook Technical Committee (CTC) creates [coded-wire-tags](#CWT) (referred to as auxiliary tags) when developing their annual [exploitation rate](#ER) analysis.  These tags are developed by regional CTC leads and represent impacts in fisheries where poor or no sampling occurred.  It may be possible to develop similar auxiliary tags for Coho to address this issue.
<br>
For freshwater areas it may be possible to develop auxiliary tags by assuming 100% of the catch comes from the local stock.  However, this assumption becomes problematic in freshwater areas that contain multiple stocks.   It would also be problematic if attempting to create auxiliaries for a marine area due to the mixed-stock nature of catch.  It is possible in mixed-stock fisheries to use [base period](#BP) [CWT](#CWT) catches to develop contemporary auxiliaries, assuming that the [base period](#BP) CWTs are representative of contemporary CWTs. However, at least two adjustments would be required.  

* An adjustment must be made to [base period](#BP) [CWT](#CWT) estimates to account for differences in fishery magnitude ([base period](#BP) vs. contemporary).  This adjustment may be achieved via post-season fishery scalers, which represent the magnitude of the contemporary fishery relative to the [base period](#BP).

* An adjustment must be made to [base period](#BP) [CWT](#CWT) estimates to account for differences in abundance ([base period](#BP) vs. contemporary).  It would be most appropriate to compare base period cohort size to contemporary cohort size, but the contemporary cohort size would be dependent upon the auxiliaries imputed. Therefore, this step represents a mathematical challenge as the cohort size is necessary to develop auxiliary CWTs, which are necessary to develop cohort size, an iterative approach is need for a solution.  
<br>

**2.)	Issue: How to aggregate annual data to produce "average" base [exploitation rate](#ER) estimates**
<br>
Description: There are two potential approaches to producing [base period](#BP) [exploitation rates](#ER).  One approach could be to calculate exploitation rates for each stock-age-fishery-[time step](#TS) annually, subsequently averaging annual rates across years to produce a [base period](#BP) exploitation rate.  An alternative approach could be to aggregate [CWT](#CWT) catch and cohort sizes across years to produce a [base period](#BP) exploitation rate.  See an example of each approach below:

```{r results ='asis'}
tibble(
  `Year` = c(2017, 2018, 2019),
  `Cohort Size` =c(500, 1000, 750),
  `CWT Catch`=c(50, 70, 150)  ) |> 
  gt::gt() |> 
  gt::cols_label(.list = list(`Year` = gt::md("**Year**"), 
                              `Cohort Size` = gt::md("**Cohort Size**"),
                              `CWT Catch`=gt::md("**CWT Catch**")
                              )) |> 
  gt::tab_header(title = "Example")
```
$$Approach 1 BPER = \frac{\frac{50}{500}+\frac{70}{1000}+\frac{150}{750}}{3}$$

$$Approach 2 BPER = \frac{50+70+150}{500+1000+750}$$

In approach 1, the BPER is calculated as 12.3%.  In approach 2, the BPER is calculated as 12%.
The advantage to approach 1 is that each year receives an equal weight in determining the [exploitation rate](#ER), irrespective of fishery catch or abundance. In approach 2, there is a greater weighting for years with a higher abundance.
The advantage to approach 2 is that many fisheries have a small number of [CWTs](#CWT) recovered.  If we use approach 1, it is possible that many zeroes will be averaged along with years with recoveries and an aggregate may better capture small sample size scenarios.<br>

**3.)	Issue: Need for additional QAQC following out-of-base procedures**<br>
Description: There were a number of anomalous results produced for out-of-base stocks relative to in-base stocks during the development of the Chinook [base period](#BP) that warranted additional investigation and QAQC (see [FRAM
documentation](https://framverse.github.io/fram_doc/)). As out-of-base expansions are applied to fisheries recoveries during the out-of-base years, additional expansions to the estimated [CWTs](#CWT) should be examined.<br>

**4.)	Issue: In the contemporary period, when fishing has not occurred in a particular area-time period, how can we develop [base period](#BP) [exploitation rates](#ER) in case that fishery does occur in the future?** <br>
Description: Contemporary fisheries are generally of a smaller magnitude than during the Coho [base period](#BP) years, with increasing spatio-temporal restrictions.  In some cases, fishery-[time step](#TS) combinations that previously had catch no longer have catches in the contemporary period.  This is problematic if there were a need to model the fishery in the future but solely use contemporary [coded-wire-tags](#CWT).
<br>
Potential Solutions: Develop methods to use old base period [coded-wire-tags](#CWT) to represent the contemporary period.  If solely using contemporary tags (e.g., not combining with current base period tags), it may be possible to use the old [base period](#BP) information as a surrogate for these fishery-time periods.  To do so, it may be necessary to divide both fishery catch and [CWT](#CWT) recoveries from the old [base period](#BP) by a large number (e.g., 10000) prior to combining with contemporary tags.  Doing so would have minimal effect on the contemporary cohort size because added impacts would be minimal.  However, through the use of a fishery scaler, it would be possible in future runs that did include the surrogate fishery-time step to model catches.  <br>

**5.)	Issue: How to model stock-specific non-retention impacts during the out-of-base procedure**<br>
Description: Expanded [CWT](#CWT) recoveries only include landed catch, as CWTs are not recovered from released fish.  Non-retention fisheries do occur for Coho coast wide and impacts from these fisheries should be accounted for during the cohort reconstruction process.  The stock composition for non-retention impacts would typically be determined during the “all stocks” run, which is not performed during an out-of-base procedure.<br>
Potential Solutions: It may be possible to determine the percentage contribution of the stock-age to the [marked](#M) retained catch in a fishery-[time step](#TS) and use that percentage to estimate the number of encounters for the stock given an estimate of marked released encounters in the fishery.  To examine this metric, biologists would need coast-wide estimates of marked released and marked retained fish by fishery-[time step](#TS).  It may be possible to alternatively calculate an average non-retention impacts per marked fish and apply this across fisheries to produce estimates more quickly at the cost of fine-scale fishery accuracy.<br>

**6.)	Issue: During processing, it is possible for [exploitation rates](#ER) to exceed 100%**<br>
Description: During the out-of-base process, a number of expansions are occurring on [CWTs](#CWT) and it is mathematically possible to produce exploitation rates that exceed 100% for a stock-age.  
Potential Solutions: If a stock-age exceeds a 100% [exploitation rate](#ER), it may be valuable to insert an error handling procedure in the program code that alerts the modeler when this occurs.  If QAQC does not find any errors in the code or input information, exploitation rates over 100% can be reset to a maximum value (e.g., 99%).<br>

**7.)	Issue: If a fishery occurred in the contemporary period, but not the [base period](#BP), how can we integrate it into the current [base period](#BP)?**
Description: If there is a contemporary fishery that was not fished during the [base period](#BP) years, no tag recoveries will exist for the fishery.  We cannot directly add the contemporary [CWTs](#CWT) from the fishery to the base years if it did not occur because doing so will affect the cohort size and all of the other fishery [base period](#BP) [exploitation rates](#ER).<br>
Potential Solutions: It may be possible to add the contemporary impacts to the old [base period](#BP) via [CWTs](#CWT) that have been scaled by dividing by a large number (e.g., 10,000).  The resulting impacts to the cohort size would be negligible and it would be possible to scale the fishery in subsequent years using a fishery scaler.  

## [Mixed Stock Model Software Documentation](#AHB2) 

The [Mixed Stock Model](#MSM) was rewritten in Visual Basic in 2007. Although the document ["Coho FRAM Base Development"](#Pac1) describes the project, the current [MSM](#MSM) lacks code documentation and a user guide, thus requiring a great deal of initiative by a new user to research the code, the data requirements, and re-(run) analyses. [Mixed Stock Model Software Documentation](#AHB2) was developed during the early stages of this project. This documentation was developed as a separate stand-alone product submitted along with this Southern Fund project.

## Mapping RMIS Recoveries to [FRAM](#FRAM) Stocks and Fisheries

For descriptions of the procedures to map RMIS [CWT](#CWT) recoveries to [FRAM](#FRAM) stocks and fisheries see chapters [3.1](#SM) and [3.2](#FM). This project developed two key mapping files. Both files are located on the [GitHub site](https://github.com/PSC-CoTC/coho_CWT_analysis/blob/main/lu_ahb_FishMap_20220826.csv) for this project. <br>
**1. Mapping RMIS [CWT](#CWT) releases to [FRAM](#FRAM) production regions, management units, and stocks.**<br>
We mapped every [marked](#M) [CWT](#CWT) release along the West Coast starting with brood year 1994 to a [FRAM](#FRAM) stocks. The mapping file "lu_cct_MSM_Stock_CWT_20220913.csv" is in a format that can be imported into the [MSM](#MSM) database and is compatible with the existing [MSM](#MSM) mapping procedures.<br>
**2. Mapping RMIS [CWT](#CWT) recoveries to [FRAM](#FRAM) fisheries.**<br>
Fishery mapping instructions in the file "lu_ahb_FishMap_20220826.csv" map the vast majority of RMIS recoveries since 1986 to [FRAM](#FRAM) fisheries. Fisheries not mapped do not correspond to a [FRAM](#FRAM) fishery, such as certain High Seas or test fishery recoveries. The mapping instructions (see chapter [3.2](#FM)) rely on the same logic needed by the [MSM](#MSM) program to make valid [FRAM](#FRAM) fishery assignments.



\newpage

# Appendix{#A}

## List of Coho [FRAM](#FRAM) Fisheries {#F}

```{r echo = FALSE, results='asis'}
readRDS(paste0(table_dir,"framvs_doc_lu_fishery_coho.rds"))
```

## List of Coho [FRAM](#FRAM) Stocks {#S}

```{r echo = FALSE, results='asis'}
readRDS(paste0(table_dir,"framvs_doc_lu_stock_coho.rds"))
```

## Table of [CWT](#CWT) Adjustments (Expansions), [FRAM](#FRAM) Catches, RMIS Catches, and Sampled Catches by Fishery and Year {#CWTAdj}

```{r}
# see fram_samp_filla.xlsx 

#compare FRAM annual landed to RMIS
rmis_landed<-rmis_landed|>
  filter(!is.na(Total_rmis))
fram_samp_full<-rmis_landed[,c(1:2,4:5)]|>
  full_join(fram_landed,
            by=c("RunYear","FisheryID"))|>
  left_join(d_uniquefishery,
            by = c("FisheryID"))|>
  filter(FisheryID!=0)|>
  filter(FisheryID!=300)|>
  filter(RunYear!=2020)|>
  left_join(d_rmis_no_samp|>
              select(-number_sampled),
            by=c("RunYear","FisheryID","FisheryName"))
 

action <- na.omit(read.csv(paste0(table_dir,"/MR_Action.csv")))

#ColRiver fisheries are included in FRAM escapement and need to be reassigned to RMIS escapement for comparison to FRAM
ColR<-c(24,25,26,31,32)

fram_samp_fulla <- fram_samp_full |>
  select(-c(Kept_M_fram, Sampled_rmis)) |>
  ungroup() |>
  #mutate(Xpansion = NA) |> 
  pivot_longer(-c(FisheryID, FisheryName, RunYear), names_to = "Source", values_to = "Landed") |>
  group_by(FisheryID, FisheryName, Source) |>
  pivot_wider(names_from = RunYear, values_from = Landed) |>
  mutate(Source = if_else(Source == "Total_fram", "FRAM",
                          if_else(Source == "Total_rmis", "RMIS",
                                  if_else(Source == "number_caught", "Sampled", Source)))) |>
  arrange(FisheryID, Source) |>
  ungroup() |>
  rowwise(FisheryID, FisheryName, Source) |>
  mutate(Total = sum(c_across(`2010`:`2019`),
                     na.rm = TRUE)) |>
  ungroup() |>
  left_join(action, 
            select(FisheryID,Action),
            by = c("FisheryID", "FisheryName")) |> 
  bind_rows(fram_samp_full |>
              select(-c(Kept_M_fram, Sampled_rmis)) |>
              ungroup() |>
              #mutate(Xpansion = NA) |> 
              pivot_longer(-c(FisheryID, FisheryName, RunYear), names_to = "Source", values_to = "Landed") |>
              group_by(FisheryID, FisheryName, Source) |>
              pivot_wider(names_from = RunYear, values_from = Landed) |>
              mutate(Source = if_else(Source == "Total_fram", "FRAM",
                                      if_else(Source == "Total_rmis", "RMIS",
                                              if_else(Source == "number_caught", "Sampled", Source)))) |>
              arrange(FisheryID, Source) |>
              ungroup() |>
              left_join(action, 
                        select(FisheryID,Action),
                        by = c("FisheryID", "FisheryName")) |> 
              pivot_longer(`2010`:`2019`, names_to = "Year", values_to = "Value") |> 
              pivot_wider(names_from = "Source", values_from = "Value") |> 
              mutate(Expansion = case_when(
                        Action == "Remove  fishery" ~ 0,
                        FisheryID == 198 ~ 1,
                        FisheryID == 98 ~ 0,
                        FisheryID == 91 ~ FRAM / RMIS,
                       (abs(FRAM - RMIS) > 400 & abs(FRAM - RMIS) / RMIS > 0.095)|FRAM<=10|is.na(FRAM)|RMIS<=10|is.na(RMIS)|is.na(Sampled)|Sampled<=10|RMIS/Sampled>8 ~ 0,
                        is.na(RMIS) ~ 0,
                        is.na(Sampled) ~ 0,
                        TRUE ~ RMIS / Sampled
                    )
              ) |> 
              mutate(Expansion = if_else(FisheryID %in% ColR,RMIS/Sampled,Expansion))|>
              mutate(Expansion = if_else(FisheryID ==31 ,3,Expansion))|> # approximate sampling expansion when summing over all years
              mutate(Expansion = if_else(FisheryID %in% ColR & Expansion >8,8,Expansion))|># this will prevent the expansion of 10 for fishery 24 in 2010 from being excluded from Dan's code
              group_by(Year) |> 
              mutate(Expansion = if_else(FisheryID == 97, 
                                         (RMIS[FisheryID == 98] + RMIS) / (Sampled[FisheryID == 98] + Sampled),
                                         Expansion),
                     Expansion = if_else(is.na(Expansion), 
                                         RMIS / Sampled,
                                         Expansion)) |> #filter(FisheryID %in% 97:98) |> 
              ungroup() |> 
              select(-c(FRAM, RMIS, Sampled)) |> 
              pivot_wider(names_from = Year, values_from = Expansion) |> 
              mutate(Source = "Xpansion")
  ) |> 
  
  arrange(FisheryID, Source)
#write.csv(fram_samp_fulla, "fram_samp_fulla.csv", na = "")
```


```{r tab_foo,   results='asis' } 
table <- fram_samp_fulla
table$Total<-round(table$Total,0)

reactable(table |> mutate(Source = if_else(Source == "Xpansion", "Expansion", Source)),
          defaultColDef = colDef(align = "center",
                                 headerStyle = list(background = "#f7f7f8")),
          rowStyle = function(index){
            if(index %% 4 == 0){
              list(`border-bottom` = "solid")
            } else{
              FALSE
            }
          },
          columns = list(FisheryID = colDef(name = "Fishery ID",
                                            filterable = TRUE,
                                            style = function(value, index){
                                               if(index %% 4 == 0){
                                                 list(background = "#ff8c00", fontWeight = "normal")
                                               }
                                              }),
                         FisheryName = colDef(name = "Fishery Name",
                                              minWidth = 140,
                                              filterable = TRUE,
                                              style = function(value, index){
                                                 if(index %% 4 == 0){
                                                   list(background = "#ff8c00", fontWeight = "normal")
                                                 }
                                                }),
                         Source = colDef(filterable = TRUE,
                                            style = function(value, index){
                                               if(index %% 4 == 0){
                                                 list(background = "#ff8c00", fontWeight = "normal")
                                               }
                                              }),
                         `2010` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if(is.na(value) | table$`2010`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2011` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if(is.na(value) | table$`2011`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2012` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2012`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2013` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2013`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2014` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2014`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2015` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2015`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2016` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2016`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2017` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2017`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2018` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2018`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         `2019` = colDef(na = "-",
                                         cell = function(value, index){
                                             if(is.na(value)){
                                               value <- "-"
                                             } else if(value == 0){
                                               value
                                             } else if(index %% 4 == 0) {
                                               value <- format(round(value, 2), nsmall = 2)
                                             } else {
                                               value
                                             }
                                         },
                                         style = function(value, index){
                                           if(index %% 4 == 1){
                                             fram <- index
                                             rmis <- index + 1
                                             sampled <- index + 2
                                             xpansion <- index + 3
                                           } else if(index %% 4 == 2){
                                             fram <- index - 1
                                             rmis <- index
                                             sampled <- index + 1
                                             xpansion <- index + 2
                                           } else if(index %% 4 == 3){
                                             fram <- index - 2
                                             rmis <- index - 1
                                             sampled <- index
                                             xpansion <- index + 1
                                           } else{
                                             fram <- index - 3
                                             rmis <- index - 2
                                             sampled <- index  - 1
                                             xpansion <- index
                                           }
                                           if(index %% 4 == 0){
                                             list(background = "#ff8c00", fontWeight = "bold")
                                           } else if (is.na(value) | table$`2019`[xpansion] == 0){
                                             list(background = "#a9a9a9")
                                           } else{
                                             FALSE
                                           }
                                         }),
                         Total = colDef(na = "-",
                                        style = function(value, index){
                                               if(index %% 4 == 0){
                                                 list(background = "#ff8c00", fontWeight = "normal")
                                               }
                                              }),
                         Action = colDef(filterable = TRUE,
                                         align = "left",
                                         minWidth = 270,
                                         style = function(value, index){
                                               if(index %% 4 == 0){
                                                 list(background = "#ff8c00", fontWeight = "normal")
                                               }
                                              })
                     ),
          bordered = FALSE,
          outlined = TRUE,
          highlight = TRUE,
          showPageSizeOptions = TRUE,
          pageSizeOptions = c(4, 12, 24, 58),
          defaultPageSize = 12) #%>%
          #reactablefmtr::add_title("FRAM Catches, RMIS Catches, Catch Associated with Sampling and Expansions by Fishery and Year", font_size=16)


# write.csv(fram_samp_full,"fram_samp_full.csv")
```


## Coded-wire Tag Groups which were Listed as Clipped under “cwt_1st_mark”, but were Released as an Unclipped Group {#WrongClip}

```{r results='asis' }
t_unclipped|> gt::gt()
```

## Coho Salmon Production Regions and Management Units {#CoPRMU}

```{r results='asis'}
t_prmu|> gt::gt()
```

## Coho Stock Mapping by Tag Code {#StkMap}

```{r results='asis'}
stockmaptable<-as.data.frame(d_stockmap)|>
  filter(brood_year>2009)|> 
  DT::datatable(filter = 'top',
                extensions = 'Buttons',
                options = list(
                  dom = 'Bfrtip',
                  buttons = c('copy', 'csv', 'excel', 'pdf', 'print')))
stockmaptable
```

## Number of [CWT](#CWT) Recoveries 2010-2019 per Assessed, Mapped [FRAM](#FRAM) Stock {#TagRec}

```{r results='asis'}
gt::gt(rmis_n_recs_per_stock) |> 
  gt::sub_missing()
```


## Run Information {#RunInfo}

[FRAM](#FRAM) Access Database Name: 'FramVS2-PSC-Coho-Backwards-86 to 2016 BEST.mdb'

[FRAM](#FRAM) Version: 2.20March25

Base Period ID: 27

Description of Model Runs for the Sensitivity Analysis

```{r echo = FALSE, results='asis'}
tablex_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_x_section1c.xlsx"))

reactable::reactable(tablex_1c |> 
                     dplyr::rename(`Run Name` = RunName),
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       align = "center",
                                                       headerStyle = list(background = "#f7f7f8")),
                     highlight = TRUE,
                     outlined = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     columns = list(
                      `Run ID` = reactable::colDef(minWidth = 80),
                      `Run Name` = reactable::colDef(minWidth = 225),
                      Description = reactable::colDef(minWidth = 200),
                      `Run Year` = reactable::colDef(minWidth = 90)
                     ),
                     rowStyle = function(index){
                       if(tablex_1c$`Run Year`[index] == 2013){
                         list(background = "#ADD8E680")
                       } else if(tablex_1c$`Run Year`[index] == 2015){
                         list(background = "#90ee9080")
                       }
                     }
                )
```

## List of Fisheries Set to Zero for the Sensitivity Analysis {#FZero}

```{r echo = FALSE, results = 'asis'}
# **ID**|**Fishery Title**|**ID**|**Fishery Title**|**ID**|**Fishery Title**
# :---:|:---|:---:|:---|:---:|:---
# 48|Grays Harbor Sport (2.2)|75|Hoh R C&S|126|Green/Duwamish R Net
# 49|South Grays Harbor Spt|82|6D NT Net (Dungeness Bay & R)|127|Green/Duwamish R Sport
# 50|Grays Harbor Est. Net|83|6D T Net (Dungeness Bay & R)|128|Lk Wash/Samm./Tribs Spt
# 51|Humptulips R Sport|84|Elwha R Net|134|Puyallup R Net
# 52|Lower Chehalis R Net|85|West JDF Straits Trib Net|135|Puyallup R Sport
# 53|Humptulips R C&S|86|East JDF Straits Trib Net|139|A 13C NT Net (Chambers)
# 54|Chehalis R Sport|89|East JDF Straits Trib Sport|140|A 13C TNet (Chambers)
# 55|Humptulips R Net|90|West JDF Straits Trib Sport|141|A 13A NT Net (Carr Inlet)
# 56|Upper Chehalis R Net|94|Dungeness R Sport|142|A 13A T Net (Carr Inlet)
# 57|Chehalis R C&S|95|Elwha R Sport|147|Nisqually R Net
# 58|Wynochee R Sport|98|Nooksack R Net|148|McAllister Creek Net
# 59|Hoquiam R Sport|99|Nooksack R Sport|150|Nisqually R Sport
# 60|Wishkah R Sport|100|Samish R Sport|151|Deschutes R Sport
# 61|Satsop R Sport|103|Skagit R Net|155|A 9/9A NT Net
# 62|Quinault R Sport|104|Skagit River Test Net|156|A 9/9A T Net (On Res)
# 63|Quinault R Net|105|Swinomish Channel Net|157|A 12A NT Net (Quilcene Bay)
# 64|Quinault R C&S|108|Skagit R Sport|158|A 12A T Net (Quilcene Bay)
# 65|Queets R Sport|111|WA A 8D NT Net (Tulalip Bay)|161|Skokomish R Net
# 66|Clearwater R Sport|112|WA A 8D T Net (Tulalip Bay)|162|Quilcene R Net
# 67|Salmon R Sport (Queets)|113|Stillaguamish R Net|163|12, 12B Trib FW Sport
# 68|Queets R Net|114|Snohomish R Net|164|12A Trib FW Sport
# 69|Queets R C&S|116|Stillaguamish R Sport|165|12C, 12D Trib FW Sport
# 70|Quillayute R Sport|117|Snohomish R Sport|166|Skokomish R Sport
# 71|Quillayute R Net|121|WA A 10A NT Net (Elliott Bay)|167|Lower Fraser R Term Catch
# 72|Quillayute R C&S|122|WA A 10A T Net (Elliott Bay)|168|Upper Fraser R Term Catch
# 73|Hoh R Sport|123|WA A 10E NT Net (East Kitsap)|169|Lower Fraser River Sport
# 74|Hoh R Net|124|WA A 10E T Net (East Kitsap)|184|Fraser R Gill Net
# ---|---|125|WA A 10F-G T Net (Lake Union)|193|Alberni Canal Sport

table4_1c <- readxl::read_excel(paste0(table_dir, "SA_Table_4_section1c.xlsx"))

reactable::reactable(table4_1c,
                     outlined = TRUE,
                     striped = TRUE,
                     height = 500,
                     showPageSizeOptions = TRUE,
                     defaultPageSize = 20,
                     pageSizeOptions = c(10, 20, 50, 100),
                     highlight = TRUE,
                     defaultColDef = reactable::colDef(vAlign = "center", 
                                                       headerVAlign = "bottom", 
                                                       #align = "center", 
                                                       headerStyle = list(background = "#f7f7f8")),
                     columns = list(
                       ID = reactable::colDef(minWidth = 20, align = "center")
                     )
                    )
```

## Treaty/Non-treaty Net Fishery Pairs {#Pairs}

```{r table, echo= FALSE, results='asis'}
dplyr::tibble(`NT FishID` = c(80,82,87,96,101,109,111,119,121,123,130,132,137,139,141,143,145,153,155,157,159),
  `T FishID` =c(81,83,88,97,102,110,112,120,122,124,131,133,138,140,142,144,146,154,156,158,160),
  `Fishery Name` = c("A4B6CNet","Ar6D Net","A6-7ANet","A7BCDNet","Ar 8 Net","Ar8A Net","Ar8D Net","Ar10 Net","Ar10ANet","Ar10ENet","Ar11 Net","Ar11ANet","Ar13 Net","Ar13CNet","Ar13ANet","Ar13DNet","A13FKNet","1212BNet","A9-9ANet","Ar12ANet","A12CDNet")
       ) |> 
 knitr::kable(caption = "Treaty/non-treaty net fishery pairs") |>
  kableExtra::kable_styling() |>
  kableExtra::scroll_box( width="100%",height = "250px")

```

## [FRAM](#FRAM) Time Steps {#TStep}

```{r appendi2, echo= FALSE, results='asis'}
dplyr::tibble(
  `Time Step` = c(1,2,3,4,5),
  Duration =c("Jan-June","July","August","September","Oct-Dec")) |> 
  gt::gt() |> 
  gt::cols_label(.list = list(`Time Step` = gt::md("**Time Step**"), 
                              Duration = gt::md("**Duration**"))) |> 
  gt::tab_header(title = "Coho FRAM time steps")
```

## [Mark Rates](#MR) from [FRAM](#FRAM) and Sampling by Fishery and Time Step {#markrates}

```{r tableMR, echo= FALSE, results='asis'}

  MR_Table<-mr_fram_samp_rmis[,c(1,2,3:7, 9:12, 14)]|>
  left_join(d_uniquefishery %>% select(FisheryID,FisheryName), by = c("FisheryID"))
  
    MR_Table |> 
      mutate(TimeStep = factor(TimeStep),
             RunYear = factor(RunYear),
             FisheryID = factor(FisheryID),
             FisheryFlag = factor(FisheryFlag),
             across(c(5:6, 9:10), ~round(.x)),
             across(c(7, 11), ~round(.x, 2))) |>
      dplyr::rename(`Run Year` = RunYear,
             `Fish ID` = FisheryID, 
             `Fish Name` = FisheryName,
             `Time Step` = TimeStep,
             `Flag` = FisheryFlag,
             `FRAM Total` = Total,
             `FRAM Marked` = Marked,
             `FRAM Mark Rate` = MR,
             `FRAM Type` = Type.x,
             `Sampl total` = Tot_samp,
             `Sampl Marked` = Marked_samp,
             `Sampl Mark Rate` = MR_samp,
             `Sampl Type` = Type.y)|>
      DT::datatable(filter = 'top', 
                    extensions = 'Buttons', 
                    options = list(
                      dom = 'Bfrtip',
                      buttons = c('copy', 'csv', 'excel', 'pdf', 'print')))
    
```

## List of Model Runs, Fisheries, and [Time Steps](#TS) with Catch Reductions to Avoid Negative Escapements {#NegEsc}

**Model Run ID**|**ID**|**Fishery**|**Time Step**|**New Catch**
:---:|:---:|:---:|:---:|:---:
43, 44|160|A 12C-12D Treaty Net|5|2,000
69, 70|107|A 9 Sport (Admirality Inlet)|4|10,000
||118|A 10 Sport (Seattle)|4|5,000
91, 92|96|A 7B-7C-7D Non-Treaty Net|4|10,000
||97|A 7B-7C-7D Treaty Net|4|10,000
93/94, 97/98, 99/100|144|A 13D Treaty Net|4|2,000
||144|A 13D Treaty Net|5|2,000



\newpage

# Bibliography

<a id='Auer1'>Auerbach, Dan, Hagen-Breaux A., Dapp D., Bellman M., Miler O.. 2020. FRAM_VS User Guide. Portland:PFMC.

<a id='Pac1'>Cope, Robert. 2011. Causes and Effects of Bias in Anticipated Mar Rates in Mark-Selective Fisheries for Coho Salmon. Salmon Methodology Review Report, Portland: Pacific Fishery Management Council.

<a id='AHB1'>Hagen-Breaux, Angelika. 2018. BkFRAM May 2 2018.docx. Olympia: Washington Department of Fish and Wildlife.

<a id='AHB2'>Hagen-Breaux, Angelika. 2022. Mixed Stock Model Software Documentation (MSM_VS 1.2). Vancouver: Pacific Fishery Management Council.

<a id='Lapi1'>Lapi, L. et al. 1989. Information Content and Data Standards of a Coastwide Coded-Wire Tag Database. Report TCDS (89)-1. Vancouver: Pacific Salmon Commission.

<a id='MEW1'>Model Evaluation Workgroup. 2008. Fishery Regulation Assessment Model (FRAM) - An Overview for Coho and Chinook v 3.0. Portland: Pacific Fishery Management Council.

<a id='PSC1'>Pacific Salmon Commission. 2014. Pacific Salmon Treaty. Vancouver: Pacific Salmon Commission.

<a id='PSC2'>Pacific Salmon Commission. 2020. Pacific Salmon Treaty. Vancouver: Pacific Salmon Commission.

<a id='Pac1'></a>Packer, James and Cook-Tabor C. 2007. Coho FRAM Base Period Development. Vancouver: Pacific Salmon Commission. 

<a id='SFEC1'>Selective Fishery Evaluation Technical Committee (2021). ANALYSIS OF COHO SALMON DOUBLE INDEX TAG (DIT) GROUPS 
FOR BROOD YEARS 1998:2011. Vancouver (BC): Pacific Salmon Commission.

<a id='SMAWG1'>Salmon modeling and analysis workgroup. 2022. FRAM Overview. https://framverse.github.io/fram_doc/ built September 29, 2022.