Area-weighted index standardization

class: center, middle, inverse, title-slide

.title[
# Area-weighted index standardization
]
.subtitle[
## DFO TESA sdmTMB workshop
]
.author[
###
]
.date[
### January 16–19 2023
]

---

# Why model-based index standardization?

* Design-based methods (stratified means) are straightforward but assume perfect survey implementation

* Model-based approaches can account for (some) aspects of imperfect implementation

* Model-based approaches can leverage spatial correlation to improve precision

* Model-based approaches can integrate additional data sources, multiple gear types, accounting for gear calibration, etc.

---

# Design- vs. model-based estimators

* Design-based: fish are homogenously distributed within strata; noise comes from sampling

* Model-based: fish are distributed according to a mathematical model; noise from uncertainty in fish distribution and sampling

---

# Calculating an area-weighted population index

1. Fit a sensible model: .blue[`sdmTMB()`]

2. Predict on a grid covering the survey domain: .blue[`predict(..., newdata = ...)`]

3. Sum up the density multiplied by cell area to calculate total abundance or biomass: .blue[`get_index()`]

---

# Fit a sensible model

```r
mesh <- make_mesh(pcod, xy_cols = c("X", "Y"), cutoff = 10)
fit <- sdmTMB(
* density ~ 0 + as.factor(year),
  data = pcod,
  mesh = mesh,
  family = tweedie(link = "log"),
  spatial = "on",
  time = "year",
* spatiotemporal = "iid",
  silent = FALSE # show progress!
)
```

---

# Predict over the survey domain

.small[

```r
pred <- predict(
  fit, newdata = qcs_grid, return_tmb_object = TRUE
)
#> Free parallelADFun object.
select(pred$data, year, X, Y, est) |>
  head()
#>   year   X    Y      est
#> 1 2003 456 5636 1.639835
#> 2 2003 458 5636 1.682362
#> 3 2003 460 5636 1.724889
#> 4 2003 462 5636 1.767417
#> 5 2003 464 5636 1.889949
#> 6 2003 466 5636 2.361589
```
]

---

# Predict over the survey domain

---

# Sum up the density multiplied by cell area

(*and calculate standard errors)

```r
index <- get_index(pred, area = 4, bias_correct = TRUE)
head(index)
#>   year       est       lwr       upr  log_est        se
#> 1 2003  936192.2  653707.9 1340745.5 13.74958 0.1832484
#> 2 2004 1832131.6 1359016.7 2469952.2 14.42099 0.1524103
#> 3 2005 1757226.3 1224183.3 2522371.1 14.37925 0.1844231
#> 4 2007  452107.9  328781.9  621693.4 13.02168 0.1625163
#> 5 2009  722994.0  518720.1 1007711.7 13.49116 0.1694095
#> 6 2011 1357912.0 1028868.6 1792187.0 14.12146 0.1415783
```

---

# The resulting standardized index

.xsmall[

```r
ggplot(index, aes(year, est)) + geom_line() +
  geom_ribbon(aes(ymin = lwr, ymax = upr), alpha = 0.4) +
  xlab('Year') + ylab('Biomass estimate (kg)')
```

<img src="09-index-standardization_files/figure-html/pcod-st-index-plot-1.png" width="600px" style="display: block; margin: auto;" />
]

---

# Remember...

There is *one* design-based index.

A good model-based estimated depends on the modeller constructing a reasonable model!

There are *many* possible model-based indexes. Thankfully, results are often (but not always!) qualitatively similar across configurations.

.tiny[
Commander, C.J.C, L.A.K. Barnett, E.J. Ward, S.C. Anderson, T.E. Essington. 2022. The shadow model: how and why small choices in spatially explicit species distribution models affect predictions. PeerJ. 10: e12783. <https://doi.org/10.7717/peerj.12783>
]