mean_excess {qrmtools} | R Documentation |
Sample mean excess function, mean excess function of a GPD and sample mean excess plot.
mean_excess_np(x, omit = 3) mean_excess_plot(x, omit = 3, xlab = "Threshold", ylab = "Mean excess over threshold", ...) mean_excess_GPD(x, shape, scale)
x |
|
omit |
number >= 1 of unique last observations to be omitted from the sorted data (as mean excess plot becomes unreliable for these observations as thresholds). |
xlab |
x-axis label. |
ylab |
y-axis label. |
... |
additional arguments passed to the underlying
|
shape |
GPD shape parameter xi. |
scale |
GPD scale parameter beta. |
Mean excess plots can be used in the peaks-over-threshold method for choosing a threshold. To this end, one chooses the smallest threshold above which the mean excess plot is roughly linear.
mean_excess_np()
returns a two-column matrix giving
the sorted data without the omit
-largest unique values
(first column) and the corresponding values of the sample mean excess
function (second column). It is mainly used in mean_excess_plot()
.
mean_excess_plot()
returns invisible()
.
mean_excess_GPD()
returns the mean excess function of a
generalized Pareto distribution evaluated at x
.
Marius Hofert
## (Sample) mean excess function data(fire) ME <- mean_excess_np(fire) stopifnot(dim(ME) == c(2164, 2), all.equal(ME[nrow(ME),], c(65.707491, 121.066231), check.attributes = FALSE)) ## A 'manual' (sample) mean excess plot plot(ME, xlab = "Threshold", ylab = "Mean excess over threshold") ## (Sample) mean excess plot mean_excess_plot(fire) ## => Any value in [10, 20] seems reasonable here as threshold choice ## (one would probably go with 10 to benefit from a larger sample size). ## With mean excess functions of two fitted GPDs overlaid u <- c(10, 20) # thresholds fit <- lapply(u, function(u.) fit_GPD_MLE(fire[fire > u.] - u.)) q <- lapply(u, function(u.) seq(u., ME[nrow(ME),"x"], length.out = 129)) MEF.GPD <- lapply(1:2, function(k) mean_excess_GPD(q[[k]]-u[k], shape = fit[[k]]$par[["shape"]], scale = fit[[k]]$par[["scale"]])) mean_excess_plot(fire, ylim = range(ME, unlist(MEF.GPD))) col <- c("royalblue3", "maroon3") for(k in 1:2) lines(q[[k]], MEF.GPD[[k]], col = col[k]) legend("bottomright", col = rev(col), lty = rep(1, length(u)), bty = "n", legend = as.expression(sapply(rev(seq_along(u)), function(k) substitute("Threshold choice"~~u==u., list(u. = u[k])))))