---
title: "Estimating extinction dates from sighting records"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Estimating extinction dates from sighting records}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

## Overview

EDE implements methods for inference about extinction from sighting records
and fossil ranges. Some estimate an endpoint directly; others test whether a
species could plausibly have persisted to a candidate time. The tests return
frequentist p-values, not posterior probabilities that the species is extant;
the Bayesian models later in this vignette give those.
Their source papers define no point estimate, so EDE reports the one-sided
upper confidence bound obtained by inverting each test.

The column containing these p-values is currently named `chance` for
compatibility with EDE 0.1.0.

## Input and observation origin

```{r}
library(EDE)

years <- c(1900, 1902, 1903, 1905, 1907, 1908, 1910, 1912,
           1915, 1918, 1920, 1923, 1925, 1928, 1930, 1933, 1936)
counts <- c(4, 3, 5, 2, 3, 4, 2, 1, 2, 1, 1, 2, 1, 1, 1, 1, 1)
sd <- sighting_data(data.frame(year = years, sightings = counts))
```

The earliest supplied time defines the observation origin. Include an initial
zero-count row if observation began before the first sighting. Methods differ
in their treatment of counts, following their source papers: OLE, the
Weibull test, both Solow (1993) tests and Burgman's tests count every
sighting event; the Robson-Whitlock, Strauss-Sadler, McInerny,
Solow-Roberts and Jarić-Ebenhard methods use only distinct sighting times.
Strauss-Sadler, Solow-Roberts and Jarić-Ebenhard do not use the observation
origin at all.

## Endpoint estimators

`robson1964()` uses the final gap for a Quenouille-type bias-reduced point
estimate and a one-sided confidence interval. `strauss1989()` supplies an
unbiased endpoint estimate and one-sided interval under uniform occurrence.
`ole()` fits the Weibull extreme-value model to the `k` most recent sighting
events and reports the two-sided interval of Roberts & Solow (2003).

```{r}
robson1964(sd)
strauss1989(sd)
ole(sd)
ole(sd, k = 10)
```

## Persistence tests

```{r}
solow1993(sd, test_year = 2000)
solow1993b(sd, test_year = 2000)
solow2005(sd, test_year = 2000)
mcinerny2006(sd, test_year = 2000)
burgman1995(sd, test_year = 2000)
solow_roberts2003(sd, test_year = 2000)
jaric2010(sd, test_year = 2000)
```

- `solow1993()` assumes a stationary Poisson process.
- `solow1993b()` allows an exponentially declining sighting rate. Its upper
  bound may not exist: the p-value has a positive floor however late the
  candidate time.
- `solow2005()` uses the Weibull extreme-value p-value associated with OLE.
  Its bound is one-sided, so it equals the upper limit of `ole()` at twice
  the significance level.
- `mcinerny2006()` conditions on the previous binary sighting rate.
- `burgman1995()` uses discrete periods. By default it applies equation 2,
  the test for extinction based on the final run of empty periods;
  `equation = 4` tests for a change in conservation status using the longest
  empty run anywhere in the record.
- `solow_roberts2003()` uses only the two most recent distinct sightings; its
  bound coincides with that of `robson1964()`.
- `jaric2010()` uses the average interval and can adjust it for a trend in
  consecutive interval lengths. It is an index rather than a formal test.

The full p-value curve is available with `data_out = TRUE`:

```{r}
curve <- jaric2010(sd, test_year = 2000, data_out = TRUE)
plot(curve$time, curve$chance, type = "l",
     xlab = "candidate time", ylab = "index")
abline(h = 0.05, lty = 2)
```

## A Bayesian view

Under the constant-rate model, `solow1993_bayes()` gives the Bayes factor of
Solow (1993a) in favour of persistence to a candidate time, assuming a
uniform prior on the extinction time and a noninformative prior on the
sighting rate. Combined with a prior probability of persistence, it gives the
posterior probability that the species is extant, which is the quantity the
p-values above do not provide.

```{r}
solow1993_bayes(sd, test_year = 2000, prior = 0.5)
```

## Uncertain sightings

`jaric_roberts2014()` assigns a probability of validity to every sighting.
EDE requires binary counts because reliability belongs to individual
records; this is a package restriction, not one of the source paper. A
sighting at the observation origin fixes the origin and its reliability is
not used.

```{r}
uncertain <- sighting_data(data.frame(
  year = c(1900, 1910, 1920, 1930),
  sightings = 1
))

jaric_roberts2014(
  uncertain,
  reliability = c(1.0, 0.9, 0.6, 0.3)
)
```

The result includes the effective number of sightings and the
reliability-adjusted endpoint as additional components. Supply `test_year`
to add the reliability-adjusted p-value of equations 7-8. Reliabilities should
be elicited independently of the extinction analysis.

Other approaches to uncertain sightings:

- `thompson2013()` handles several observation types (for example specimens,
  expert sightings and local reports), each with its own range of detection
  and validity probabilities, in discrete time.
- `brook2019()` resamples the record, keeping each sighting with probability
  equal to its reliability, and applies any EDE estimator to every
  sub-record. Solow (2023) criticised the statistical properties of the
  resulting date; the authors present it as pre-processing rather than an
  estimator, and its help page discusses both views.

## Other Bayesian models

`rout2009()` generalises the Solow (1993a) Bayes factor to a sighting rate
that declines as a power of time, and `rout2009_absences()` turns it into a
decision rule: the number of consecutive absent surveys after which declaring
eradication minimises the expected cost. `caley_barry2014()` models yearly
detections with constant detection and extinction probabilities; its
posterior is computed exactly:

```{r}
fox <- sighting_data(data.frame(year = c(2001, 2003, 2005, 2006), count = 1))
caley_barry2014(fox, test_year = 2012)
```

## Records and surveys

`thompson2017()` implements the iterative model of Thompson et al. (2017),
which underlies the IUCN Red List tool for estimating the probability that a
taxon is extinct. Each year contributes either a record, with a probability
of being correct, or a survey that failed to find the species, described by
its extent, its chance of detecting the species and its reliability. Ranges
for these quantities give interval and Monte Carlo bounds. `thompson2019()`
casts the same data as Bayesian updating, in which each year's posterior is
the next year's prior.

## Fossil ranges

The same question arises for the last appearance of a taxon in the fossil
record. `strauss1989()` assumes fossil horizons are uniformly distributed;
`marshall1994()` drops that assumption and bounds the confidence interval
using only the distribution of gap sizes, with `marshall1994_trend()` to test
whether gap size changes through the section; `marshall1997()` accepts a
fossil recovery potential that varies through the section; and
`weiss_marshall1999()` gives a Bayesian estimate for sections sampled at
discrete horizons, with absences recorded. EDE always estimates the end-point
at the largest `time`, so ages in millions of years should be entered as
`time = -age`, and first appearances as negative stratigraphic positions.

## Comparing methods

`ede_compare()` applies every method that needs only the sighting record:

```{r}
ede_compare(sighting_data(monk_seal), test_year = 2050)
```

## Choosing a method

- Use OLE when the recent record is sufficiently long and an extreme-value
  model is defensible.
- Use `solow1993()` only when a constant sighting rate is plausible.
- Use `solow1993b()` for an approximately exponential decline.
- Use `jaric2010()` when the intervals themselves show a gradual trend.
- Use `burgman1995()` for discrete frequency data; use `equation = 4` only
  for questions about a change in status.
- Use `jaric_roberts2014()` when observations have explicit reliability
  assessments, and `thompson2013()` when they fall into types with uncertain
  detection probabilities.
- Use `thompson2017()` when unsuccessful surveys are documented, as in an
  IUCN Red List assessment.
- Use the Bayesian models when a posterior probability that the species is
  extant is needed; report the prior alongside it.

Running several scientifically defensible methods is useful, but p-values from
models whose assumptions are violated should not be combined or interpreted as
extinction probabilities.

## References

References for the methods introduced after the core estimators (Rout et
al. 2009; Thompson et al. 2013, 2017, 2019; Caley & Barry 2014; Brook et al.
2019; Marshall 1994, 1997; Weiss & Marshall 1999) are given on the help page
of each function.

Burgman, M. A., Grimson, R. C., & Ferson, S. (1995). Inferring threat from
scientific collections. *Conservation Biology*, 9(4), 923-928.

Cooke, P. (1980). Optimal linear estimation of bounds of random variables.
*Biometrika*, 67, 257-258.

Jarić, I., & Ebenhard, T. (2010). A method for inferring extinction based on
sighting records that change in frequency over time. *Wildlife Biology*,
16(3), 267-275.

Jarić, I., & Roberts, D. L. (2014). Accounting for observation reliability
when inferring extinction based on sighting records. *Biodiversity and
Conservation*, 23(11), 2801-2815.

McInerny, G. J., Roberts, D. L., Davy, A. J., & Cribb, P. J. (2006).
Significance of sighting rate in inferring extinction and threat.
*Conservation Biology*, 20(2), 562-567.

Roberts, D. L., & Solow, A. R. (2003). When did the dodo become extinct?
*Nature*, 426, 245.

Robson, D. S., & Whitlock, J. H. (1964). Estimation of a truncation point.
*Biometrika*, 51, 33-39.

Solow, A. R. (1993a). Inferring extinction from sighting data. *Ecology*, 74,
962-964.

Solow, A. R. (1993b). Inferring extinction in a declining population.
*Journal of Mathematical Biology*, 32, 79-82.

Solow, A. R. (2005). Inferring extinction from a sighting record.
*Mathematical Biosciences*, 195, 47-55.

Solow, A. R., & Roberts, D. L. (2003). A nonparametric test for extinction
based on a sighting record. *Ecology*, 84, 1329-1332.

Strauss, D., & Sadler, P. M. (1989). Classical confidence intervals and
Bayesian probability estimates for ends of local taxon ranges.
*Mathematical Geology*, 21, 411-427.