Calculates maximum likelihood estimates of the extremal index \(\theta\) based on the \(K\)-gaps model for threshold inter-exceedances times of Suveges and Davison (2010).

`kgaps(data, u, k = 1, inc_cens = TRUE)`

- data
A numeric vector or numeric matrix of raw data. If

`data`

is a matrix then the log-likelihood is constructed as the sum of (independent) contributions from different columns. A common situation is where each column relates to a different year.If

`data`

contains missing values then`split_by_NAs`

is used to divide the data further into sequences of non-missing values, stored in different columns in a matrix. Again, the log-likelihood is constructed as a sum of contributions from different columns.- u
A numeric scalar. Extreme value threshold applied to data.

- k
A non-negative numeric scalar. Run parameter \(K\), as defined in Suveges and Davison (2010). Threshold inter-exceedances times that are not larger than

`k`

units are assigned to the same cluster, resulting in a \(K\)-gap equal to zero. Specifically, the \(K\)-gap \(S\) corresponding to an inter-exceedance time of \(T\) is given by \(S = \max(T - K, 0)\). In practice, \(k\) should be no smaller than 1, because when \(k\) is less than 1 the estimate of \(\theta\) is always equal to 1.- inc_cens
A logical scalar indicating whether or not to include contributions from right-censored inter-exceedance times, relating to the first and last observations. It is known that these times are greater than or equal to the time observed. See Attalides (2015) for details. If

`data`

has multiple columns then there will be right-censored first and last inter-exceedance times for each column.

An object (a list) of class `c("kgaps", "exdex")`

containing

`theta`

The maximum likelihood estimate (MLE) of \(\theta\).

`se`

The estimated standard error of the MLE, calculated using an algebraic expression for the observed information. If

`k = 0`

then`se`

is returned as`0`

.`se_exp`

The estimated standard error of the MLE, calculated using an algebraic expression for the expected information. If the estimate of \(\theta\) is 0 or 1 then

`se_exp`

is`NA`

.`ss`

The list of summary statistics returned from

`kgaps_stat`

.`k, u, inc_cens`

The input values of

`k`

,`u`

and`inc_cens`

.`max_loglik`

The value of the log-likelihood at the MLE.

`call`

The call to

`kgaps`

.

If `inc_cens = FALSE`

then the maximum likelihood estimate of
the extremal index \(\theta\) under the \(K\)-gaps model of
Suveges and Davison (2010) is calculated.

If `inc_cens = TRUE`

then information from right-censored
first and last inter-exceedance times is also included in the likelihood
to be maximized, following Attalides (2015). The form of the
log-likelihood is given in the **Details** section of
`kgaps_stat`

.

It is possible that the estimate of \(\theta\) is equal to 1, and also
possible that it is equal to 0. `kgaps_stat`

explains the
respective properties of the data that cause these events to occur.

Suveges, M. and Davison, A. C. (2010) Model
misspecification in peaks over threshold analysis, *Annals of
Applied Statistics*, **4**(1), 203-221.
doi:10.1214/09-AOAS292

Attalides, N. (2015) Threshold-based extreme value modelling, PhD thesis, University College London. https://discovery.ucl.ac.uk/1471121/1/Nicolas_Attalides_Thesis.pdf

`kgaps_confint`

to estimate confidence intervals
for \(\theta\).

`kgaps_methods`

for S3 methods for `"kgaps"`

objects.

`kgaps_imt`

for the information matrix test, which
may be used to inform the choice of the pair (`u, k`

).

`choose_uk`

for a diagnostic plot based on
`kgaps_imt`

.

`kgaps_stat`

for the calculation of sufficient
statistics for the \(K\)-gaps model.

`kgaps_post`

in the
`revdbayes`

package for Bayesian inference
about \(\theta\) using the \(K\)-gaps model.

```
### S&P 500 index
u <- quantile(sp500, probs = 0.60)
theta <- kgaps(sp500, u)
theta
#>
#> Call:
#> kgaps(data = sp500, u = u)
#>
#> Estimate of the extremal index theta:
#> theta
#> 0.6953
summary(theta)
#>
#> Call:
#> kgaps(data = sp500, u = u)
#>
#> Estimate Std. Error
#> theta 0.6953 0.007234
coef(theta)
#> theta
#> 0.6953391
nobs(theta)
#> [1] 2901
vcov(theta)
#> theta
#> theta 5.232578e-05
logLik(theta)
#> 'log Lik.' -3811.894 (df=1)
### Newlyn sea surges
u <- quantile(newlyn, probs = 0.60)
theta <- kgaps(newlyn, u, k = 2)
theta
#>
#> Call:
#> kgaps(data = newlyn, u = u, k = 2)
#>
#> Estimate of the extremal index theta:
#> theta
#> 0.1758
summary(theta)
#>
#> Call:
#> kgaps(data = newlyn, u = u, k = 2)
#>
#> Estimate Std. Error
#> theta 0.1758 0.009211
### Cheeseboro wind gusts
theta <- kgaps(cheeseboro, 45, k = 3)
theta
#>
#> Call:
#> kgaps(data = cheeseboro, u = 45, k = 3)
#>
#> Estimate of the extremal index theta:
#> theta
#> 0.2405
summary(theta)
#>
#> Call:
#> kgaps(data = cheeseboro, u = 45, k = 3)
#>
#> Estimate Std. Error
#> theta 0.2405 0.02336
```