---
title: "Variance Methods: Taylor, Bootstrap, and Jackknife"
output:
  rmarkdown::html_vignette:
    highlight: null
vignette: >
  %\VignetteIndexEntry{Variance Methods: Taylor, Bootstrap, and Jackknife}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
```

Every estimator in tidycreel takes a `variance` argument with three choices:
`"taylor"` (the default), `"bootstrap"` and `"jackknife"`. This page explains
what the three have in common, how each one works, and when they can
disagree.

```{r load, message = FALSE}
library(tidycreel)
```

## What all three estimate

All three estimate the same thing: how much the estimate would vary if a
different set of days had been sampled under the same design. That is a
**design-based** variance. None of the three assumes the counts or the catch
follow any particular distribution; they differ only in how they compute it.

- **Taylor linearization** uses formulas. For a total these are the
  between-day and within-day terms of `vignette("effort-pipeline")`; for a
  rate or a product of estimates, the estimate is approximated by a linear
  function of totals (a first-order Taylor series) and the variance of that
  is computed.
- **Jackknife** drops one sampled day from a stratum, re-computes the
  estimate, puts the day back, and repeats for every sampled day. The
  variance comes from how much the re-computed estimates spread. tidycreel
  uses the stratified delete-one jackknife (JKn) from the survey package.
- **Bootstrap** draws sampled days with replacement within each stratum, 500
  times, re-computing the estimate each time, using Canty and Davison's (1999)
  bootstrap from the survey package. It is random, so set a seed for a
  reproducible standard error.

All three resample or vary **days**, the primary sampling unit, never single
counts or interviews.

## Same season, three methods

The season from `vignette("tidycreel")`: 120 days, 36 of them sampled, 28
weekdays and 8 weekend days.

```{r design, message = FALSE}
design <- creel_design(example_season_calendar, date = date, strata = day_type)
design <- add_counts(design, example_season_counts,
  count_col = total_anglers,
  period_length_col = daylight_hours
)
design <- add_interviews(design, example_season_interviews,
  catch = catch_total,
  effort = hours_fished,
  trip_status = trip_status,
  n_anglers = n_anglers
)
```

```{r compare, message = FALSE}
set.seed(2031)
methods <- c("taylor", "bootstrap", "jackknife")
data.frame(
  variance = methods,
  effort_se = sapply(methods, function(v) {
    estimate_effort(design, target = "period_total", variance = v)$estimates$se
  }),
  rate_se = sapply(methods, function(v) {
    estimate_catch_rate(design, variance = v)$estimates$se
  }),
  total_catch_se = sapply(methods, function(v) {
    estimate_total_catch(design, target = "period_total", variance = v)$estimates$se
  }),
  row.names = NULL
)
```

The point estimates are identical; only the standard errors change, and here
by no more than a few percent. For effort, a total, the jackknife reproduces
the Taylor value exactly. The bootstrap differs a little, and would differ a
little again with another seed.

## When they can disagree

Taylor linearization for a rate or a product rests on a first-order
approximation. Wolter (2007, chapter 6) reports that with enough sampled units
it is usually reliable, but that it may not be for highly skewed populations,
and catch data often are: most parties catch nothing and a few catch many.
When the three methods disagree noticeably on your data, treat that as a
finding: the standard error is less certain than any one of them suggests.

There is no rule that one method is always better. Taylor is the default
because it is fast and deterministic; the replicate methods are a check on it.

## One sampled day in a stratum

With one sampled day in a stratum there is no day-to-day variation to measure,
so no method can estimate that stratum's variance. All three refuse:

```{r single-psu, error = TRUE}
one_day <- creel_design(
  data.frame(
    date     = as.Date(c("2024-06-03", "2024-06-04", "2024-06-08")),
    day_type = c("weekday", "weekday", "weekend")
  ),
  date = date, strata = day_type
)
one_day <- add_counts(one_day,
  data.frame(
    date      = as.Date(c("2024-06-03", "2024-06-04", "2024-06-08")),
    day_type  = c("weekday", "weekday", "weekend"),
    n_anglers = c(10L, 14L, 15L),
    hours     = 10
  ),
  period_length_col = hours
)

estimate_effort(one_day, variance = "taylor")
```

`add_counts()` already warned, when the counts were attached, that both
strata have fewer than three sampled days; the estimate is where the weekend's
single day becomes a refusal. The remedies are the ones the error gives: sample more days in that stratum,
or merge it with a similar stratum before estimating. (The jackknife and
bootstrap errors currently suggest switching to Taylor, which refuses the
same design; that is #503.)

## References

- Canty, A. J., and Davison, A. C. (1999). Resampling-based variance
  estimation for labour force surveys. *The Statistician*, 48, 379–391.
- Wolter, K. M. (2007). *Introduction to Variance Estimation*, 2nd edition.
  Springer, New York.