which_species_occur.Rmd

---
title: "Which species occur in a given area over a given time?"
author: "Andrew Edwards"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Which species occur in a given area over a given time?}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 10,
  fig.height = 8
)
```

This vignette shows how to find out which species were caught in a given area
during a given time.

For this example, we are specifically asking what species were seen in the Strait of
Georgia in 2018, the only year in which the IPHC survey went into the
Strait. Dick Beamish was asking, based on a hallway discussion originally
motivated by divers filming Bluntnose Sixgill Shark in Port Alberni.

This vignette requires the data for each species to already be saved
locally, by running the 'data for all species' vignette.

```{r setup, results = FALSE}
load_all()
# library(gfiphc)
library(dplyr)

path_data = "all-species-data"          # folder where have already saved
                                        # extracted raw data
path_results = "all-species-results"    # folder where have already saved
                                        # results (think just time series though)
```

## Define species of interest (all in this case)

Want to check all species used in the groundfish synopsis report:
```{r species}
sp_vec <- gfsynopsis::get_spp_names() %>%
  dplyr::arrange(species_code) %>%
  dplyr::pull(species_common_name)
# sp_vec
```

## Select only stations in Strait of Georgia

Easiest to look for new stations in 2018 that weren't in 2017, then just keep
the desired ones (using the saved `yelloweye_rockfish` data that includes all sets).

```{r stations}
sog_lon_range <- c(-126, -122)
sog_lat_range <- c(48.5, 50.6)

stations_in_2017 <- filter(yelloweye_rockfish$set_counts,
                           year == 2017)

stations_in_2018 <- filter(yelloweye_rockfish$set_counts,
                           year == 2018)

# Reduce to required stations, keeping yelloweye data for map plotting
stations_in_sog_yelloweye <- filter(stations_in_2018,
                                       !(station %in% stations_in_2017$station),
                                       lon > sog_lon_range[1],
                                       lat > sog_lat_range[1],
                                       lat < sog_lat_range[2],
                                       !(station %in% c("2206",
                                                        "2207",
                                                        "2209",
                                                        "2213",
                                                        "2217")))

stations_in_sog <- pull(stations_in_sog_yelloweye,
                           station)
stations_in_sog

# To originally figure out the above stations to exclude:
# stations_in_sog_yelloweye %>% arrange(lat) %>% select(station, lat, lon) %>% as.data.frame()

# Confirm they are all usable:
expect_equal(stations_in_sog_yelloweye$usable,
             rep("Y", nrow(stations_in_sog_yelloweye)))
```

Now plot the Strait stations on a map

```{r map}
plot_BC(xlim = sog_lon_range,
        ylim = sog_lat_range,
        main = "All 2018 stations in Strait of Georgia")

plot_iphc_map(stations_in_sog_yelloweye,
              sp_short_name = NULL,
              years = 2018,
              add_to_existing = TRUE,
              indicate_standard = FALSE,
              include_legend = FALSE)
```

The effective skate number for each station is essentially equivalent to the number of skates each
of 100 hooks. For these `r nrow(stations_in_sog_yelloweye)`
stations ranges from
`r round(min(stations_in_sog_yelloweye$E_it), 3)` to
`r round(max(stations_in_sog_yelloweye$E_it))`
so the sampling was fairly consistent.

## Load cached data and extract numbers for just the Strait

Need to use raw data:

```{r extract}
species_in_sog_2018 <- tibble("species" = character(),
                           "total_number_caught" = numeric(),
                           "average_catch_rate" = numeric())

species_not_in_sog_2018 <- character()

for(sp in sp_vec){      # increment these to test, done up to number shown

  sp_set_counts <- readRDS(paste0(path_data, "/", sp_hyphenate(sp)))

  sp_counts_in_sog_2018 <- filter(sp_set_counts$set_counts,
                                  year == 2018,
                                  station %in% stations_in_sog,
                                  usable == "Y") %>%
    select(station,
           N_it,
           C_it)

  if(nrow(sp_counts_in_sog_2018) == 0){
    species_not_in_sog_2018 <- c(species_not_in_sog_2018,
                                 sp)
  } else {
    if(max(sp_counts_in_sog_2018$N_it > 0)){
      species_in_sog_2018 <- add_row(species_in_sog_2018,
                                  "species" = sp,
                                  "total_number_caught" = sum(sp_counts_in_sog_2018$N_it),
                                  "average_catch_rate" = mean(sp_counts_in_sog_2018$C_it))
    } else {
      species_not_in_sog_2018 <- c(species_not_in_sog_2018,
                                   sp)
    }
  }
}

species_in_sog_2018 <- arrange(species_in_sog_2018,
                               desc(total_number_caught))
```

## Results

So the species that were caught by the IPHC survey in 2018 in the Strait of
Georgia were (ordered by numbers caught)
```{r speciescaught}
knitr::kable(species_in_sog_2018, digits = 3)
```
where the columns give the total number caught across all
`r nrow(stations_in_sog_yelloweye)`
stations, and the average catch rate
(calculated as the catch rate per effective skate for each station, then
averaged over stations). These calculations do not account for hook competition,
competition between species for the bait on the finite number of hooks
(but can be if really needed).

Clearly dogfish were caught in abundance, and Bluntnose Sixgill Shark were the
second highest catch. And only three halibut were caught, which presumably
explains why the IPHC haven't returned to the Strait since.

The species that were not caught, are
```{r speciesnotcaught}
knitr::kable(species_not_in_sog_2018)
```

This is checking all the species listed in the groundfish synopsis report, so includes some that we
obviously don't expect to catch (such as Basking Shark).

## Save as .csv files

```{r save}
write.csv(species_in_sog_2018,
          "species-in-sog-2018.csv",
          row.names = FALSE,
          quote = FALSE)

write.csv(species_not_in_sog_2018,
          "species-not-in-sog-2018.csv",
          row.names = FALSE,
          quote = FALSE)
```

## Plot locations of catch for a specific species

```{r locationsonespecies}
sp <- "bluntnose sixgill shark"
sp_short_name <- "Bluntnose"

sp_set_counts <- readRDS(paste0(path_data, "/", sp_hyphenate(sp)))

sp_counts_in_sog_2018 <- filter(sp_set_counts$set_counts,
                                year == 2018,
                                station %in% stations_in_sog,
                                usable == "Y")

plot_BC(xlim = sog_lon_range,
        ylim = sog_lat_range,
        main = "2018 stations in SoG that caught Bluntnose")

plot_iphc_map(sp_counts_in_sog_2018,
              sp_short_name = sp_short_name,
              years = 2018,
              add_to_existing = TRUE,
              indicate_standard = FALSE)
```