| Title: | Extreme Values in R |
|---|---|
| Description: | Functions for extreme value theory, which may be divided into the following groups; exploratory data analysis, block maxima, peaks over thresholds (univariate and bivariate), point processes, gev/gpd distributions. |
| Authors: | Bernhard Pfaff [aut, cre], Alexander McNeil [aut] (S original (EVIS)), Alec Stephenson [trl] (R port of EVIS) |
| Maintainer: | Bernhard Pfaff <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 1.7-4 |
| Built: | 2024-11-01 11:16:31 UTC |
| Source: | https://github.com/bpfaff/evir |
These data are the daily log returns on BMW share price
from Tuesday 2nd January 1973 until Tuesday 23rd July 1996.
The data are contained in a numeric vector.
The dates of each observation are contained in a times
attribute, which is an object of class "POSIXct" (see
DateTimeClasses).
Note that these data form an irregular time series because
no trading takes place at the weekend.
data(bmw)data(bmw)
A numeric vector containing 6146 observations, with a
times attribute which is a POSIXct object
of the same length.
These data describe large fire insurance claims
in Denmark from Thursday 3rd January 1980 until
Monday 31st December 1990.
The data are contained in a numeric vector.
The dates of each observation are contained in a times
attribute, which is an object of class "POSIXct" (see
DateTimeClasses).
They were supplied by Mette Rytgaard of Copenhagen Re.
Note that these data form an irregular time series.
data(danish)data(danish)
A numeric vector containing 2167 observations, with a
times attribute which is a POSIXct object
of the same length.
Declusters clustered point process data so that Poisson assumption is more tenable over a high threshold.
decluster(series, run = NA, picture = TRUE)decluster(series, run = NA, picture = TRUE)
series |
a numeric vector of threshold exceedances
with a |
run |
parameter to be used in the runs method; any two consecutive threshold exceedances separated by more than this number of observations/days are considered to belong to different clusters |
picture |
whether or not a picture of declustering should be drawn |
The declustered object.
Embrechts, P., Klueppelberg, C., Mikosch, T. (1997). Modelling Extremal Events. Springer. Chapter 8, 413–429.
# decluster the 200 exceedances of a particular threshold in # the negative BMW data data(bmw) out <- pot(-bmw, ne = 200) decluster(out$data, 30)# decluster the 200 exceedances of a particular threshold in # the negative BMW data data(bmw) out <- pot(-bmw, ne = 200) decluster(out$data, 30)
Cumulative probability, quantiles, density and random generation from the generalized extreme value distribution.
pgev(q, xi = 1, mu = 0, sigma = 1) qgev(p, xi = 1, mu = 0, sigma = 1) dgev(x, xi = 1, mu = 0, sigma = 1) rgev(n, xi = 1, mu = 0, sigma = 1)pgev(q, xi = 1, mu = 0, sigma = 1) qgev(p, xi = 1, mu = 0, sigma = 1) dgev(x, xi = 1, mu = 0, sigma = 1) rgev(n, xi = 1, mu = 0, sigma = 1)
q |
vector of quantiles |
p |
vector of probabilities |
x |
vector of values at which to evaluate density |
n |
sample size |
xi |
shape parameter |
mu |
location parameter |
sigma |
scale parameter |
Probability (pgev), quantile (qgev), density
(dgev) or random sample (rgev) for the GEV
distribution with shape xi.
Cumulative probability, quantiles, density and random generation from the generalized Pareto distribution.
pgpd(q, xi, mu = 0, beta = 1) qgpd(p, xi, mu = 0, beta = 1) dgpd(x, xi, mu = 0, beta = 1) rgpd(n, xi, mu = 0, beta = 1)pgpd(q, xi, mu = 0, beta = 1) qgpd(p, xi, mu = 0, beta = 1) dgpd(x, xi, mu = 0, beta = 1) rgpd(n, xi, mu = 0, beta = 1)
q |
vector of quantiles. |
p |
vector of probabilities. |
x |
vector of values at which to evaluate density |
n |
sample size |
xi |
shape parameter. |
mu |
location parameter. |
beta |
scale parameter |
Probability (pgpd), quantile (qgpd), density
(dgpd) or random sample (rgpd) for the GPD
distribution with shape xi.
Plots empirical distribution function of a sample.
emplot(data, alog = "x", labels = TRUE, ...)emplot(data, alog = "x", labels = TRUE, ...)
data |
data vector |
alog |
whether axes are to be logged: |
labels |
whether or not axes should be labelled |
... |
other graphics parameters |
This is a simple explanatory function. A straight line on the double log scale indicates Pareto tail behaviour.
## Not run: data(danish) ## Not run: emplot(danish) # Danish fire insurance data show Pareto tail behaviour## Not run: data(danish) ## Not run: emplot(danish) # Danish fire insurance data show Pareto tail behaviour
Plot estimates of extremal index using the blocks method.
exindex(data, block, start = 5, end = NA, reverse = FALSE, auto.scale = TRUE, labels = TRUE, ...)exindex(data, block, start = 5, end = NA, reverse = FALSE, auto.scale = TRUE, labels = TRUE, ...)
data |
data vector (raw values not block maxima). |
block |
the block size. A numeric value is interpreted as
the number of data values in each successive block.
All the data is used, so the last block may not contain
|
start |
lowest value of K at which to plot a point; K is the number of blocks in which a specified threshold is exceeded |
end |
highest value of K at which to plot a point |
reverse |
whether plot is to be by increasing threshold
( |
auto.scale |
whether or not plot should be automatically
scaled; if not, |
labels |
whether or not axes should be labelled |
... |
other graphics parameters |
A table of results is returned invisibly.
Embrechts, P., Klueppelberg, C., Mikosch, T. (1997). Modelling Extremal Events. Springer. Chapter 8, 413-429.
## Not run: data(bmw) ## Not run: exindex(bmw, 100) ## Not run: exindex(-bmw, 100) # calculate extremal index for the right and left tails of the BMW # log returns## Not run: data(bmw) ## Not run: exindex(bmw, 100) ## Not run: exindex(-bmw, 100) # calculate extremal index for the right and left tails of the BMW # log returns
Finds a threshold so that a given number of extremes lie above.
findthresh(data, ne)findthresh(data, ne)
data |
data vector |
ne |
vector giving number of extremes above the threshold |
When the data are tied a threshold is found so that at least the specified number of extremes lie above.
Vector of suitable thresholds.
# Find threshold giving (at least) fifty exceedances for Danish data data(danish) findthresh(danish, 50)# Find threshold giving (at least) fifty exceedances for Danish data data(danish) findthresh(danish, 50)
Fits generalized extreme value distribution (GEV) to block maxima data.
gev(data, block = NA, ...)gev(data, block = NA, ...)
data |
data vector. Interpretation depends on value of block: if no block size is specified then data are interpreted as block maxima; if block size is set, then data are interpreted as raw data and block maxima are calculated. |
block |
the block size. A numeric value is interpreted as the
number of data values in each successive block. All the data is
used, so the last block may not contain |
... |
arguments passed to |
An object of class gev describing the fit and including
parameter estimates and standard errors.
Fitting is carried out using maximum likelihood.
plot.gev, gumbel,
optim, as.POSIXct
# Fit GEV to monthly maxima data(bmw) out <- gev(bmw, "month") # Fit GEV to maxima of blocks of 100 observations out <- gev(bmw, 100) # Fit GEV to the data in nidd.annual, the annual maximum water # levels of the River Nidd, using the "BFGS" optimization method data(nidd.annual) out <- gev(nidd.annual, method = "BFGS", control = list(maxit = 500))# Fit GEV to monthly maxima data(bmw) out <- gev(bmw, "month") # Fit GEV to maxima of blocks of 100 observations out <- gev(bmw, 100) # Fit GEV to the data in nidd.annual, the annual maximum water # levels of the River Nidd, using the "BFGS" optimization method data(nidd.annual) out <- gev(nidd.annual, method = "BFGS", control = list(maxit = 500))
Returns an object of class "gpd" representing the fit of
a generalized Pareto model to excesses over a high threshold.
gpd(data, threshold = NA, nextremes = NA, method = c("ml", "pwm"), information = c("observed", "expected"), ...)gpd(data, threshold = NA, nextremes = NA, method = c("ml", "pwm"), information = c("observed", "expected"), ...)
data |
data vector |
threshold |
a threshold value (either this or
|
nextremes |
the number of upper extremes to be used (either
this or |
method |
whether parameters should be estimated by the maximum
likelihood method |
information |
whether standard errors should be calculated with
|
... |
arguments passed to |
The function uses the general purpose optimization function
optim when method = "ml" is chosen.
An object of class "gpd" describing the fit and including
parameter estimates and standard errors.
Parameter and quantile estimation for the generalized Pareto distribution, JRM Hosking and JR Wallis, Technometrics 29(3), pages 339-349, 1987.
data(danish) out <- gpd(danish, 10) # Fits GPD to excess losses over 10 for the Danish # fire insurance datadata(danish) out <- gpd(danish, 10) # Fits GPD to excess losses over 10 for the Danish # fire insurance data
Calculates quantile estimates and confidence intervals for high quantiles above the threshold in a GPD analysis, and adds a graphical representation to an existing plot.
gpd.q(x, pp, ci.type = c("likelihood", "wald"), ci.p = 0.95, like.num = 50)gpd.q(x, pp, ci.type = c("likelihood", "wald"), ci.p = 0.95, like.num = 50)
x |
a list object returned by |
pp |
the desired probability for quantile estimate (e.g. 0.99 for the 99th percentile) |
ci.type |
method for calculating a confidence interval:
|
ci.p |
probability for confidence interval (must be less than 0.999) |
like.num |
number of times to evaluate profile likelihood |
The GPD approximation in the tail is used to estimate quantile.
The "wald" method uses the observed Fisher information
matrix to calculate confidence interval. The "likelihood"
method reparametrizes the likelihood in terms of the unknown
quantile and uses profile likelihood arguments to construct a
confidence interval.
gpd, plot.gpd,
gpd.sfall, tailplot
## Not run: data(danish) ## Not run: out <- gpd(danish, 10) ## Not run: tp <- tailplot(out) ## Not run: gpd.q(tp, 0.999) # Estimates 99.9th percentile of Danish fire losses## Not run: data(danish) ## Not run: out <- gpd(danish, 10) ## Not run: tp <- tailplot(out) ## Not run: gpd.q(tp, 0.999) # Estimates 99.9th percentile of Danish fire losses
Calculates expected shortfall (tail conditional expectation) estimates and confidence intervals for high quantiles above the threshold in a GPD analysis, and adds a graphical representation to an existing plot.
gpd.sfall(x, pp, ci.p = 0.95, like.num = 50)gpd.sfall(x, pp, ci.p = 0.95, like.num = 50)
x |
a list object returned by |
pp |
the desired probability for expected shortfall estimate (e.g. 0.99 for the 99th percentile) |
ci.p |
probability for confidence interval (must be less than 0.999) |
like.num |
number of times to evaluate profile likelihood |
Expected shortfall is the expected size of the loss, given that a particular quantile of the loss distribution is exceeded. The GPD approximation in the tail is used to estimate expected shortfall. The likelihood is reparametrised in terms of the unknown expected shortfall and profile likelihood arguments are used to construct a confidence interval.
gpd, plot.gpd,
tailplot, gpd.q
## Not run: data(danish) ## Not run: out <- gpd(danish, 10) ## Not run: tp <- tailplot(out) ## Not run: gpd.q(tp, 0.999) # Estimates 99.9th percentile of Danish fire losses ## Not run: gpd.sfall(tp, 0.999) # Estimates associated expected shortfall for Danish fire losses## Not run: data(danish) ## Not run: out <- gpd(danish, 10) ## Not run: tp <- tailplot(out) ## Not run: gpd.q(tp, 0.999) # Estimates 99.9th percentile of Danish fire losses ## Not run: gpd.sfall(tp, 0.999) # Estimates associated expected shortfall for Danish fire losses
Returns an object of class "gpdbiv" representing the fit
of a bivariate POT (peaks over thresholds) model for joint
excesses over thresholds.
gpdbiv(data1 = NA, data2 = NA, u1 = NA, u2 = NA, ne1 = NA, ne2 = NA, global = FALSE, method = "BFGS", ...)gpdbiv(data1 = NA, data2 = NA, u1 = NA, u2 = NA, ne1 = NA, ne2 = NA, global = FALSE, method = "BFGS", ...)
data1 |
first data vector |
data2 |
second data vector |
u1 |
threshold for |
u2 |
threshold for |
ne1 |
number of upper extremes to be used for |
ne2 |
number of upper extremes to be used for |
global |
should a global maximisation of the likelihood with respect to marginal and dependence parameters be undertaken. The default alternative is a two-stage local fit where first the marginal parameters are estimated and then the dependence parameter. This is much faster than a global fit. |
method |
the optimization method (see |
... |
other arguments passed to |
This function implements a model suggested by Richard Smith (see references below). The marginal excess distributions are GPD distributions, as suggested by univariate EVT and implemented in gpd. The dependence specification is known as the logistic or Gumbel dependence structure, but it would be easy to program alternatives.
An object of class "gpdbiv" representing the fit and including
parameter estimates and standard errors.
Multivariate Threshold Methods, Richard L. Smith, in Extreme Value Theory and Applications, ed. J. Galambos, published by Kluwer, pages 225-248, 1994.
Markov Chain Models for Threshold Exceedances, R.L. Smith, J.A. Tawn, S.G. Coles, Biometrika 84, 249-268, 1997.
gpd, plot.gpdbiv,
interpret.gpdbiv
data(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) interpret.gpdbiv(out, 0.05, 0.05) ## Not run: plot(out)data(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) interpret.gpdbiv(out, 0.05, 0.05) ## Not run: plot(out)
Fits gumbel distribution (GEV with xi = 0) to
block maxima data.
gumbel(data, block = NA, ...)gumbel(data, block = NA, ...)
data |
data vector. Interpretation depends on value of block: if no block size is specified then data are interpreted as block maxima; if block size is set, then data are interpreted as raw data and block maxima are calculated. |
block |
the block size. A numeric value is interpreted as the
number of data values in each successive block. All the data is
used, so the last block may not contain |
... |
arguments passed to |
This function is primarily intended for comparison with GEV for assessing the need for a heavy-tailed Frechet (or short-tailed Weibull) to model block maxima.
An object of class "gev" describing the fit and including
parameter estimates and standard errors.
Fitting is carried out using maximum likelihood.
plot.gev, gev,
optim, as.POSIXct
# Fit Gumbel to maxima of blocks of 100 observations data(bmw) out <- gumbel(bmw, 100) # Fit Gumbel to the data in nidd.annual, the annual maximum water # levels of the River Nidd, using the "BFGS" optimization method data(nidd.annual) out <- gumbel(nidd.annual, method = "BFGS", control = list(maxit = 500))# Fit Gumbel to maxima of blocks of 100 observations data(bmw) out <- gumbel(bmw, 100) # Fit Gumbel to the data in nidd.annual, the annual maximum water # levels of the River Nidd, using the "BFGS" optimization method data(nidd.annual) out <- gumbel(nidd.annual, method = "BFGS", control = list(maxit = 500))
Plot the Hill estimate of the tail index of heavy-tailed data, or of an associated quantile estimate.
hill(data, option = c("alpha","xi","quantile"), start = 15, end = NA, reverse = FALSE, p = NA, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)hill(data, option = c("alpha","xi","quantile"), start = 15, end = NA, reverse = FALSE, p = NA, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)
data |
data vector |
option |
whether |
start |
lowest number of order statistics at which to plot a point |
end |
highest number of order statistics at which to plot a point |
reverse |
whether plot is to be by increasing threshold
( |
p |
probability required when option |
ci |
probability for asymptotic confidence band; for no
confidence band set |
auto.scale |
whether or not plot should be automatically
scaled; if not, |
labels |
whether or not axes should be labelled |
... |
other graphics parameters |
This plot is usually calculated from the alpha perspective.
For a generalized Pareto analysis of heavy-tailed data using
the gpd function, it helps to plot the Hill estimates
for xi.
## Not run: data(danish) ## Not run: hill(danish) # Hill plot of heavy-tailed Danish fire insurance data ## Not run: hill(danish, option = "quantile", end = 500, p = 0.999) # Hill plot of estimated 0.999 quantile of Danish fire insurance data## Not run: data(danish) ## Not run: hill(danish) # Hill plot of heavy-tailed Danish fire insurance data ## Not run: hill(danish, option = "quantile", end = 500, p = 0.999) # Hill plot of estimated 0.999 quantile of Danish fire insurance data
Interprets the results of a bivariate GPD model fitted using the bivariate POT method.
interpret.gpdbiv(out, x, y)interpret.gpdbiv(out, x, y)
out |
a |
x |
a scalar value greater than first threshold |
y |
a scalar value greater than second threshold |
First marginal probabilities of exceeding the points x and y are calculated, and then joint and conditional probabilities.
A vector of probabilities is invisibly returned, in printed order.
A simple interpretation of the fit in terms of exceedance probabilities for the point (x,y) is printed.
data(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) interpret.gpdbiv(out, 0.05, 0.05) # probabilities of 5% falls in BMW and Siemens stock pricesdata(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) interpret.gpdbiv(out, 0.05, 0.05) # probabilities of 5% falls in BMW and Siemens stock prices
Plots sample mean excesses over increasing thresholds.
meplot(data, omit = 3, labels = TRUE, ...)meplot(data, omit = 3, labels = TRUE, ...)
data |
data vector |
omit |
number of upper plotting points to be omitted |
labels |
whether or not axes are to be labelled |
... |
other graphics parameters |
An upward trend in plot shows heavy-tailed behaviour. In particular, a straight line with positive gradient above some threshold is a sign of Pareto behaviour in tail. A downward trend shows thin-tailed behaviour whereas a line with zero gradient shows an exponential tail. Because upper plotting points are the average of a handful of extreme excesses, these may be omitted for a prettier plot.
## Not run: data(danish) ## Not run: meplot(danish) # Sample mean excess plot of heavy-tailed Danish fire insurance data## Not run: data(danish) ## Not run: meplot(danish) # Sample mean excess plot of heavy-tailed Danish fire insurance data
These data represent annual maximal levels of the
River Nidd in Yorkshire.
These data are suitable for analysis with gev.
data(nidd.annual)data(nidd.annual)
A numeric vector containing 35 observations.
These data represent high river levels of the
River Nidd in Yorkshire above a threshold value
of 65. These data are suitable for analysis with
gpd.
data(nidd.thresh)data(nidd.thresh)
A numeric vector containing 154 observations.
The plot method plot.gev provides two different residual plots for assessing fitted GEV model. The user selects the plot type from a menu. See the examples below.
## S3 method for class 'gev' plot(x, ...)## S3 method for class 'gev' plot(x, ...)
x |
a |
... |
other graphics parameters |
Data are converted to unit exponentially distributed residuals under null hypothesis that GEV fits. Two diagnostics for iid exponential data are offered.
data(bmw) out <- gev(bmw, 100) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Scatterplot of Residuals ## Not run: 2: plot: QQplot of Residualsdata(bmw) out <- gev(bmw, 100) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Scatterplot of Residuals ## Not run: 2: plot: QQplot of Residuals
The plot method plot.gpd provides four different plots for assessing fitted GPD model. The user selects the plot type from a menu. See the examples below.
## S3 method for class 'gpd' plot(x, optlog = NA, extend = 1.5, labels = TRUE, ...)## S3 method for class 'gpd' plot(x, optlog = NA, extend = 1.5, labels = TRUE, ...)
x |
a |
optlog |
optional argument for plots 1 and 2 giving a
particular choice of logarithmic axes: |
extend |
optional argument for plots 1 and 2 expressing how far x-axis should extend as a multiple of the largest data value. This argument must take values greater than 1 and is useful for showing estimated quantiles beyond data. |
labels |
optional argument for plots 1 and 2 specifying whether or not axes should be labelled |
... |
other graphics parameters |
If plot 1 or 2 is selected as the final plot, a list object containing
details of the plot is returned invisibly.
This object should be used as the first argument of gpd.q
or gpd.sfall to add quantile estimates or expected shortfall
estimates to the plot.
data(danish) out <- gpd(danish, 10) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Excess Distribution ## Not run: 2: plot: Tail of Underlying Distribution ## Not run: 3: plot: Scatterplot of Residuals ## Not run: 4: plot: QQplot of Residualsdata(danish) out <- gpd(danish, 10) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Excess Distribution ## Not run: 2: plot: Tail of Underlying Distribution ## Not run: 3: plot: Scatterplot of Residuals ## Not run: 4: plot: QQplot of Residuals
Provides a number of plots summarising a bivariate GPD model fitted using the bivariate POT method. See the examples below.
## S3 method for class 'gpdbiv' plot(x, extend = 1.1, n.contours = 15, ...)## S3 method for class 'gpdbiv' plot(x, extend = 1.1, n.contours = 15, ...)
x |
a |
extend |
optional argument expressing how far x-axis should extend as a multiple of the largest data value. |
n.contours |
number of contours in bivariate contour plots |
... |
other graphics parameters |
Option 1 plots the threshold exceedance data; option 2 plots contours of the fitted bivariate distribution function in the joint upper tail (above both thresholds); option 3 plots corresponding contours of the fitted joint survival function; plots 4 and 5 show the fitted tails of the marginal distributions.
gpd, gpdbiv,
tailplot, interpret.gpdbiv,
plot.gpd
data(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Exceedance data ## Not run: 2: plot: Contours of Bivariate Distribution Function ## Not run: 3: plot: Contours of Bivariate Survival Function ## Not run: 4: plot: Tail of Marginal 1 ## Not run: 5: plot: Tail of Marginal 2data(bmw) ; data(siemens) out <- gpdbiv(-bmw, -siemens, ne1 = 100, ne2 = 100) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Exceedance data ## Not run: 2: plot: Contours of Bivariate Distribution Function ## Not run: 3: plot: Contours of Bivariate Survival Function ## Not run: 4: plot: Tail of Marginal 1 ## Not run: 5: plot: Tail of Marginal 2
The plot method plot.potd provides seven different plots for assessing fitted POT model. The user selects the plot type from a menu. See the examples below.
## S3 method for class 'potd' plot(x, ...)## S3 method for class 'potd' plot(x, ...)
x |
an object returned by the function |
... |
other graphics parameters |
Plot 1 displays the exceedance process of the chosen threshold. Plots 2-4 assess the Poisson nature of the exceedance process by looking at the scaled gaps between exceedances, which should be iid unit exponentially distributed. Plots 5-6 assess the GPD nature of the excesses by looking at suitably defined residuals, which should again be iid unit exponentially distributed. Option 8 allows the user to call GPD plotting functions.
If plot 1 or 2 from the GPD plots is selected as the final plot
(i.e. option 8 is selected, followed by option 1 or 2), a list
object containing details of the plot is returned invisibly.
This object should be used as the first argument of gpd.q
or gpd.sfall to add quantile estimates or expected shortfall
estimates to the plot.
data(danish) out <- pot(danish,10) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Point Process of Exceedances ## Not run: 2: plot: Scatterplot of Gaps ## Not run: 3: plot: Qplot of Gaps ## Not run: 4: plot: ACF of Gaps ## Not run: 5: plot: Scatterplot of Residuals ## Not run: 6: plot: Qplot of Residuals ## Not run: 7: plot: ACF of Residuals ## Not run: 8: plot: Go to GPD Plotsdata(danish) out <- pot(danish,10) ## Not run: plot(out) ## Not run: Make a plot selection (or 0 to exit): ## Not run: 1: plot: Point Process of Exceedances ## Not run: 2: plot: Scatterplot of Gaps ## Not run: 3: plot: Qplot of Gaps ## Not run: 4: plot: ACF of Gaps ## Not run: 5: plot: Scatterplot of Residuals ## Not run: 6: plot: Qplot of Residuals ## Not run: 7: plot: ACF of Residuals ## Not run: 8: plot: Go to GPD Plots
Fits a Poisson point process to the data, an approach
sometimes known as peaks over thresholds (POT), and
returns an object of class "potd".
pot(data, threshold = NA, nextremes = NA, run = NA, picture = TRUE, ...)pot(data, threshold = NA, nextremes = NA, run = NA, picture = TRUE, ...)
data |
numeric vector of data, which may have a |
threshold |
a threshold value (either this or |
nextremes |
the number of upper extremes to be used (either
this or |
run |
if the data are to be declustered the run length
parameter for the runs method (see |
picture |
whether or not a picture should be drawn if declustering is performed |
... |
arguments passed to |
Uses optim for point process likelihood maximization.
An object of class "potd" describing the fit and including
parameter estimates and standard errors.
gpd, plot.potd,
plot.gpd, decluster,
optim, as.POSIXct
data(danish) out <- pot(danish, 10) # Fits POT model to Danish fire insurance lossesdata(danish) out <- pot(danish, 10) # Fits POT model to Danish fire insurance losses
Creates a QQplot for threshold data against the exponential distribution or the generalized Pareto distribution.
qplot(data, xi = 0, trim = NA, threshold = NA, line = TRUE, labels = TRUE, ...)qplot(data, xi = 0, trim = NA, threshold = NA, line = TRUE, labels = TRUE, ...)
data |
data vector |
xi |
the xi value of a generalized Pareto distribution |
trim |
value at which data are to be right-truncated |
threshold |
value at which data are to be left-truncated |
line |
whether or not a straight line is to be added |
labels |
whether or not the axes are to be labelled |
... |
other graphics parameters |
If xi is zero the reference distribution is the exponential; if xi is non-zero the reference distribution is the generalized Pareto with that value of xi. In the case of the exponential, the plot is interpreted as follows. Concave departures from a straight line are a sign of heavy-tailed behaviour. Convex departures show thin-tailed behaviour.
## Not run: data(danish) ## Not run: qplot(danish) # QQplot of heavy-tailed Danish fire insurance data## Not run: data(danish) ## Not run: qplot(danish) # QQplot of heavy-tailed Danish fire insurance data
Creates a plot showing how the estimate of a high quantile in the tail of a dataset based on the GPD approximation varies with threshold or number of extremes.
quant(data, p = 0.99, models = 30, start = 15, end = 500, reverse = TRUE, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)quant(data, p = 0.99, models = 30, start = 15, end = 500, reverse = TRUE, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)
data |
numeric vector of data |
p |
desired probability for quantile estimate (e.g. 0.99 gives 99th percentile) |
models |
number of consecutive gpd models to be fitted |
start |
lowest number of exceedances to be considered |
end |
maximum number of exceedances to be considered |
reverse |
should plot be by increasing threshold
( |
ci |
probability for asymptotic confidence band; for no confidence band set to zero |
auto.scale |
whether or not plot should be automatically scaled; if not, xlim and ylim graphical parameters may be entered |
labels |
whether or not axes should be labelled |
... |
other graphics parameters |
For every model gpd is called. Evaluation may be slow.
Confidence intervals by the Wald method (which is fastest).
A table of results is returned invisibly.
## Not run: data(danish) ## Not run: quant(danish, 0.999) # Estimates of the 99.9th percentile of the Danish losses using # the GPD model with various thresholds## Not run: data(danish) ## Not run: quant(danish, 0.999) # Estimates of the 99.9th percentile of the Danish losses using # the GPD model with various thresholds
Creates a data frame showing the development of records in a dataset and calculating the expected behaviour for iid data.
records(data, do.plot = TRUE, conf.level = 0.95, ...)records(data, do.plot = TRUE, conf.level = 0.95, ...)
data |
data vector |
do.plot |
whether a plot of record development should be created |
conf.level |
confidence level for record development plot |
... |
graphics parameters |
Records are counted and the observations at which they occur recorded. This is compared with the expected behaviour for iid data.
A data frame.
## Not run: data(danish) ## Not run: records(danish) # Record fire insurance losses in Denmark## Not run: data(danish) ## Not run: records(danish) # Record fire insurance losses in Denmark
Makes a rapid calculation of point estimates of prescribed
quantiles and expected shortfalls using the output of the
function gpd.
riskmeasures(x, p)riskmeasures(x, p)
x |
results of a |
p |
a vector of probability levels |
This function simply calculates point estimates and (at present)
makes no attempt to calculate confidence intervals for the risk
measures. If confidence levels are required use gpd.q and
gpd.sfall which interact with graphs of the tail of a loss
distribution and are much slower.
A matrix with three columns: probability level, quantile estimate, shortfall estimate.
gpd, tailplot,
gpd.q, gpd.sfall
data(danish) out <- gpd(danish, 10) riskmeasures(out, c(0.999, 0.9999)) # gives estimates of 0.999 and 0.9999 quantiles of Danish loss # distribution as well as the associated expected shortfall estimatesdata(danish) out <- gpd(danish, 10) riskmeasures(out, c(0.999, 0.9999)) # gives estimates of 0.999 and 0.9999 quantiles of Danish loss # distribution as well as the associated expected shortfall estimates
Calculates the k-block return level and 95% confidence interval based on a GEV model for block maxima, where k is specified by the user. The k-block return level is that level exceeded once every k blocks, on average.
rlevel.gev(out, k.blocks = 20, add = FALSE, ...)rlevel.gev(out, k.blocks = 20, add = FALSE, ...)
out |
an object returned by the function |
k.blocks |
specifies the particular return level to be estimated; default set arbitrarily to 20 |
add |
whether the return level should be added graphically
to a time series plot; if |
... |
other graphics parameters |
The GEV likelihood is reparameterized in terms of the unknown return level and profile likelihood arguments are used to construct a confidence interval.
Vector containing lower 95% bound of confidence interval, estimated return level and upper 95% bound.
data(bmw) out <- gev(bmw, "month") # Fit GEV to monthly maxima of daily returns on BMW share price ## Not run: rlevel.gev(out, 40) # Calculate the 40 month return leveldata(bmw) out <- gev(bmw, "month") # Fit GEV to monthly maxima of daily returns on BMW share price ## Not run: rlevel.gev(out, 40) # Calculate the 40 month return level
Creates a plot showing how the estimate of shape varies with threshold or number of extremes.
shape(data, models = 30, start = 15, end = 500, reverse = TRUE, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)shape(data, models = 30, start = 15, end = 500, reverse = TRUE, ci = 0.95, auto.scale = TRUE, labels = TRUE, ...)
data |
numeric vector of data |
models |
number of consecutive gpd models to be fitted |
start |
lowest number of exceedances to be considered |
end |
maximum number of exceedances to be considered |
reverse |
should plot be by increasing threshold
( |
ci |
probability for asymptotic confidence band; for no confidence band set to zero |
auto.scale |
whether or not plot should be automatically
scaled; if not, |
labels |
whether or not axes should be labelled |
... |
other graphics parameters |
For every model gpd is called. Evaluation may be slow.
A table of results is returned invisibly.
## Not run: data(danish) ## Not run: shape(danish) # Shape plot of heavy-tailed Danish fire insurance data## Not run: data(danish) ## Not run: shape(danish) # Shape plot of heavy-tailed Danish fire insurance data
These data are the daily log returns on Siemens share price
from Tuesday 2nd January 1973 until Tuesday 23rd July 1996.
The data are contained in a numeric vector.
The dates of each observation are contained in a times
attribute, which is an object of class "POSIXct" (see
DateTimeClasses).
Note that these data form an irregular time series because
no trading takes place at the weekend.
data(siemens)data(siemens)
A numeric vector containing 6146 observations, with a
times attribute which is a POSIXct object
of the same length.
The daily closing values of the S&P index from
Monday 4th January 1960 until Friday 11th June 1993.
The data are contained in a numeric vector.
The dates of each observation are contained in a times
attribute, which is an object of class "POSIXct" (see
DateTimeClasses).
data(sp.raw)data(sp.raw)
A numeric vector containing 8415 observations, with a
times attribute which is a POSIXct object
of the same length.
The daily log returns on the S&P index value from
Tuesday 5th January 1960 until Friday 16 October 1987.
The data are contained in a numeric vector.
The dates of each observation are contained in a times
attribute, which is an object of class "POSIXct" (see
DateTimeClasses).
data(spto87)data(spto87)
A numeric vector containing 6985 observations, with a
times attribute which is a POSIXct object
of the same length.
Interacts with the output of gpd to produce a plot of the tail
of the underlying distribution of the data.
This is one of the options of plot.gpd, but tailplot enables
the user to bypass the menu of the former.
tailplot(x, optlog = NA, extend = 1.5, labels = TRUE, ...)tailplot(x, optlog = NA, extend = 1.5, labels = TRUE, ...)
x |
a |
optlog |
optional argument giving a particular choice
of logarithmic axes: |
extend |
optional argument expressing how far x-axis should extend as a multiple of the largest data value. This argument must take values greater than 1 and is useful for showing estimated quantiles beyond data. |
labels |
optional argument specifying whether or not axes should be labelled |
... |
other graphics parameters |
A list object containing details of the plot is returned invisibly.
This object should be used as the first argument of gpd.q
or gpd.sfall to add quantile estimates or expected shortfall
estimates to the plot.
data(danish) out <- gpd(danish, 10) ## Not run: tailplot(out)data(danish) out <- gpd(danish, 10) ## Not run: tailplot(out)