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)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.
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.
A numeric scalar. Extreme value threshold applied to data.
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.
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.
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. tools:::Rd_expr_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
summary(theta)
coef(theta)
nobs(theta)
vcov(theta)
logLik(theta)
### Newlyn sea surges
u <- quantile(newlyn, probs = 0.60)
theta <- kgaps(newlyn, u, k = 2)
theta
summary(theta)
### Cheeseboro wind gusts
theta <- kgaps(cheeseboro, 45, k = 3)
theta
summary(theta)
Run the code above in your browser using DataLab