diff --git a/Rscripts/logbook_analys.R b/Rscripts/logbook_analys.R new file mode 100644 index 0000000..b394aa5 --- /dev/null +++ b/Rscripts/logbook_analys.R @@ -0,0 +1,87 @@ +library(dplyr) + +observer <- read.csv('data/confidential/wcgop_logbook/ward_oken_observer_data_2024-09-10.csv') |> + as_tibble() |> + filter(AVG_LAT > 40+1/6, + sector == 'Midwater Rockfish', + DATATYPE == 'Analysis Data') + +observer_yt <- observer|> + group_by(HAUL_ID) |> + summarise(has_yt = 'Yellowtail Rockfish' %in% species, + across(AVG_DEPTH:YEAR, first)) |> + filter(!has_yt) |> + mutate(species = 'Yellowtail Rockfish', MT = 0) |> + select(-has_yt) |> + bind_rows( + filter(observer, species == 'Yellowtail Rockfish') + ) |> + mutate(DRVID = as.integer(DRVID)) + +em <- read.csv('data/confidential/wcgop_logbook/ward_oken_em_data_2024-12-03.csv') |> + as_tibble() |> + filter(AVG_LAT > 40+1/6, sector == 'Midwater Rockfish EM') + +em_yt <- em |> + group_by(HAUL_ID) |> + summarise(has_yt = 'Yellowtail Rockfish' %in% species, + across(AVG_DEPTH:YEAR, first)) |> + filter(!has_yt) |> + mutate(species = 'Yellowtail Rockfish', MT = 0) |> + select(-has_yt) |> + bind_rows( + filter(em, species == 'Yellowtail Rockfish') + ) + +all_hauls_yt <- bind_rows(list(obs = observer_yt, em = em_yt), .id = 'program') |> + filter(HAUL_DURATION > 0, HAUL_DURATION < 10) |> + mutate(depth_std = (AVG_DEPTH - mean(AVG_DEPTH))/sd(AVG_DEPTH), + fyear = factor(YEAR), log_duration = log(HAUL_DURATION), + cpue = MT/HAUL_DURATION, month = lubridate::month(as.Date(SET_DATE)), + fDRVID = factor(DRVID), fport = factor(R_PORT), + program = factor(program)) # note D_PORT is NA in some cases, R_PORT is always defined (fish ticket) + +nrow(filter(all_hauls_yt, MT>0)) / nrow(all_hauls_yt) +# good to monitor % positive hauls. This is pretty high (70%). + +yt_mesh <- sdmTMB::make_mesh(all_hauls_yt, xy_cols = c('AVG_LONG', 'AVG_LAT'), n_knots = 50) +plot(yt_mesh) +# this does not look like a st model is a good idea. very clear "hotspots" of fishing effort. +# but we need it for sdmTMB to work (not actually used) + +mod_simple <- sdmTMB::sdmTMB(formula = cpue ~ fyear + s(depth_std) + s(month, bs = 'cc') + program + (1|fDRVID) + (1|fport) - 1, family = sdmTMB::delta_lognormal(), + data = all_hauls_yt, spatial = 'off', spatiotemporal = 'off', time = 'fyear', mesh = yt_mesh) + +summary(mod_simple) + + +resids <- residuals(mod_simple, model = 2) # the presence/absence model looks consistently good +qqnorm(resids) +qqline(resids) + +sdmTMB::visreg_delta(mod_simple, xvar = "depth_std", model = 1, gg = TRUE) +sdmTMB::visreg_delta(mod_simple, xvar = "depth_std", model = 2, gg = TRUE) +sdmTMB::visreg_delta(mod_simple, xvar = "month", model = 1, gg = TRUE) +sdmTMB::visreg_delta(mod_simple, xvar = "month", model = 2, gg = TRUE) + +preds <- predict(mod_simple, + newdata = tibble(fyear = factor(2012:2023), depth_std = 0, month = 6, program = 'obs', fDRVID = '546053', fport = 'NEWPORT'), + return_tmb_object = TRUE, re_form = NA, re_form_iid = NA) + +ind <- sdmTMB::get_index(preds, bias_correct = TRUE) + +ind |> + ggplot(aes(x = as.numeric(fyear), y = est, ymin = lwr, ymax = upr)) + + geom_line() + + geom_errorbar() + + +all_hauls_yt |> + ggplot() + + geom_point(aes(x = HAUL_DURATION, y = MT)) + geom_density(aes(x = HAUL_DURATION, col = YEAR, group = YEAR)) + +### general observation: as I have added terms to the model, the signal in the index went from well-defined to noisy and meaningless. +### the port (8 levels) and vessel (46 levels) effects really did it in. should think more carefully about what terms to include. +### note 3962 hauls, 2762 positive hauls over 12 years. ~330 hauls/yr, 230 positive hauls/yr. # of hauls varies a lot by year. +### (increases 2012-2018) \ No newline at end of file diff --git a/docs/data_summary_doc.qmd b/docs/data_summary_doc.qmd index 10bea23..17e95e5 100644 --- a/docs/data_summary_doc.qmd +++ b/docs/data_summary_doc.qmd @@ -25,12 +25,12 @@ This documents data streams for use in the Northern Yellowtail Rockfish stock as Some summary figures: ```{r pacfin-catch, cache=TRUE} -load(here('data/confidential/commercial/pacfin_12-09-2024/PacFIN.YTRK.CompFT.09.Dec.2024.RData')) +load(here('data/confidential/commercial/pacfin_12-11-2024/PacFIN.YTRK.CompFT.11.Dec.2024.RData')) # I have checked and LANDING_YEAR = PACFIN_YEAR for all rows. yt_n_catch_pacfin <- catch.pacfin |> filter( - AGENCY_CODE == 'W' | AGENCY_CODE == 'O' | COUNTY_NAME == 'DEL NORTE' | COUNTY_NAME == 'HUMBOLDT', + AGENCY_CODE == 'W' | AGENCY_CODE == 'O' | COUNTY_NAME == 'DEL NORTE' | COUNTY_NAME == 'HUMBOLDT', # OR, WA, N. CA LANDING_YEAR < 2024 ) @@ -38,24 +38,29 @@ yt_n_catch_pacfin <- catch.pacfin |> or_comm_catch <- read.csv(here("Data/Confidential/commercial/Oregon Commercial landings_433_2023.csv")) |> tibble::as_tibble() |> group_by(YEAR) |> - summarise(catch_mt = sum(TOTAL)) + summarise(catch_mt = sum(TOTAL)) |> + mutate(AGENCY_CODE = 'O') -# WA URCK/UPOP reconstruction +# WA URCK/UPOP reconstruction, use prior to 2000 wa_urck <- read.csv(here('Data/Raw_not_confidential/WA_URCK_UPOP_recon_1981-2016.csv')) |> - filter(Species_AfterComps %in% c('YTRK', 'YTR1')) |> + filter(Species_AfterComps %in% c('YTRK', 'YTR1'), + Year < 2001) |> group_by(Year) |> summarise(catch_mt = sum(RoundPounds_AfterComps) * 0.453592 / 1000) |> - mutate(AGENCY_CODE = 'W', - model = 'URCK/UPOP') + mutate(AGENCY_CODE = 'W') |> + rename(YEAR = Year) # Put it all together comm_draft_landings <- yt_n_catch_pacfin |> - filter(AGENCY_CODE != 'O') |> + filter(AGENCY_CODE != 'O', + AGENCY_CODE != 'W' | PACFIN_YEAR >= 2001 + ) |> group_by(PACFIN_YEAR, AGENCY_CODE)|> summarise(catch_mt = sum(ROUND_WEIGHT_MTONS)) |> rename(YEAR = PACFIN_YEAR) |> bind_rows( - mutate(or_comm_catch, AGENCY_CODE = 'O') + or_comm_catch, + wa_urck )|> group_by(YEAR, AGENCY_CODE) |> summarise(catch_mt = sum(catch_mt)) |> @@ -109,23 +114,7 @@ comm_draft_landings |> geom_line(aes(x = Year, y = catch_mt, col = AGENCY_CODE, linetype = model), linewidth = 1) + viridis::scale_color_viridis(discrete = TRUE) + labs(title = 'Comparing commercial landings in 2017 model vs current', x = 'YEAR') - -# compare 2017 landings to current landings, WA only -comm_draft_landings |> - rename(Year = YEAR) |> - mutate(model = 'PacFIN') |> - filter(AGENCY_CODE == 'W') |> - bind_rows( - filter(comm2017, AGENCY_CODE == 'W') - )|> - bind_rows(wa_urck) |> - ggplot() + - geom_line(aes(x = Year, y = catch_mt, col = model), linewidth = 1) + - viridis::scale_color_viridis(discrete = TRUE) + - labs(title = 'Possible sources for WA landings (vert. line = 1994)', x = 'YEAR') + - geom_vline(xintercept = 1994) ``` -URCK/UPOP = sum of RoundPounds_AfterComps for Species_AfterComps = YTR1 or YTRK, converted to mt. Table of commercial landings by state, 1981 onward. @@ -139,9 +128,11 @@ comm_draft_landings |> knitr::kable(align = 'l', digits = 1) ``` -WDFW has alerted the STAT these do not include all tribal catches in recent years. +WDFW has alerted the STAT these do not include all tribal catches in recent years. 1981-2000 are from the UPOP/URCK reconstruction, summing catches for YTRK and YTR1. 2001 onward comes from PacFIN -Note this is for catch landed into Del Norte and Humboldt counties only. +Oregon landings were provided directly to the STAT from ODFW and come from PacFIN and a UPOP/URCK reconstruction. + +California landings are for catch landed into Del Norte and Humboldt counties only. ### At-Sea landings @@ -159,14 +150,15 @@ ashop_catch |> ```{r rec-catch-modern} # WASHINGTON -wa_modern <- read.csv(here('Data/Raw_not_confidential/RecFIN_WA_catch_to_2023.csv')) |> - filter(RECFIN_WATER_AREA_NAME != 'Canada', RECFIN_YEAR <= 2023) |> +wa_modern <- read.csv(here('Data/Confidential/Rec/CTE501-WASHINGTON-1990---2023.csv')) |> + filter(RECFIN_WATER_AREA_NAME != 'CANADA', RECFIN_YEAR <= 2023, + RECFIN_CATCH_AREA_NAME != 'SEKIU AND PILLAR POINT') |> tibble::as_tibble() |> group_by(RECFIN_YEAR) |> - summarise(Dead_Catch = sum(SUM_TOTAL_MORTALITY_MT)) |> + summarise(Dead_Catch = sum(TOTAL_MORTALITY_MT)) |> mutate(State = 'WA (mt)') -wa_historical <- read.csv(here('Data/Raw_not_confidential/WA_historical_rec.csv')) |> +wa_historical <- read.csv(here('Data/Raw_not_confidential/WA_historical_rec.csv')) |> tibble::as_tibble() |> filter(AREA %in% 1:4) |> # these are marine areas, per T. Tsou on 12/10/24 in meeting notes group_by(RECFIN_YEAR) |> @@ -247,7 +239,8 @@ ca_mrfss_catch <- read.csv(here('Data/Confidential/Rec/RecFIN_CA_MRFSS.csv')) |> n_ca_catch_mt * weighted_ratio), State = 'CA (mt)') |> select(-n_ca_catch_mt, -weighted_ratio) |> - rename(RECFIN_YEAR = YEAR) + rename(RECFIN_YEAR = YEAR) |> + filter(RECFIN_YEAR != 1980) # Interpolate between 3-year average before and after MRFSS break ca_mrfss_interpolate <- tibble( @@ -298,7 +291,7 @@ For CA rec catch estimates for the assessment: Initial length sample sizes after running `PacFIN.Utilities::cleanPacFIN()`: ```{r comm-bio, cache=TRUE} -load(here('Data/Confidential/Commercial/pacfin_12-09-2024/PacFIN.YTRK.bds.09.Dec.2024.RData')) +load(here('Data/Confidential/Commercial/pacfin_12-11-2024/PacFIN.YTRK.bds.11.Dec.2024.RData')) bds_clean <- PacFIN.Utilities::cleanPacFIN(bds.pacfin, verbose = FALSE) bds_clean |> filter(state == 'WA' | state == 'OR' | PACFIN_GROUP_PORT_CODE == 'CCA' | diff --git a/docs/data_summary_doc_files/figure-html/ashop-bio-fig-1.png b/docs/data_summary_doc_files/figure-html/ashop-bio-fig-1.png index 9c3550e..ddfcd34 100644 Binary files a/docs/data_summary_doc_files/figure-html/ashop-bio-fig-1.png and b/docs/data_summary_doc_files/figure-html/ashop-bio-fig-1.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-catch-1.png b/docs/data_summary_doc_files/figure-html/pacfin-catch-1.png index 41b89ff..dd5bd8f 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-catch-1.png and b/docs/data_summary_doc_files/figure-html/pacfin-catch-1.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-catch-2.png b/docs/data_summary_doc_files/figure-html/pacfin-catch-2.png index c6835ad..e847795 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-catch-2.png and b/docs/data_summary_doc_files/figure-html/pacfin-catch-2.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-catch-3.png b/docs/data_summary_doc_files/figure-html/pacfin-catch-3.png index a994b63..e37dcd8 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-catch-3.png and b/docs/data_summary_doc_files/figure-html/pacfin-catch-3.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-catch-4.png b/docs/data_summary_doc_files/figure-html/pacfin-catch-4.png index 0c11745..1cc3916 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-catch-4.png and b/docs/data_summary_doc_files/figure-html/pacfin-catch-4.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-catch-5.png b/docs/data_summary_doc_files/figure-html/pacfin-catch-5.png index 3a245ce..ad7a574 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-catch-5.png and b/docs/data_summary_doc_files/figure-html/pacfin-catch-5.png differ diff --git a/docs/data_summary_doc_files/figure-html/pacfin-midwater-vs-bottom-figs-1.png b/docs/data_summary_doc_files/figure-html/pacfin-midwater-vs-bottom-figs-1.png index 27b56af..fc5562f 100644 Binary files a/docs/data_summary_doc_files/figure-html/pacfin-midwater-vs-bottom-figs-1.png and b/docs/data_summary_doc_files/figure-html/pacfin-midwater-vs-bottom-figs-1.png differ diff --git a/docs/data_summary_doc_files/figure-html/recfin-bio-1.png b/docs/data_summary_doc_files/figure-html/recfin-bio-1.png index aeaff32..d6fe25d 100644 Binary files a/docs/data_summary_doc_files/figure-html/recfin-bio-1.png and b/docs/data_summary_doc_files/figure-html/recfin-bio-1.png differ diff --git a/docs/index.html b/docs/index.html index 735e7c7..17c1b54 100644 --- a/docs/index.html +++ b/docs/index.html @@ -6,7 +6,7 @@ - +