Estimator of small exceedance probabilities and large return periods using censored Hill
Estimator of small exceedance probabilities and large return periods using censored Hill
Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the estimates for the EVI obtained from the Hill estimator adapted for right censoring.
cProb(data, censored, gamma1, q, plot =FALSE, add =FALSE, main ="Estimates of small exceedance probability",...)cReturn(data, censored, gamma1, q, plot =FALSE, add =FALSE, main ="Estimates of large return period",...)
Arguments
data: Vector of n observations.
censored: A logical vector of length n indicating if an observation is censored.
gamma1: Vector of n−1 estimates for the EVI obtained from cHill.
q: The used large quantile (we estimate P(X>q) or 1/P(X>q) for q large).
plot: Logical indicating if the estimates should be plotted as a function of k, default is FALSE.
add: Logical indicating if the estimates should be added to an existing plot, default is FALSE.
main: Title for the plot, default is "Estimates of small exceedance probability" for cProb
and "Estimates of large return period" for cReturn.
...: Additional arguments for the plot function, see plot for more details.
Details
The probability is estimated as
P^(X>q)=(1−km)×(q/Zn−k,n)−1/Hk,nc
with Zi,n the i-th order statistic of the data, Hk,nc
the Hill estimator adapted for right censoring and km the Kaplan-Meier estimator for the CDF evaluated in Zn−k,n.
Returns
A list with following components: - k: Vector of the values of the tail parameter k.
P: Vector of the corresponding probability estimates, only returned for cProb.
R: Vector of the corresponding estimates for the return period, only returned for cReturn.
q: The used large quantile.
References
Beirlant, J., Guillou, A., Dierckx, G. and Fils-Villetard, A. (2007). "Estimation of the Extreme Value Index and Extreme Quantiles Under Random Censoring." Extremes, 10, 151--174.
Author(s)
Tom Reynkens
See Also
cHill, cQuant, Prob, KaplanMeier
Examples
# Set seedset.seed(29072016)# Pareto random sampleX <- rpareto(500, shape=2)# Censoring variableY <- rpareto(500, shape=1)# Observed sampleZ <- pmin(X, Y)# Censoring indicatorcensored <-(X>Y)# Hill estimator adapted for right censoringchill <- cHill(Z, censored=censored, plot=TRUE)# Small exceedance probabilityq <-10cProb(Z, censored=censored, gamma1=chill$gamma1, q=q, plot=TRUE)# Return periodcReturn(Z, censored=censored, gamma1=chill$gamma1, q=q, plot=TRUE)