p.exceed: probability of exceeding the threshold u; for the Smith estimator, this is mean(x > threshold) for x
being the data.
shape: GPD shape parameter xi (a real number).
scale: GPD scale parameter beta (a positive number).
lower.tail: logical; if TRUE (default) probabilities are P(X<=x) otherwise, P(X>x).
log, log.p: logical; if TRUE, probabilities p are given as log(p).
Returns
dGPDtail() computes the density, pGPDtail() the distribution function, qGPDtail() the quantile function and rGPDtail() random variates of the GPD-based tail distribution in the POT method.
Details
Let u denote the threshold (threshold), pu the exceedance probability (p.exceed) and FGPD the GPD distribution function. Then the distribution function of the GPD-based tail distribution is given by
F(q)=1−pu(1−FGPD(q−u))
. The quantile function is
F−1(p)=u+FGPD−1(1−(1−p)/pu)
and the density is
f(x)=pufGPD(x−u)
, where fGPD denotes the GPD density.
Note that the distribution function has a jumpt of height P(X<=u) (1-p.exceed) at u.
Author(s)
Marius Hofert
References
McNeil, A. J., Frey, R., and Embrechts, P. (2015). Quantitative Risk Management: Concepts, Techniques, Tools. Princeton University Press.
Examples
## Generate data to work withset.seed(271)X <- rt(1000, df =3.5)# in MDA(H_{1/df}); see MFE (2015, Section 16.1.1)## Determine thresholds for POT methodmean_excess_plot(X[X >0])abline(v =1.5)u <-1.5# threshold## Fit GPD to the excesses (per margin)fit <- fit_GPD_MLE(X[X > u]- u)fit$par
1/fit$par["shape"]# => close to df## Estimate threshold exceedance probabilitiesp.exceed <- mean(X > u)## Define corresponding densities, distribution function and RNGdF <-function(x) dGPDtail(x, threshold = u, p.exceed = p.exceed, shape = fit$par["shape"], scale = fit$par["scale"])pF <-function(q) pGPDtail(q, threshold = u, p.exceed = p.exceed, shape = fit$par["shape"], scale = fit$par["scale"])rF <-function(n) rGPDtail(n, threshold = u, p.exceed = p.exceed, shape = fit$par["shape"], scale = fit$par["scale"])## Basic check of dF()curve(dF, from = u -1, to = u +5)## Basic check of pF()curve(pF, from = u, to = u +5, ylim =0:1)# quite flat hereabline(v = u, h =1-p.exceed, lty =2)# mass at u is 1-p.exceed (see 'Details')## Basic check of rF()set.seed(271)X. <- rF(1000)plot(X., ylab ="Losses generated from the fitted GPD-based tail distribution")stopifnot(all.equal(mean(X. == u),1-p.exceed, tol =7e-3))# confirms the above## Pick out 'continuous part'X.. <- X.[X. > u]plot(pF(X..), ylab ="Probability-transformed tail losses")# should be U[1-p.exceed, 1]