Version:	1.0.15
Date:	2024-11-30
Title:	Functions from "Reinsurance: Actuarial and Statistical Aspects"
Description:	Functions from the book "Reinsurance: Actuarial and Statistical Aspects" (2017) by Hansjoerg Albrecher, Jan Beirlant and Jef Teugels https://www.wiley.com/en-us/Reinsurance%3A+Actuarial+and+Statistical+Aspects-p-9780470772683.
Maintainer:	Tom Reynkens <tomreynkens.r@gmail.com>
Depends:	R (≥ 3.4)
Imports:	stats, graphics, methods, survival, utils, parallel, foreach (≥ 1.4.1), doParallel (≥ 1.0.10), Rcpp (≥ 0.12.12)
Suggests:	testthat, knitr, rmarkdown, interval, Icens
LinkingTo:	Rcpp
License:	GPL-2 | GPL-3 [expanded from: GPL (≥ 2)]
URL:	https://github.com/TReynkens/ReIns
BugReports:	https://github.com/TReynkens/ReIns/issues
VignetteBuilder:	knitr
Encoding:	UTF-8
NeedsCompilation:	yes
ByteCompile:	yes
Packaged:	2024-11-30 13:02:34 UTC; tom
Author:	Tom Reynkens [aut, cre], Roel Verbelen [aut] (R code for Mixed Erlang distribution), Anastasios Bardoutsos [ctb] (Original R code for cEPD estimator), Dries Cornilly [ctb] (Original R code for EVT estimators for truncated data), Yuri Goegebeur [ctb] (Original S-Plus code for basic EVT estimators), Klaus Herrmann [ctb] (Original R code for GPD estimator)
Repository:	CRAN
Date/Publication:	2024-12-02 09:30:05 UTC

The Burr distribution

Description

Density, distribution function, quantile function and random generation for the Burr distribution (type XII).

Usage

dburr(x, alpha, rho, eta = 1, log = FALSE)
pburr(x, alpha, rho, eta = 1, lower.tail = TRUE, log.p = FALSE)
qburr(p, alpha, rho, eta = 1, lower.tail = TRUE, log.p = FALSE)
rburr(n, alpha, rho, eta = 1)

Arguments

Details

The Cumulative Distribution Function (CDF) of the Burr distribution is equal to F(x) = 1-((\eta+x^{-\rho\times\alpha})/\eta)^{1/\rho} for all x \ge 0 and F(x)=0 otherwise. We need that \alpha>0, \rho<0 and \eta>0.

Beirlant et al. (2004) uses parameters \eta, \tau, \lambda which correspond to \eta, \tau=-\rho\times\alpha and \lambda=-1/\rho.

Value

dburr gives the density function evaluated in x, pburr the CDF evaluated in x and qburr the quantile function evaluated in p. The length of the result is equal to the length of x or p.

rburr returns a random sample of length n.

Author(s)

Tom Reynkens.

References

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

tBurr, Distributions

Examples

# Plot of the PDF
x <- seq(0, 10, 0.01)
plot(x, dburr(x, alpha=2, rho=-1), xlab="x", ylab="PDF", type="l")

# Plot of the CDF
x <- seq(0, 10, 0.01)
plot(x, pburr(x, alpha=2, rho=-1), xlab="x", ylab="CDF", type="l")

Conditional Tail Expectation

Description

Compute Conditional Tail Expectation (CTE) CTE_{1-p} of the fitted spliced distribution.

Usage

CTE(p, splicefit)

ES(p, splicefit)

Arguments

Details

The Conditional Tail Expectation is defined as

CTE_{1-p} = E(X | X>Q(1-p)) = E(X | X>VaR_{1-p}) = VaR_{1-p} + \Pi(VaR_{1-p})/p,

where \Pi(u)=E((X-u)_+) is the premium of the excess-loss insurance with retention u.

If the CDF is continuous in p, we have CTE_{1-p}=TVaR_{1-p}= 1/p \int_0^p VaR_{1-s} ds with TVaR the Tail Value-at-Risk.

See Reynkens et al. (2017) and Section 4.6 of Albrecher et al. (2017) for more details.

The ES function is the same function as CTE but is deprecated.

Value

Vector with the CTE corresponding to each element of p.

Author(s)

Tom Reynkens with R code from Roel Verbelen for the mixed Erlang quantiles.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

qSplice, ExcessSplice, SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)

p <- seq(0.01, 0.99, 0.01)
# Plot of CTE
plot(p, CTE(p, splicefit), type="l", xlab="p", ylab=bquote(CTE[1-p]))

## End(Not run)

EPD estimator

Description

Fit the Extended Pareto Distribution (GPD) to the exceedances (peaks) over a threshold. Optionally, these estimates are plotted as a function of k.

Usage

EPD(data, rho = -1, start = NULL, direct = FALSE, warnings = FALSE, 
    logk = FALSE, plot = FALSE, add = FALSE, main = "EPD estimates of the EVI", ...)

Arguments

Details

We fit the Extended Pareto distribution to the relative excesses over a threshold (X/u). The EPD has distribution function F(x) = 1-(x(1+\kappa-\kappa x^{\tau}))^{-1/\gamma} with \tau = \rho/\gamma <0<\gamma and \kappa>\max(-1,1/\tau).

The parameters are determined using MLE and there are two possible approaches: maximise the log-likelihood directly (direct=TRUE) or follow the approach detailed in Beirlant et al. (2009) (direct=FALSE). The latter approach uses the score functions of the log-likelihood.

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Joossens, E. and Segers, J. (2009). "Second-Order Refined Peaks-Over-Threshold Modelling for Heavy-Tailed Distributions." Journal of Statistical Planning and Inference, 139, 2800–2815.

Fraga Alves, M.I. , Gomes, M.I. and de Haan, L. (2003). "A New Class of Semi-parametric Estimators of the Second Order Parameter." Portugaliae Mathematica, 60, 193–214.

See Also

GPDmle, ProbEPD

Examples

data(secura)

# EPD estimates for the EVI
epd <- EPD(secura$size, plot=TRUE)

# Compute return periods
ReturnEPD(secura$size, 10^10, gamma=epd$gamma, kappa=epd$kappa, 
          tau=epd$tau, plot=TRUE)

Fit EPD using MLE

Description

Fit the Extended Pareto Distribution (EPD) to data using Maximum Likelihood Estimation (MLE).

Usage

EPDfit(data, tau, start = c(0.1, 1), warnings = FALSE)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A vector with the MLE estimate for the \gamma parameter of the EPD as the first component and the MLE estimate for the \kappa parameter of the EPD as the second component.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Joossens, E. and Segers, J. (2009). "Second-Order Refined Peaks-Over-Threshold Modelling for Heavy-Tailed Distributions." Journal of Statistical Planning and Inference, 139, 2800–2815.

See Also

EPD, GPDfit

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Fit EPD to last 500 observations
res <- EPDfit(SOAdata/sort(soa$size)[500], tau=-1)

EVT fit

Description

Create an S3 object using an EVT (Extreme Value Theory) fit.

Usage

EVTfit(gamma, endpoint = NULL, sigma = NULL)

Arguments

Details

See Reynkens et al. (2017) and Section 4.3 of Albrecher et al. (2017) for details.

Value

An S3 object containing the above input arguments.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

See Also

SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

# Create MEfit object
mefit <- MEfit(p=c(0.65,0.35), shape=c(39,58), theta=16.19, M=2)

# Create EVTfit object
evtfit <- EVTfit(gamma=c(0.76,0.64), endpoint=c(39096, Inf))

# Create SpliceFit object
splicefit <- SpliceFit(const=c(0.5,0.996), trunclower=0, t=c(1020,39096), type=c("ME","TPa","Pa"),
                       MEfit=mefit, EVTfit=evtfit)

# Show summary
summary(splicefit)

Estimates for excess-loss premiums using EPD estimates

Description

Estimate premiums of excess-loss reinsurance with retention R and limit L using EPD estimates.

Usage

ExcessEPD(data, gamma, kappa, tau, R, L = Inf, warnings = TRUE, plot = TRUE, add = FALSE, 
          main = "Estimates for premium of excess-loss insurance", ...)

Arguments

Details

We need that u \ge X_{n-k,n}, the (k+1)-th largest observation. If this is not the case, we return NA for the premium. A warning will be issued in that case if warnings=TRUE.

The premium for the excess-loss insurance with retention R and limit L is given by

E(\min{(X-R)_+, L}) = \Pi(R) - \Pi(R+L)

where \Pi(u)=E((X-u)_+)=\int_u^{\infty} (1-F(z)) dz is the premium of the excess-loss insurance with retention u. When L=\infty, the premium is equal to \Pi(R).

We estimate \Pi by

\hat{\Pi}(u) = (k+1)/(n+1) \times (X_{n-k,n})^{1/\hat{\gamma}} \times ((1-\hat{\kappa}/\hat{\gamma})(1/\hat{\gamma}-1)^{-1}u^{1-1/\hat{\gamma}} + \hat{\kappa}/(\hat{\gamma}X_{n-k,n}^{\hat{\tau}})(1/\hat{\gamma}-\hat{\tau}-1)^{-1}u^{1+\hat{\tau}-1/\hat{\gamma}})

with \hat{\gamma}, \hat{\kappa} and \hat{\tau} the estimates for the parameters of the EPD.

See Section 4.6 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

EPD, ExcessHill, ExcessGPD

Examples

data(secura)

# EPD estimator
epd <- EPD(secura$size)

# Premium of excess-loss insurance with retention R
R <- 10^7
ExcessEPD(secura$size, gamma=epd$gamma, kappa=epd$kappa, tau=epd$tau, R=R, ylim=c(0,2*10^4))

Estimates for excess-loss premiums using GPD-MLE estimates

Description

Estimate premiums of excess-loss reinsurance with retention R and limit L using GPD-MLE estimates.

Usage

ExcessGPD(data, gamma, sigma, R, L = Inf, warnings = TRUE, plot = TRUE, add = FALSE, 
          main = "Estimates for premium of excess-loss insurance", ...)

Arguments

Details

We need that u \ge X_{n-k,n}, the (k+1)-th largest observation. If this is not the case, we return NA for the premium. A warning will be issued in that case if warnings=TRUE. One should then use global fits: ExcessSplice.

The premium for the excess-loss insurance with retention R and limit L is given by

E(\min{(X-R)_+, L}) = \Pi(R) - \Pi(R+L)

where \Pi(u)=E((X-u)_+)=\int_u^{\infty} (1-F(z)) dz is the premium of the excess-loss insurance with retention u. When L=\infty, the premium is equal to \Pi(R).

We estimate \Pi by

\hat{\Pi}(u) = (k+1)/(n+1) \times \hat{\sigma}_k/ (1-\hat{\gamma}_k) \times (1+\hat{\gamma}_k/\hat{\sigma}_k (u-X_{n-k,n}))^{1-1/\hat{\gamma}_k},

with \hat{\gamma}_k and \hat{\sigma}_k the estimates for the parameters of the GPD.

See Section 4.6 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

GPDmle, ExcessHill, ExcessEPD

Examples

data(secura)

# GPDmle estimator
mle <- GPDmle(secura$size)

# Premium of excess-loss insurance with retention R
R <- 10^7
ExcessGPD(secura$size, gamma=mle$gamma, sigma=mle$sigma, R=R, ylim=c(0,2*10^4))

Estimates for excess-loss premiums using a Pareto model

Description

Estimate premiums of excess-loss reinsurance with retention R and limit L using a (truncated) Pareto model.

Usage

ExcessPareto(data, gamma, R, L = Inf, endpoint = Inf, warnings = TRUE, plot = TRUE, 
        add = FALSE, main = "Estimates for premium of excess-loss insurance", ...)
        
ExcessHill(data, gamma, R, L = Inf, endpoint = Inf, warnings = TRUE, plot = TRUE, 
        add = FALSE, main = "Estimates for premium of excess-loss insurance", ...)

Arguments

Details

We need that u \ge X_{n-k,n}, the (k+1)-th largest observation. If this is not the case, we return NA for the premium. A warning will be issued in that case if warnings=TRUE. One should then use global fits: ExcessSplice.

The premium for the excess-loss insurance with retention R and limit L is given by

E(\min{(X-R)_+, L}) = \Pi(R) - \Pi(R+L)

where \Pi(u)=E((X-u)_+)=\int_u^{\infty} (1-F(z)) dz is the premium of the excess-loss insurance with retention u. When L=\infty, the premium is equal to \Pi(R).

We estimate \Pi (for the untruncated Pareto distribution) by

\hat{\Pi}(u) = (k+1)/(n+1) / (1/H_{k,n}-1) \times (X_{n-k,n}^{1/H_{k,n}} u^{1-1/H_{k,n}}),

with H_{k,n} the Hill estimator.

The ExcessHill function is the same function but with a different name for compatibility with old versions of the package.

See Section 4.6 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

Hill, ExcessEPD, ExcessGPD, ExcessSplice

Examples

data(secura)

# Hill estimator
H <- Hill(secura$size)

# Premium of excess-loss insurance with retention R
R <- 10^7
ExcessPareto(secura$size, H$gamma, R=R)

Estimates for excess-loss premiums using splicing

Description

Estimate premiums of excess-loss reinsurance with retention R and limit L using fitted spliced distribution.

Usage

ExcessSplice(R, L=Inf, splicefit)

Arguments

Details

The premium for the excess-loss insurance with retention R and limit L is given by

E(\min{(X-R)_+, L}) = \Pi(R) - \Pi(R+L)

where \Pi(u)=E((X-u)_+)=\int_u^{\infty} (1-F(z)) dz is the premium of the excess-loss insurance with retention u. When L=\infty, the premium is equal to \Pi(R).

See Reynkens et al. (2017) and Section 4.6 of Albrecher et al. (2017) for more details.

Value

An estimate for the premium is returned (for every value of R).

Author(s)

Tom Reynkens with R code from Roel Verbelen for the estimates for the excess-loss premiums using the mixed Erlang distribution.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.8)

# Excess-loss premium 
ExcessSplice(R=2, splicefit=splicefit)

## End(Not run)

Exponential quantile plot

Description

Computes the empirical quantiles of a data vector and the theoretical quantiles of the standard exponential distribution. These quantiles are then plotted in an exponential QQ-plot with the theoretical quantiles on the x-axis and the empirical quantiles on the y-axis.

Usage

ExpQQ(data, plot = TRUE, main = "Exponential QQ-plot", ...)

Arguments

Details

The exponential QQ-plot is defined as

( -\log(1-i/(n+1)), X_{i,n} )

for i=1,...,n, with X_{i,n} the i-th order statistic of the data.

Note that the mean excess plot is the derivative plot of the Exponential QQ-plot.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

MeanExcess, LognormalQQ, ParetoQQ, WeibullQQ

Examples

data(norwegianfire)

# Exponential QQ-plot for Norwegian Fire Insurance data for claims in 1976.
ExpQQ(norwegianfire$size[norwegianfire$year==76])

# Pareto QQ-plot for Norwegian Fire Insurance data for claims in 1976.
ParetoQQ(norwegianfire$size[norwegianfire$year==76])

The Extended Pareto Distribution

Description

Density, distribution function, quantile function and random generation for the Extended Pareto Distribution (EPD).

Usage

depd(x, gamma, kappa, tau = -1, log = FALSE)
pepd(x, gamma, kappa, tau = -1, lower.tail = TRUE, log.p = FALSE)
qepd(p, gamma, kappa, tau = -1, lower.tail = TRUE, log.p = FALSE)
repd(n, gamma, kappa, tau = -1)

Arguments

Details

The Cumulative Distribution Function (CDF) of the EPD is equal to F(x) = 1-(x(1+\kappa-\kappa x^{\tau}))^{-1/\gamma} for all x > 1 and F(x)=0 otherwise.

Note that an EPD random variable with \tau=-1 and \kappa=\gamma/\sigma-1 is GPD distributed with \mu=1, \gamma and \sigma.

Value

depd gives the density function evaluated in x, pepd the CDF evaluated in x and qepd the quantile function evaluated in p. The length of the result is equal to the length of x or p.

repd returns a random sample of length n.

Author(s)

Tom Reynkens.

References

Beirlant, J., Joossens, E. and Segers, J. (2009). "Second-Order Refined Peaks-Over-Threshold Modelling for Heavy-Tailed Distributions." Journal of Statistical Planning and Inference, 139, 2800–2815.

See Also

Pareto, GPD, Distributions

Examples

# Plot of the PDF
x <- seq(0, 10, 0.01)
plot(x, depd(x, gamma=1/2, kappa=1, tau=-1), xlab="x", ylab="PDF", type="l")

# Plot of the CDF
x <- seq(0, 10, 0.01)
plot(x, pepd(x, gamma=1/2, kappa=1, tau=-1), xlab="x", ylab="CDF", type="l")

The Frechet distribution

Description

Density, distribution function, quantile function and random generation for the Fréchet distribution (inverse Weibull distribution).

Usage

dfrechet(x, shape, loc = 0, scale = 1, log = FALSE)
pfrechet(x, shape, loc = 0, scale = 1, lower.tail = TRUE, log.p = FALSE)
qfrechet(p, shape, loc = 0, scale = 1, lower.tail = TRUE, log.p = FALSE)
rfrechet(n, shape, loc = 0, scale = 1)

Arguments

Details

The Cumulative Distribution Function (CDF) of the Fréchet distribution is equal to F(x) = \exp(-((x-loc)/scale)^{-shape}) for all x \ge loc and F(x)=0 otherwise. Both shape and scale need to be strictly positive.

Value

dfrechet gives the density function evaluated in x, pfrechet the CDF evaluated in x and qfrechet the quantile function evaluated in p. The length of the result is equal to the length of x or p.

rfrechet returns a random sample of length n.

Author(s)

Tom Reynkens.

See Also

tFréchet, Distributions

Examples

# Plot of the PDF
x <- seq(1,10,0.01)
plot(x, dfrechet(x, shape=2), xlab="x", ylab="PDF", type="l")

# Plot of the CDF
x <- seq(1,10,0.01)
plot(x, pfrechet(x, shape=2), xlab="x", ylab="CDF", type="l")

The generalised Pareto distribution

Description

Density, distribution function, quantile function and random generation for the Generalised Pareto Distribution (GPD).

Usage

dgpd(x, gamma, mu = 0, sigma, log = FALSE)
pgpd(x, gamma, mu = 0, sigma, lower.tail = TRUE, log.p = FALSE)
qgpd(p, gamma, mu = 0, sigma, lower.tail = TRUE, log.p = FALSE)
rgpd(n, gamma, mu = 0, sigma)

Arguments

Details

The Cumulative Distribution Function (CDF) of the GPD for \gamma \neq 0 is equal to F(x) = 1-(1+\gamma (x-\mu)/\sigma)^{-1/\gamma} for all x \ge \mu and F(x)=0 otherwise. When \gamma=0, the CDF is given by F(x) = 1-\exp((x-\mu)/\sigma) for all x \ge \mu and F(x)=0 otherwise.

Value

dgpd gives the density function evaluated in x, pgpd the CDF evaluated in x and qgpd the quantile function evaluated in p. The length of the result is equal to the length of x or p.

rgpd returns a random sample of length n.

Author(s)

Tom Reynkens.

References

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

tGPD, Pareto, EPD, Distributions

Examples

# Plot of the PDF
x <- seq(0, 10, 0.01)
plot(x, dgpd(x, gamma=1/2, sigma=5), xlab="x", ylab="PDF", type="l")

# Plot of the CDF
x <- seq(0, 10, 0.01)
plot(x, pgpd(x, gamma=1/2, sigma=5), xlab="x", ylab="CDF", type="l")

Fit GPD using MLE

Description

Fit the Generalised Pareto Distribution (GPD) to data using Maximum Likelihood Estimation (MLE).

Usage

GPDfit(data, start = c(0.1, 1), warnings = FALSE)

Arguments

Details

See Section 4.2.2 in Albrecher et al. (2017) for more details.

Value

A vector with the MLE estimate for the \gamma parameter of the GPD as the first component and the MLE estimate for the \sigma parameter of the GPD as the second component.

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur and R code from Klaus Herrmann.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

GPDmle, EPDfit

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Fit GPD to last 500 observations
res <- GPDfit(SOAdata-sort(soa$size)[500])

GPD-ML estimator

Description

Fit the Generalised Pareto Distribution (GPD) to the exceedances (peaks) over a threshold using Maximum Likelihood Estimation (MLE). Optionally, these estimates are plotted as a function of k.

Usage

GPDmle(data, start = c(0.1,1), warnings = FALSE, logk = FALSE, 
       plot = FALSE, add = FALSE, main = "POT estimates of the EVI", ...)

POT(data, start = c(0.1,1), warnings = FALSE, logk = FALSE, 
    plot = FALSE, add = FALSE, main = "POT estimates of the EVI", ...)

Arguments

Details

The POT function is the same function but with a different name for compatibility with the old S-Plus code.

For each value of k, we look at the exceedances over the (k+1)th largest observation: X_{n-k+j,n}-X_{n-k,n} for j=1,...,k, with X_{j,n} the jth largest observation and n the sample size. The GPD is then fitted to these k exceedances using MLE which yields estimates for the parameters of the GPD: \gamma and \sigma.

See Section 4.2.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur and R code from Klaus Herrmann.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

GPDfit, GPDresiduals, EPD

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Plot GPD-ML estimates as a function of k
GPDmle(SOAdata, plot=TRUE)

GPD residual plot

Description

Residual plot to check GPD fit for peaks over a threshold.

Usage

GPDresiduals(data, t, gamma, sigma, plot = TRUE, 
             main = "GPD residual plot", ...)

Arguments

Details

Consider the POT values Y=X-t and the transformed variable

R= 1/\gamma \log(1+\gamma/\sigma Y),

when \gamma \neq 0 and

R = Y/\sigma,

otherwise. We can assess the goodness-of-fit of the GPD when modelling POT values Y=X-t by constructing an exponential QQ-plot of the transformed variable R since R is standard exponentially distributed if Y follows the GPD.

See Section 4.2.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

GPDfit, ExpQQ

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Plot POT-MLE estimates as a function of k
pot <- GPDmle(SOAdata, plot=TRUE)

# Residual plot
k <- 200
GPDresiduals(SOAdata, sort(SOAdata)[length(SOAdata)-k], pot$gamma[k], pot$sigma[k])

Hill estimator

Description

Computes the Hill estimator for positive extreme value indices (Hill, 1975) as a function of the tail parameter k. Optionally, these estimates are plotted as a function of k.

Usage

Hill(data, k = TRUE, logk = FALSE, plot = FALSE, add = FALSE, 
     main = "Hill estimates of the EVI", ...)

Arguments

Details

The Hill estimator can be seen as the estimator of slope in the upper right corner (k last points) of the Pareto QQ-plot when using constrained least squares (the regression line has to pass through the point (-\log((k+1)/(n+1)),\log X_{n-k})). It is given by

H_{k,n}=1/k\sum_{j=1}^k \log X_{n-j+1,n}- \log X_{n-k,n}.

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Hill, B. M. (1975). "A simple general approach to inference about the tail of a distribution." Annals of Statistics, 3, 1163–1173.

See Also

ParetoQQ, Hill.2oQV, genHill

Examples

data(norwegianfire)

# Plot Hill estimates as a function of k
Hill(norwegianfire$size[norwegianfire$year==76],plot=TRUE)

Bias-reduced MLE (Quantile view)

Description

Computes bias-reduced ML estimates of gamma based on the quantile view.

Usage

Hill.2oQV(data, start = c(1,1,1), warnings = FALSE, logk = FALSE, 
          plot = FALSE, add = FALSE, main = "Estimates of the EVI", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur and R code from Klaus Herrmann.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Dierckx, G., Goegebeur Y. and Matthys, G. (1999). "Tail Index Estimation and an Exponential Regression Model." Extremes, 2, 177–200.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Examples

data(norwegianfire)

# Plot bias-reduced MLE (QV) as a function of k
Hill.2oQV(norwegianfire$size[norwegianfire$year==76],plot=TRUE)

Select optimal threshold for Hill estimator

Description

Select optimal threshold for the Hill estimator by minimising the Asymptotic Mean Squared Error (AMSE).

Usage

Hill.kopt(data, start = c(1, 1, 1), warnings = FALSE, logk = FALSE, 
          plot = FALSE, add = FALSE, main = "AMSE plot", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Beirlant J., Vynckier, P. and Teugels, J. (1996). "Tail Index Estimation, Pareto Quantile Plots, and Regression Diagnostics." Journal of the American Statistical Association, 91, 1659–1667.

See Also

Hill, Hill.2oQV

Examples

data(norwegianfire)

# Plot Hill estimator as a function of k
Hill(norwegianfire$size[norwegianfire$year==76],plot=TRUE)

# Add optimal value of k
Hill.kopt(norwegianfire$size[norwegianfire$year==76],add=TRUE)

Kaplan-Meier estimator

Description

Computes the Kaplan-Meier estimator for the survival function of right censored data.

Usage

KaplanMeier(x, data, censored, conf.type="plain", conf.int = 0.95)

Arguments

Details

We consider the random right censoring model where one observes Z = \min(X,C) where X is the variable of interest and C is the censoring variable.

This function is merely a wrapper for survfit.formula from survival.

This estimator is only suitable for right censored data. When the data are interval censored, one can use the Turnbull estimator implemented in Turnbull.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Kaplan, E. L. and Meier, P. (1958). "Nonparametric Estimation from Incomplete Observations." Journal of the American Statistical Association, 53, 457–481.

See Also

survfit.formula, Turnbull

Examples

data <- c(1, 2.5, 3, 4, 5.5, 6, 7.5, 8.25, 9, 10.5)
censored <- c(0, 1, 0, 0, 1, 0, 1, 1, 0, 1)

x <- seq(0, 12, 0.1)

# Kaplan-Meier estimator
plot(x, KaplanMeier(x, data, censored)$surv, type="s", ylab="Kaplan-Meier estimator")

Least Squares tail estimator

Description

Computes the Least Squares (LS) estimates of the EVI based on the last k observations of the generalised QQ-plot.

Usage

LStail(data, rho = -1, lambda = 0.5, logk = FALSE, plot = FALSE, add = FALSE, 
       main = "LS estimates of the EVI", ...)
            
TSfraction(data, rho = -1, lambda = 0.5, logk = FALSE, plot = FALSE, add = FALSE, 
           main = "LS estimates of the EVI", ...)

Arguments

Details

We estimate \gamma (EVI) and b using least squares on the following regression model (Beirlant et al., 2005): Z_j = \gamma + b(n/k) (j/k)^{-\rho} + \epsilon_j with Z_j = (j+1) \log(UH_{j,n}/UH_{j+1,n}) and UH_{j,n}=X_{n-j,n}H_{j,n}, where H_{j,n} is the Hill estimator with threshold X_{n-j,n}.

See Section 5.8 of Beirlant et al. (2004) for more details.

The function TSfraction is included for compatibility with the old S-Plus code.

Value

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Beirlant, J., Dierckx, G. and Guillou, A. (2005). "Estimation of the Extreme Value Index and Regression on Generalized Quantile Plots." Bernoulli, 11, 949–970.

Beirlant, J., Dierckx, G., Guillou, A. and Starica, C. (2002). "On Exponential Representations of Log-spacing of Extreme Order Statistics." Extremes, 5, 157–180.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

genQQ

Examples

data(soa)

# LS tail estimator
LStail(soa$size, plot=TRUE, ylim=c(0,0.5))

Log-normal quantile plot

Description

Computes the empirical quantiles of the log-transform of a data vector and the theoretical quantiles of the standard normal distribution. These quantiles are then plotted in a log-normal QQ-plot with the theoretical quantiles on the x-axis and the empirical quantiles on the y-axis.

Usage

LognormalQQ(data, plot = TRUE, main = "Log-normal QQ-plot", ...)

Arguments

Details

By definition, a log-transformed log-normal random variable is normally distributed. We can thus obtain a log-normal QQ-plot from a normal QQ-plot by replacing the empirical quantiles of the data vector by the empirical quantiles from the log-transformed data. We hence plot

(\Phi^{-1}(i/(n+1)), \log(X_{i,n}) )

for i=1,\ldots,n, where \Phi is the standard normal CDF.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

ExpQQ, ParetoQQ, WeibullQQ

Examples

data(norwegianfire)

# Log-normal QQ-plot for Norwegian Fire Insurance data for claims in 1976.
LognormalQQ(norwegianfire$size[norwegianfire$year==76])

Derivative plot of the log-normal QQ-plot

Description

Computes the derivative plot of the log-normal QQ-plot. These values can be plotted as a function of the data or as a function of the tail parameter k.

Usage

LognormalQQ_der(data, k = FALSE, plot = TRUE, 
                main = "Derivative plot of log-normal QQ-plot", ...)

Arguments

Details

The derivative plot of a log-normal QQ-plot is

(k, H_{k,n}/N_{k,n})

or

(\log X_{n-k,n}, H_{k,n}/N_{k,n})

with H_{k,n} the Hill estimates and

N_{k,n} = (n+1)/(k+1) \phi(\Phi^{-1}(a)) - \Phi^{-1}(a).

Here is a=1-(k+1)/(n+1), \phi the standard normal PDF and \Phi the standard normal CDF.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

LognormalQQ, Hill, MeanExcess, ParetoQQ_der, WeibullQQ_der

Examples

data(norwegianfire)

# Log-normal QQ-plot for Norwegian Fire Insurance data for claims in 1976.
LognormalQQ(norwegianfire$size[norwegianfire$year==76])

# Derivate plot
LognormalQQ_der(norwegianfire$size[norwegianfire$year==76])

Mixed Erlang fit

Description

Create an S3 object using a Mixed Erlang (ME) fit.

Usage

MEfit(p, shape, theta, M, M_initial = NULL)

Arguments

Details

The rate parameter \lambda used in Albrecher et al. (2017) is equal to 1/\theta.

See Reynkens et al. (2017) and Section 4.3 of Albrecher et al. (2017) for more details

Value

An S3 object which contains the input arguments in a list.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

# Create MEfit object
mefit <- MEfit(p=c(0.65,0.35), shape=c(39,58), theta=16.19, M=2)

# Create EVTfit object
evtfit <- EVTfit(gamma=c(0.76,0.64), endpoint=c(39096, Inf))

# Create SpliceFit object
splicefit <- SpliceFit(const=c(0.5,0.996), trunclower=0, t=c(1020,39096), type=c("ME","TPa","Pa"),
                       MEfit=mefit, EVTfit=evtfit)

# Show summary
summary(splicefit)

Mean excess function

Description

Computes the mean excess values for a vector of observations. These mean excess values can then be plotted as a function of the data or as a function of the tail parameter k.

Usage

MeanExcess(data, plot = TRUE, k = FALSE, main = "Mean excess plot", ...)

Arguments

Details

The mean excess plot is

(k,e_{k,n})

or

(X_{n-k,n}, e_{k,n})

with

e_{k,n}=1/k\sum_{j=1}^k X_{n-j+1,n}-X_{n-k,n}.

Note that the mean excess plot is the derivative plot of the Exponential QQ-plot.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

ExpQQ, LognormalQQ_der, ParetoQQ_der, WeibullQQ_der

Examples

data(norwegianfire)

# Mean excess plots for Norwegian Fire Insurance data for claims in 1976.

# Mean excess values as a function of k
MeanExcess(norwegianfire$size[norwegianfire$year==76], k=TRUE)

# Mean excess values as a function of the data
MeanExcess(norwegianfire$size[norwegianfire$year==76], k=FALSE)

Mean excess function using Turnbull estimator

Description

Computes mean excess values using the Turnbull estimator. These mean excess values can then be plotted as a function of the empirical quantiles (computed using the Turnbull estimator) or as a function of the tail parameter k.

Usage

MeanExcess_TB(L, U = L, censored, trunclower = 0, truncupper = Inf, 
              plot = TRUE, k = FALSE, intervalpkg = TRUE, 
              main = "Mean excess plot", ...)

Arguments

Details

The mean excess values are given by

\hat{e}^{TB}(v)=(\int_v^{\infty} 1-\hat{F}^{TB}(u) du)/(1-\hat{F}^{TB}(v))

where \hat{F}^{TB} is the Turnbull estimator for the CDF. More specifically, we use the values v=\hat{Q}^{TB}(p) for p=1/(n+1), \ldots, (n-1)/(n+1) where \hat{Q}^{TB}(p) is the empirical quantile function corresponding to the Turnbull estimator.

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u.

If the interval package is installed and intervalpkg=TRUE, the icfit function is used to compute the Turnbull estimator. Otherwise, survfit.formula from survival is used.

Use MeanExcess for non-censored data.

See Section 4.3 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

MeanExcess, Turnbull, icfit

Examples

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Mean excess plot
MeanExcess_TB(Z, U, censored, k=FALSE)

Moment estimator

Description

Compute the moment estimates for real extreme value indices as a function of the tail parameter k. Optionally, these estimates are plotted as a function of k.

Usage

Moment(data, logk = FALSE, plot = FALSE, add = FALSE, 
       main = "Moment estimates of the EVI", ...)

Arguments

Details

The moment estimator for the EVI is introduced by Dekkers et al. (1989) and is a generalisation of the Hill estimator.

See Section 4.2.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Dekkers, A.L.M, Einmahl, J.H.J. and de Haan, L. (1989). "A Moment Estimator for the Index of an Extreme-value Distribution." Annals of Statistics, 17, 1833–1855.

See Also

Hill, genHill

Examples

data(soa)

# Hill estimator
H <- Hill(soa$size, plot=FALSE)
# Moment estimator
M <- Moment(soa$size)
# Generalised Hill estimator
gH <- genHill(soa$size, gamma=H$gamma)

# Plot estimates
plot(H$k[1:5000], M$gamma[1:5000], xlab="k", ylab=expression(gamma), type="l", ylim=c(0.2,0.5))
lines(H$k[1:5000], gH$gamma[1:5000], lty=2)
legend("topright", c("Moment", "Generalised Hill"), lty=1:2)

The Pareto distribution

Description

Density, distribution function, quantile function and random generation for the Pareto distribution (type I).

Usage

dpareto(x, shape, scale = 1, log = FALSE)
ppareto(x, shape, scale = 1, lower.tail = TRUE, log.p = FALSE)
qpareto(p, shape, scale = 1, lower.tail = TRUE, log.p = FALSE)
rpareto(n, shape, scale = 1)

Arguments

Details

The Cumulative Distribution Function (CDF) of the Pareto distribution is equal to F(x) = 1-(x/scale)^{-shape} for all x \ge scale and F(x)=0 otherwise. Both shape and scale need to be strictly positive.

Value

dpareto gives the density function evaluated in x, ppareto the CDF evaluated in x and qpareto the quantile function evaluated in p. The length of the result is equal to the length of x or p.

rpareto returns a random sample of length n.

Author(s)

Tom Reynkens.

See Also

tPareto, GPD, Distributions

Examples

# Plot of the PDF
x <- seq(1, 10, 0.01)
plot(x, dpareto(x, shape=2), xlab="x", ylab="PDF", type="l")

# Plot of the CDF
x <- seq(1, 10, 0.01)
plot(x, ppareto(x, shape=2), xlab="x", ylab="CDF", type="l")

Pareto quantile plot

Description

Computes the empirical quantiles of the log-transform of a data vector and the theoretical quantiles of the standard exponential distribution. These quantiles are then plotted in a Pareto QQ-plot with the theoretical quantiles on the x-axis and the empirical quantiles on the y-axis.

Usage

ParetoQQ(data, plot = TRUE, main = "Pareto QQ-plot", ...)

Arguments

Details

It can be easily seen that a log-transformed Pareto random variable is exponentially distributed. We can hence obtain a Pareto QQ-plot from an exponential QQ-plot by replacing the empirical quantiles from the data vector by the empirical quantiles from the log-transformed data. We hence plot

( -\log(1-i/(n+1)), \log X_{i,n} )

for i=1,...,n, with X_{i,n} the i-th order statistic of the data.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

ParetoQQ_der, ExpQQ, genQQ, LognormalQQ, WeibullQQ

Examples

data(norwegianfire)

# Exponential QQ-plot for Norwegian Fire Insurance data for claims in 1976.
ExpQQ(norwegianfire$size[norwegianfire$year==76])

# Pareto QQ-plot for Norwegian Fire Insurance data for claims in 1976.
ParetoQQ(norwegianfire$size[norwegianfire$year==76])

Derivative plot of the Pareto QQ-plot

Description

Computes the derivative plot of the Pareto QQ-plot. These values can be plotted as a function of the data or as a function of the tail parameter k.

Usage

ParetoQQ_der(data, k = FALSE, plot = TRUE, 
             main = "Derivative plot of Pareto QQ-plot", ...)

Arguments

Details

The derivative plot of a Pareto QQ-plot is

(k,H_{k,n})

or

(\log X_{n-k,n}, H_{k,n})

with H_{k,n} the Hill estimates.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

ParetoQQ, Hill, MeanExcess, LognormalQQ_der, WeibullQQ_der

Examples

data(norwegianfire)

# Pareto QQ-plot for Norwegian Fire Insurance data for claims in 1976.
ParetoQQ(norwegianfire$size[norwegianfire$year==76])

# Derivate plot
ParetoQQ_der(norwegianfire$size[norwegianfire$year==76])

Weissman estimator of small exceedance probabilities and large return periods

Description

Compute estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the approach of Weissman (1978).

Usage

Prob(data, gamma, q, plot = FALSE, add = FALSE, 
     main = "Estimates of small exceedance probability", ...)
          
Return(data, gamma, q, plot = FALSE, add = FALSE, 
       main = "Estimates of large return period", ...)
        
Weissman.p(data, gamma, q, plot = FALSE, add = FALSE, 
           main = "Estimates of small exceedance probability", ...)
          
Weissman.r(data, gamma, q, plot = FALSE, add = FALSE, 
           main = "Estimates of large return period", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Weissman.p and Weissman.r are the same functions as Prob and Return but with a different name for compatibility with the old S-Plus code.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Weissman, I. (1978). "Estimation of Parameters and Large Quantiles Based on the k Largest Observations." Journal of the American Statistical Association, 73, 812–815.

See Also

Quant

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Hill estimator
H <- Hill(SOAdata)

# Exceedance probability
q <- 10^6
# Weissman estimator
Prob(SOAdata,gamma=H$gamma,q=q,plot=TRUE)

# Return period
q <- 10^6
# Weissman estimator
Return(SOAdata,gamma=H$gamma,q=q,plot=TRUE)

Estimator of small exceedance probabilities and large return periods using EPD

Description

Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the parameters from the EPD fit.

Usage

ProbEPD(data, q, gamma, kappa, tau, plot = FALSE, add = FALSE, 
        main = "Estimates of small exceedance probability", ...)

ReturnEPD(data, q, gamma, kappa, tau, plot = FALSE, add = FALSE, 
          main = "Estimates of large return period", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Joossens, E. and Segers, J. (2009). "Second-Order Refined Peaks-Over-Threshold Modelling for Heavy-Tailed Distributions." Journal of Statistical Planning and Inference, 139, 2800–2815.

See Also

EPD, Prob

Examples

data(secura)

# EPD estimates for the EVI
epd <- EPD(secura$size, plot=TRUE)

# Compute exceedance probabilities
q <- 10^7
ProbEPD(secura$size, q=q, gamma=epd$gamma, kappa=epd$kappa, tau=epd$tau, plot=TRUE)

# Compute return periods
ReturnEPD(secura$size, q=q, gamma=epd$gamma, kappa=epd$kappa, tau=epd$tau, 
          plot=TRUE, ylim=c(0,10^4))

Estimator of small exceedance probabilities and large return periods using generalised Hill

Description

Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the generalised Hill estimates for the EVI.

Usage

ProbGH(data, gamma, q, plot = FALSE, add = FALSE, 
       main = "Estimates of small exceedance probability", ...)

ReturnGH(data, gamma, q, plot = FALSE, add = FALSE, 
         main = "Estimates of large return period", ...)

Arguments

Details

See Section 4.2.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Beirlant, J., Vynckier, P. and Teugels, J.L. (1996). "Excess Function and Estimation of the Extreme-value Index". Bernoulli, 2, 293–318.

See Also

QuantGH, genHill, ProbMOM, Prob

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Hill estimator
H <- Hill(SOAdata)
# Generalised Hill estimator
gH <- genHill(SOAdata, H$gamma)

# Exceedance probability
q <- 10^7
ProbGH(SOAdata, gamma=gH$gamma, q=q, plot=TRUE)

# Return period
q <- 10^7
ReturnGH(SOAdata, gamma=gH$gamma, q=q, plot=TRUE)

Estimator of small exceedance probabilities and large return periods using GPD-MLE

Description

Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the GPD fit for the peaks over a threshold.

Usage

ProbGPD(data, gamma, sigma, q, plot = FALSE, add = FALSE, 
        main = "Estimates of small exceedance probability", ...)

ReturnGPD(data, gamma, sigma, q, plot = FALSE, add = FALSE, 
          main = "Estimates of large return period", ...)

Arguments

Details

See Section 4.2.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

QuantGPD, GPDmle, Prob

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# GPD-ML estimator
pot <- GPDmle(SOAdata)

# Exceedance probability
q <- 10^7
ProbGPD(SOAdata, gamma=pot$gamma, sigma=pot$sigma, q=q, plot=TRUE)

# Return period
q <- 10^7
ReturnGPD(SOAdata, gamma=pot$gamma, sigma=pot$sigma, q=q, plot=TRUE)

Estimator of small exceedance probabilities and large return periods using MOM

Description

Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the Method of Moments estimates for the EVI.

Usage

ProbMOM(data, gamma, q, plot = FALSE, add = FALSE, 
        main = "Estimates of small exceedance probability", ...)

ReturnMOM(data, gamma, q, plot = FALSE, add = FALSE, 
          main = "Estimates of large return period", ...)

Arguments

Details

See Section 4.2.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Dekkers, A.L.M, Einmahl, J.H.J. and de Haan, L. (1989). "A Moment Estimator for the Index of an Extreme-value Distribution." Annals of Statistics, 17, 1833–1855.

See Also

QuantMOM, Moment, ProbGH, Prob

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# MOM estimator
M <- Moment(SOAdata)

# Exceedance probability
q <- 10^7
ProbMOM(SOAdata, gamma=M$gamma, q=q, plot=TRUE)

# Return period
q <- 10^7
ReturnMOM(SOAdata, gamma=M$gamma, q=q, plot=TRUE)

Estimator of small tail probability in regression

Description

Estimator of small tail probability 1-F_i(q) in the regression case where \gamma is constant and the regression modelling is thus only solely placed on the scale parameter.

Usage

ProbReg(Z, A, q, plot = FALSE, add = FALSE, 
        main = "Estimates of small exceedance probability", ...)

Arguments

Details

The estimator is defined as

1-\hat{F}_i(q) = \hat{A}(i/n) (k+1)/(n+1) (q/Z_{n-k,n})^{-1/H_{k,n}},

with H_{k,n} the Hill estimator. Here, it is assumed that we have equidistant covariates x_i=i/n.

See Section 4.4.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

QuantReg, ScaleReg, Prob

Examples

data(norwegianfire)

Z <- norwegianfire$size[norwegianfire$year==76]

i <- 100
n <- length(Z)

# Scale estimator in i/n
A <- ScaleReg(i/n, Z, h=0.5, kernel = "epanechnikov")$A

# Small exceedance probability
q <- 10^6
ProbReg(Z, A, q, plot=TRUE)

# Large quantile
p <- 10^(-5)
QuantReg(Z, A, p, plot=TRUE)

Weissman estimator of extreme quantiles

Description

Compute estimates of an extreme quantile Q(1-p) using the approach of Weissman (1978).

Usage

Quant(data, gamma, p, plot = FALSE, add = FALSE, 
      main = "Estimates of extreme quantile", ...)
            
Weissman.q(data, gamma, p, plot = FALSE, add = FALSE, 
           main = "Estimates of extreme quantile", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Weissman.q is the same function but with a different name for compatibility with the old S-Plus code.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Weissman, I. (1978). "Estimation of Parameters and Large Quantiles Based on the k Largest Observations." Journal of the American Statistical Association, 73, 812–815.

See Also

Prob, Quant.2oQV

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Hill estimator
H <- Hill(SOAdata)
# Bias-reduced estimator (QV)
H_QV <- Hill.2oQV(SOAdata)

# Exceedance probability
p <- 10^(-5)
# Weissman estimator
Quant(SOAdata, gamma=H$gamma, p=p, plot=TRUE)

# Second order Weissman estimator (QV)
Quant.2oQV(SOAdata, gamma=H_QV$gamma, beta=H_QV$beta, b=H_QV$b, p=p,
           add=TRUE, lty=2)

Second order refined Weissman estimator of extreme quantiles (QV)

Description

Compute second order refined Weissman estimator of extreme quantiles Q(1-p) using the quantile view.

Usage

Quant.2oQV(data, gamma, b, beta, p, plot = FALSE, add = FALSE, 
           main = "Estimates of extreme quantile", ...)
                
Weissman.q.2oQV(data, gamma, b, beta, p, plot = FALSE, add = FALSE, 
                main = "Estimates of extreme quantile", ...)

Arguments

Details

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Weissman.q.2oQV is the same function but with a different name for compatibility with the old S-Plus code.

Value

A list with following components:

Author(s)

Tom Reynkens based on S-Plus code from Yuri Goegebeur.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

Quant, Hill.2oQV

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Hill estimator
H <- Hill(SOAdata)
# Bias-reduced estimator (QV)
H_QV <- Hill.2oQV(SOAdata)

# Exceedance probability
p <- 10^(-5)
# Weissman estimator
Quant(SOAdata, gamma=H$gamma, p=p, plot=TRUE)

# Second order Weissman estimator (QV)
Quant.2oQV(SOAdata, gamma=H_QV$gamma, beta=H_QV$beta, b=H_QV$b, p=p,
           add=TRUE, lty=2)

Estimator of extreme quantiles using generalised Hill

Description

Compute estimates of an extreme quantile Q(1-p) using generalised Hill estimates of the EVI.

Usage

QuantGH(data, gamma, p, plot = FALSE, add = FALSE, 
        main = "Estimates of extreme quantile", ...)

Arguments

Details

See Section 4.2.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Beirlant, J., Vynckier, P. and Teugels, J.L. (1996). "Excess Function and Estimation of the Extreme-value Index". Bernoulli, 2, 293–318.

See Also

ProbGH, genHill, QuantMOM, Quant

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# Hill estimator
H <- Hill(SOAdata)
# Generalised Hill estimator
gH <- genHill(SOAdata, H$gamma)

# Large quantile
p <- 10^(-5)
QuantGH(SOAdata, p=p, gamma=gH$gamma, plot=TRUE)

Estimator of extreme quantiles using GPD-MLE

Description

Computes estimates of an extreme quantile Q(1-p) using the GPD fit for the peaks over a threshold.

Usage

QuantGPD(data, gamma, sigma, p, plot = FALSE, add = FALSE, 
         main = "Estimates of extreme quantile", ...)

Arguments

Details

See Section 4.2.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

See Also

ProbGPD, GPDmle, Quant

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# GPD-ML estimator
pot <- GPDmle(SOAdata)

# Large quantile
p <- 10^(-5)
QuantGPD(SOAdata, p=p, gamma=pot$gamma, sigma=pot$sigma, plot=TRUE)

Estimator of extreme quantiles using MOM

Description

Compute estimates of an extreme quantile Q(1-p) using the Method of Moments estimates of the EVI.

Usage

QuantMOM(data, gamma, p, plot = FALSE, add = FALSE, 
         main = "Estimates of extreme quantile", ...)

Arguments

Details

See Section 4.2.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant J., Goegebeur Y., Segers, J. and Teugels, J. (2004). Statistics of Extremes: Theory and Applications, Wiley Series in Probability, Wiley, Chichester.

Dekkers, A.L.M, Einmahl, J.H.J. and de Haan, L. (1989). "A Moment Estimator for the Index of an Extreme-value Distribution." Annals of Statistics, 17, 1833–1855.

See Also

ProbMOM, Moment, QuantGH, Quant

Examples

data(soa)

# Look at last 500 observations of SOA data
SOAdata <- sort(soa$size)[length(soa$size)-(0:499)]

# MOM estimator
M <- Moment(SOAdata)

# Large quantile
p <- 10^(-5)
QuantMOM(SOAdata, p=p, gamma=M$gamma, plot=TRUE)

Estimator of extreme quantiles in regression

Description

Estimator of extreme quantile Q_i(1-p) in the regression case where \gamma is constant and the regression modelling is thus only solely placed on the scale parameter.

Usage

QuantReg(Z, A, p, plot = FALSE, add = FALSE, 
         main = "Estimates of extreme quantile", ...)

Arguments

Details

The estimator is defined as

\hat{Q}_i(1-p) = Z_{n-k,n} ((k+1)/((n+1)\times p) \hat{A}(i/n))^{H_{k,n}},

with H_{k,n} the Hill estimator. Here, it is assumed that we have equidistant covariates x_i=i/n.

See Section 4.4.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

ProbReg, ScaleReg, Quant

Examples

data(norwegianfire)

Z <- norwegianfire$size[norwegianfire$year==76]

i <- 100
n <- length(Z)

# Scale estimator in i/n
A <- ScaleReg(i/n, Z, h=0.5, kernel = "epanechnikov")$A

# Small exceedance probability
q <- 10^6
ProbReg(Z, A, q, plot=TRUE)

# Large quantile
p <- 10^(-5)
QuantReg(Z, A, p, plot=TRUE)

Scale estimator

Description

Computes the estimator for the scale parameter as described in Beirlant et al. (2016).

Usage

Scale(data, gamma = NULL, logk = FALSE, plot = FALSE, add = FALSE, 
      main = "Estimates of scale parameter", ...)

Arguments

Details

The scale estimates are computed based on the following model for the CDF: 1-F(x) = A x^{-1/\gamma}, where A:= C^{1/\gamma} is the scale parameter:

\hat{A}_{k,n}=(k+1)/(n+1) X_{n-k,n}^{1/H_{k,n}}

where H_{k,n} are the Hill estimates.

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Schoutens, W., De Spiegeleer, J., Reynkens, T. and Herrmann, K. (2016). "Hunting for Black Swans in the European Banking Sector Using Extreme Value Analysis." In: Jan Kallsen and Antonis Papapantoleon (eds.), Advanced Modelling in Mathematical Finance, Springer International Publishing, Switzerland, pp. 147–166.

See Also

ScaleEPD, Scale.2o, Hill

Examples

data(secura)

# Hill estimator
H <- Hill(secura$size)

# Scale estimator
S <- Scale(secura$size, gamma=H$gamma, plot=FALSE)

# Plot logarithm of scale          
plot(S$k,log(S$A), xlab="k", ylab="log(Scale)", type="l")

Bias-reduced scale estimator using second order Hill estimator

Description

Computes the bias-reduced estimator for the scale parameter using the second-order Hill estimator.

Usage

Scale.2o(data, gamma, b, beta, logk = FALSE, plot = FALSE, add = FALSE, 
         main = "Estimates of scale parameter", ...)

Arguments

Details

The scale estimates are computed based on the following model for the CDF: 1-F(x) = A x^{-1/\gamma} ( 1+ bx^{-\beta}(1+o(1)) ), where A:= C^{1/\gamma} is the scale parameter.

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Schoutens, W., De Spiegeleer, J., Reynkens, T. and Herrmann, K. (2016). "Hunting for Black Swans in the European Banking Sector Using Extreme Value Analysis." In: Jan Kallsen and Antonis Papapantoleon (eds.), Advanced Modelling in Mathematical Finance, Springer International Publishing, Switzerland, pp. 147–166.

See Also

Scale, ScaleEPD, Hill.2oQV

Examples

data(secura)

# Hill estimator
H <- Hill(secura$size)
# Bias-reduced Hill estimator
H2o <- Hill.2oQV(secura$size)

# Scale estimator
S <- Scale(secura$size, gamma=H$gamma, plot=FALSE)
# Bias-reduced scale estimator
S2o <- Scale.2o(secura$size, gamma=H2o$gamma, b=H2o$b, 
          beta=H2o$beta, plot=FALSE)

# Plot logarithm of scale             
plot(S$k,log(S$A), xlab="k", ylab="log(Scale)", type="l")
lines(S2o$k,log(S2o$A), lty=2)

Bias-reduced scale estimator using EPD estimator

Description

Computes the bias-reduced estimator for the scale parameter using the EPD estimator (Beirlant et al., 2016).

Usage

ScaleEPD(data, gamma, kappa, logk = FALSE, plot = FALSE, add = FALSE, 
         main = "Estimates of scale parameter", ...)

Arguments

Details

The scale estimates are computed based on the following model for the CDF: 1-F(x) = A x^{-1/\gamma} ( 1+ bx^{-\beta}(1+o(1)) ), where A:= C^{1/\gamma} is the scale parameter. Using the EPD approach we replace bx^{-\beta} by \kappa/\gamma.

See Section 4.2.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Schoutens, W., De Spiegeleer, J., Reynkens, T. and Herrmann, K. (2016). "Hunting for Black Swans in the European Banking Sector Using Extreme Value Analysis." In: Jan Kallsen and Antonis Papapantoleon (eds.), Advanced Modelling in Mathematical Finance, Springer International Publishing, Switzerland, pp. 147–166.

See Also

Scale, Scale.2o, EPD

Examples

data(secura)

# Hill estimator
H <- Hill(secura$size)
# EPD estimator
epd <- EPD(secura$size)

# Scale estimator
S <- Scale(secura$size, gamma=H$gamma, plot=FALSE)
# Bias-reduced scale estimator
Sepd <- ScaleEPD(secura$size, gamma=epd$gamma, kappa=epd$kappa, plot=FALSE)

# Plot logarithm of scale             
plot(S$k,log(S$A), xlab="k", ylab="log(Scale)", type="l")
lines(Sepd$k,log(Sepd$A), lty=2)

Scale estimator in regression

Description

Estimator of the scale parameter in the regression case where \gamma is constant and the regression modelling is thus placed solely on the scale parameter.

Usage

ScaleReg(s, Z, kernel = c("normal", "uniform", "triangular", "epanechnikov", "biweight"), 
         h, plot = TRUE, add = FALSE, main = "Estimates of scale parameter", ...)

Arguments

Details

The scale estimator is computed as

\hat{A}(s) = 1/(k+1) \sum_{i=1}^n 1_{Z_i>Z_{n-k,n}} K_h(s-i/n)

with K_h(x)=K(x/h)/h, K the kernel function and h the bandwidth. Here, it is assumed that we have equidistant covariates x_i=i/n.

See Section 4.4.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

ProbReg, QuantReg, scale, Hill

Examples

data(norwegianfire)

Z <- norwegianfire$size[norwegianfire$year==76]

i <- 100
n <- length(Z)

# Scale estimator in i/n
A <- ScaleReg(i/n, Z, h=0.5, kernel = "epanechnikov")$A

# Small exceedance probability
q <- 10^6
ProbReg(Z, A, q, plot=TRUE)

# Large quantile
p <- 10^(-5)
QuantReg(Z, A, p, plot=TRUE)

Spliced distribution

Description

Density, distribution function, quantile function and random generation for the fitted spliced distribution.

Usage

dSplice(x, splicefit, log = FALSE)

pSplice(x, splicefit, lower.tail = TRUE, log.p = FALSE)

qSplice(p, splicefit, lower.tail = TRUE, log.p = FALSE)

rSplice(n, splicefit)

Arguments

Details

See Reynkens et al. (2017) and Section 4.3 in Albrecher et al. (2017) for details.

Value

dSplice gives the density function evaluated in x, pSplice the CDF evaluated in x and qSplice the quantile function evaluated in p. The length of the result is equal to the length of x or p.

rSplice returns a random sample of length n.

Author(s)

Tom Reynkens with R code from Roel Verbelen for the mixed Erlang PDF, CDF and quantiles.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758.

See Also

VaR, SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD, SpliceECDF, SpliceLL, SplicePP

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



p <- seq(0, 1, 0.01)

# Plot of splicing quantiles
plot(p, qSplice(p, splicefit), type="l", xlab="p", ylab="Q(p)")

# Plot of VaR
plot(p, VaR(p, splicefit), type="l", xlab="p", ylab=bquote(VaR[p]))



# Random sample from spliced distribution
x <- rSplice(1000, splicefit)


## End(Not run)

Plot of fitted and empirical survival function

Description

This function plots the fitted survival function of the spliced distribution together with the empirical survival function (determined using the Empirical CDF (ECDF)). Moreover, 100(1-\alpha)\% confidence bands are added.

Usage

SpliceECDF(x, X, splicefit, alpha = 0.05, ...)

Arguments

Details

Use SpliceTB for censored data.

Confidence bands are determined using the Dvoretzky-Kiefer-Wolfowitz inequality (Massart, 1990).

See Reynkens et al. (2017) and Section 4.3.1 in Albrecher et al. (2017) for more details.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Massart, P. (1990). The Tight Constant in the Dvoretzky-Kiefer-Wolfowitz Inequality. Annals of Probability, 18, 1269–1283.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758.

See Also

SpliceTB, pSplice, ecdf, SpliceFitPareto, SpliceFitGPD, SpliceLL, SplicePP, SpliceQQ

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

Splicing fit

Description

Create an S3 object using ME-Pa or ME-GPD splicing fit obtained from SpliceFitPareto, SpliceFiticPareto or SpliceFitGPD.

Usage

SpliceFit(const, trunclower, t, type, MEfit, EVTfit, loglik = NULL, IC = NULL)

Arguments

Details

See Reynkens et al. (2017) and Section 4.3 in Albrecher et al. (2017) for details.

Value

An S3 object containing the above input arguments and values for \pi, the splicing weights. These splicing weights are equal to

\pi_1=const_1, \pi_2=const_2-const_1, ...,\pi_{l+1}=1-const_l=1-(\pi_1+...+\pi_l)

when l\ge 2 and

\pi_1=const_1, \pi_2=1-const_1=1-\pi_1

when l=1, where l is the length of const.

A summary method is available.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

MEfit, EVTfit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

# Create MEfit object
mefit <- MEfit(p=c(0.65,0.35), shape=c(39,58), theta=16.19, M=2)

# Create EVTfit object
evtfit <- EVTfit(gamma=c(0.76,0.64), endpoint=c(39096, Inf))

# Create SpliceFit object
splicefit <- SpliceFit(const=c(0.5,0.996), trunclower=0, t=c(1020,39096), type=c("ME","TPa","Pa"),
                       MEfit=mefit, EVTfit=evtfit)

# Show summary
summary(splicefit)

Splicing of mixed Erlang and GPD using POT-MLE

Description

Fit spliced distribution of a mixed Erlang distribution and a Generalised Pareto Distribution (GPD). The parameters of the GPD are determined using the POT-MLE approach.

Usage

SpliceFitGPD(X, const = NULL, tsplice = NULL, M = 3, s = 1:10, trunclower = 0, 
             ncores = NULL, criterium = c("BIC","AIC"), reduceM = TRUE,
             eps = 10^(-3), beta_tol = 10^(-5), maxiter = Inf)

Arguments

Details

See Reynkens et al. (2017), Section 4.3.1 of Albrecher et al. (2017) and Verbelen et al. (2015) for details. The code follows the notation of the latter. Initial values follow from Verbelen et al. (2016).

Value

A SpliceFit object.

Author(s)

Tom Reynkens with R code from Roel Verbelen for fitting the mixed Erlang distribution.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758.

Verbelen, R., Antonio, K. and Claeskens, G. (2016). "Multivariate Mixtures of Erlangs for Density Estimation Under Censoring." Lifetime Data Analysis, 22, 429–455.

See Also

SpliceFitPareto, SpliceFiticPareto, Splice, GPDfit

Examples

## Not run: 

# GPD random sample
X <- rgpd(1000, gamma = 0.5, sigma = 2)

# Splice ME and GPD
splicefit <- SpliceFitGPD(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

Splicing of mixed Erlang and Pareto

Description

Fit spliced distribution of a mixed Erlang distribution and Pareto distribution(s). The shape parameter(s) of the Pareto distribution(s) is determined using the Hill estimator.

Usage

SpliceFitPareto(X, const = NULL, tsplice = NULL, M = 3, s = 1:10, trunclower = 0, 
                truncupper = Inf, EVTtruncation = FALSE, ncores = NULL, 
                criterium = c("BIC","AIC"), reduceM = TRUE,
                eps = 10^(-3), beta_tol = 10^(-5), maxiter = Inf)
              
SpliceFitHill(X, const = NULL, tsplice = NULL, M = 3, s = 1:10, trunclower = 0, 
              truncupper = Inf, EVTtruncation = FALSE, ncores = NULL,
              criterium = c("BIC","AIC"), reduceM = TRUE,
              eps = 10^(-3), beta_tol = 10^(-5), maxiter = Inf)

Arguments

Details

See Reynkens et al. (2017), Section 4.3.1 of Albrecher et al. (2017) and Verbelen et al. (2015) for details. The code follows the notation of the latter. Initial values follow from Verbelen et al. (2016).

The SpliceFitHill function is the same function but with a different name for compatibility with old versions of the package.

Use SpliceFiticPareto when censoring is present.

Value

A SpliceFit object.

Author(s)

Tom Reynkens with R code from Roel Verbelen for fitting the mixed Erlang distribution.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Beirlant, J., Fraga Alves, M.I. and Gomes, M.I. (2016). "Tail fitting for Truncated and Non-truncated Pareto-type Distributions." Extremes, 19, 429–462.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758.

Verbelen, R., Antonio, K. and Claeskens, G. (2016). "Multivariate Mixtures of Erlangs for Density Estimation Under Censoring." Lifetime Data Analysis, 22, 429–455.

See Also

SpliceFiticPareto, SpliceFitGPD, Splice, Hill, trHill

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

Splicing of mixed Erlang and Pareto for interval censored data

Description

Fit spliced distribution of a mixed Erlang distribution and a Pareto distribution adapted for interval censoring and truncation.

Usage

SpliceFiticPareto(L, U, censored, tsplice, M = 3, s = 1:10, trunclower = 0, 
                  truncupper = Inf, ncores = NULL, criterium = c("BIC", "AIC"), 
                  reduceM = TRUE, eps = 10^(-3), beta_tol = 10^(-5), maxiter = Inf, 
                  cpp = FALSE)

Arguments

Details

See Reynkens et al. (2017), Section 4.3.2 of Albrecher et al. (2017) and Verbelen et al. (2015) for details. The code follows the notation of the latter. Initial values follow from Verbelen et al. (2016).

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u.

Value

A SpliceFit object.

Author(s)

Tom Reynkens based on R code from Roel Verbelen for fitting the mixed Erlang distribution (without splicing).

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758.

Verbelen, R., Antonio, K. and Claeskens, G. (2016). "Multivariate Mixtures of Erlangs for Density Estimation Under Censoring." Lifetime Data Analysis, 22, 429–455.

See Also

SpliceFitPareto, SpliceFitGPD, Splice

Examples

## Not run: 

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X,Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Splice ME and Pareto
splicefit <- SpliceFiticPareto(L=Z, U=U, censored=censored, tsplice=quantile(Z,0.9))



x <- seq(0,20,0.1)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")


# Fitted survival function and Turnbull survival function 
SpliceTB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# Log-log plot with Turnbull survival function and fitted survival function
SpliceLL_TB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# PP-plot of Turnbull survival function and fitted survival function
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

# PP-plot of Turnbull survival function and 
# fitted survival function with log-scales
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit, log=TRUE)

# QQ-plot using Turnbull survival function and fitted survival function
SpliceQQ_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

## End(Not run)

LL-plot with fitted and empirical survival function

Description

This function plots the logarithm of the empirical survival function (determined using the Empirical CDF (ECDF)) versus the logarithm of the data. Moreover, the logarithm of the fitted survival function of the spliced distribution is added.

Usage

SpliceLL(x = sort(X), X, splicefit, plot = TRUE, main = "Splicing LL-plot", ...)

Arguments

Details

The LL-plot consists of the points

(\log(x_{i,n}), \log(1-\hat{F}(x_{i,n})))

for i=1,\ldots,n with n the length of the data, x_{i,n} the i-th smallest observation and \hat{F} the empirical distribution function. Then, the line

(\log(x), \log(1-\hat{F}_{spliced}(x))),

with \hat{F}_{spliced} the fitted spliced distribution function, is added.

Use SpliceLL_TB for censored data.

See Reynkens et al. (2017) and Section 4.3.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceLL_TB, pSplice, ecdf, SpliceFitPareto, SpliceFitGPD, SpliceECDF, SplicePP, SpliceQQ

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

LL-plot with fitted and Turnbull survival function

Description

This function plots the logarithm of the Turnbull survival function (which is suitable for interval censored data) versus the logarithm of the data. Moreover, the logarithm of the fitted survival function of the spliced distribution is added.

Usage

SpliceLL_TB(x = sort(L), L, U = L, censored, splicefit, plot = TRUE,
            main = "Splicing LL-plot", ...)

Arguments

Details

The LL-plot consists of the points

(\log(L_{i,n}), \log(1-\hat{F}^{TB}(L_{i,n})))

for i=1,\ldots,n with n the length of the data, x_{i,n} the i-th smallest observation and \hat{F}^{TB} the Turnbull estimator for the distribution function. Then, the line

(\log(x), \log(1-\hat{F}_{spliced}(x))),

with \hat{F}_{spliced} the fitted spliced distribution function, is added.

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u. The limits trunclower and truncupper are obtained from the SpliceFit object.

If the interval package is installed, the icfit function is used to compute the Turnbull estimator. Otherwise, survfit.formula from survival is used.

Use SpliceLL for non-censored data.

See Reynkens et al. (2017) and Section 4.3.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceLL, pSplice, Turnbull, icfit, SpliceFiticPareto, SpliceTB, SplicePP_TB, SpliceQQ_TB

Examples

## Not run: 

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X,Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Splice ME and Pareto
splicefit <- SpliceFiticPareto(L=Z, U=U, censored=censored, tsplice=quantile(Z,0.9))



x <- seq(0,20,0.1)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")


# Fitted survival function and Turnbull survival function 
SpliceTB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# Log-log plot with Turnbull survival function and fitted survival function
SpliceLL_TB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# PP-plot of Turnbull survival function and fitted survival function
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

# PP-plot of Turnbull survival function and 
# fitted survival function with log-scales
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit, log=TRUE)

# QQ-plot using Turnbull survival function and fitted survival function
SpliceQQ_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

## End(Not run)

PP-plot with fitted and empirical survival function

Description

This function plots the fitted survival function of the spliced distribution versus the empirical survival function (determined using the Empirical CDF (ECDF)).

Usage

SplicePP(X, splicefit, x = sort(X), log = FALSE, plot = TRUE, 
         main = "Splicing PP-plot", ...)

Arguments

Details

The PP-plot consists of the points

(1-\hat{F}(x_{i,n}), 1-\hat{F}_{spliced}(x_{i,n})))

for i=1,\ldots,n with n the length of the data, x_{i,n} the i-th smallest observation, \hat{F} the empirical distribution function and \hat{F}_{spliced} the fitted spliced distribution function. The minus-log version of the PP-plot consists of

(-\log(1-\hat{F}(x_{i,n})), -\log(1-\hat{F}_{spliced}(x_{i,n})))).

Use SplicePP_TB for censored data.

See Reynkens et al. (2017) and Section 4.3.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SplicePP_TB, pSplice, ecdf, SpliceFitPareto, SpliceFitGPD, SpliceECDF, SpliceLL, SpliceQQ

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

PP-plot with fitted and Turnbull survival function

Description

This function plots the fitted survival function of the spliced distribution versus the Turnbull survival function (which is suitable for interval censored data).

Usage

SplicePP_TB(L, U = L, censored, splicefit, x = NULL, log = FALSE, plot = TRUE,
            main = "Splicing PP-plot", ...)

Arguments

Details

The PP-plot consists of the points

(1-\hat{F}^{TB}(x_i), 1-\hat{F}_{spliced}(x_i)))

for i=1,\ldots,n with n the length of the data, x_i=\hat{Q}^{TB}(p_i) where p_i=i/(n+1), \hat{Q}^{TB} is the quantile function obtained using the Turnbull estimator, \hat{F}^{TB} the Turnbull estimator for the distribution function and \hat{F}_{spliced} the fitted spliced distribution function. The minus-log version of the PP-plot consists of

(-\log(1-\hat{F}^{TB}(x_i)), -\log(1-\hat{F}_{spliced}(x_i))).

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u. The limits trunclower and truncupper are obtained from the SpliceFit object.

If the interval package is installed, the icfit function is used to compute the Turnbull estimator. Otherwise, survfit.formula from survival is used.

Use SplicePP for non-censored data.

See Reynkens et al. (2017) and Section 4.3.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SplicePP, pSplice, Turnbull, icfit, SpliceFiticPareto, SpliceTB, SpliceLL_TB, SpliceQQ_TB

Examples

## Not run: 

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X,Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Splice ME and Pareto
splicefit <- SpliceFiticPareto(L=Z, U=U, censored=censored, tsplice=quantile(Z,0.9))



x <- seq(0,20,0.1)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")


# Fitted survival function and Turnbull survival function 
SpliceTB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# Log-log plot with Turnbull survival function and fitted survival function
SpliceLL_TB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# PP-plot of Turnbull survival function and fitted survival function
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

# PP-plot of Turnbull survival function and 
# fitted survival function with log-scales
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit, log=TRUE)

# QQ-plot using Turnbull survival function and fitted survival function
SpliceQQ_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

## End(Not run)

Splicing quantile plot

Description

Computes the empirical quantiles of a data vector and the theoretical quantiles of the fitted spliced distribution. These quantiles are then plotted in a splicing QQ-plot with the theoretical quantiles on the x-axis and the empirical quantiles on the y-axis.

Usage

SpliceQQ(X, splicefit, p = NULL, plot = TRUE, main = "Splicing QQ-plot", ...)

Arguments

Details

This QQ-plot is given by

(Q(p_j), \hat{Q}(p_j)),

for j=1,\ldots,n where Q is the quantile function of the fitted splicing model and \hat{Q} is the empirical quantile function and p_j=j/(n+1).

See Reynkens et al. (2017) and Section 4.3.1 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceQQ_TB, qSplice, SpliceFitPareto, SpliceFitGPD, SpliceECDF, SpliceLL, SplicePP

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



x <- seq(0, 20, 0.01)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")



# Fitted survival function and empirical survival function 
SpliceECDF(x, X, splicefit)

# Log-log plot with empirical survival function and fitted survival function
SpliceLL(x, X, splicefit)

# PP-plot of empirical survival function and fitted survival function
SplicePP(X, splicefit)

# PP-plot of empirical survival function and 
# fitted survival function with log-scales
SplicePP(X, splicefit, log=TRUE)

# Splicing QQ-plot
SpliceQQ(X, splicefit)

## End(Not run)

Splicing quantile plot using Turnbull estimator

Description

This function plots the fitted quantile function of the spliced distribution versus quantiles based on the Turnbull survival function (which is suitable for interval censored data).

Usage

SpliceQQ_TB(L, U = L, censored, splicefit, p = NULL,
            plot = TRUE, main = "Splicing QQ-plot", ...)

Arguments

Details

This QQ-plot is given by

(Q(p_j), \hat{Q}^{TB}(p_j)),

for j=1,\ldots,n where Q is the quantile function of the fitted splicing model, \hat{Q}^{TB} the quantile function obtained using the Turnbull estimator and p_j=j/(n+1).

If the interval package is installed, the icfit function is used to compute the Turnbull estimator. Otherwise, survfit.formula from survival is used.

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u. The limits trunclower and truncupper are obtained from the SpliceFit object.

Use SpliceQQ for non-censored data.

See Reynkens et al. (2017) and Section 4.3.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceQQ, qSplice, Turnbull, icfit, SpliceFiticPareto, SpliceTB, SplicePP_TB, SpliceLL_TB

Examples

## Not run: 

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X,Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Splice ME and Pareto
splicefit <- SpliceFiticPareto(L=Z, U=U, censored=censored, tsplice=quantile(Z,0.9))



x <- seq(0,20,0.1)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")


# Fitted survival function and Turnbull survival function 
SpliceTB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# Log-log plot with Turnbull survival function and fitted survival function
SpliceLL_TB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# PP-plot of Turnbull survival function and fitted survival function
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

# PP-plot of Turnbull survival function and 
# fitted survival function with log-scales
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit, log=TRUE)

# QQ-plot using Turnbull survival function and fitted survival function
SpliceQQ_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

## End(Not run)

Plot of fitted and Turnbull survival function

Description

This function plots the fitted survival function of the spliced distribution together with the Turnbull survival function (which is suitable for interval censored data). Moreover, 100(1-\alpha)\% confidence intervals are added.

Usage

SpliceTB(x = sort(L), L, U = L, censored, splicefit, alpha = 0.05, ...)

Arguments

Details

Right censored data should be entered as L=l and U=truncupper, and left censored data should be entered as L=trunclower and U=u. The limits trunclower and truncupper are obtained from the SpliceFit object.

Use SpliceECDF for non-censored data.

See Reynkens et al. (2017) and Section 4.3.2 in Albrecher et al. (2017) for more details.

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

SpliceECDF, pSplice, Turnbull, SpliceFiticPareto, SpliceLL_TB, SplicePP_TB, SpliceQQ_TB

Examples

## Not run: 

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X,Y)

# Censoring indicator
censored <- (X>Y)

# Right boundary
U <- Z
U[censored] <- Inf

# Splice ME and Pareto
splicefit <- SpliceFiticPareto(L=Z, U=U, censored=censored, tsplice=quantile(Z,0.9))



x <- seq(0,20,0.1)

# Plot of spliced CDF
plot(x, pSplice(x, splicefit), type="l", xlab="x", ylab="F(x)")

# Plot of spliced PDF
plot(x, dSplice(x, splicefit), type="l", xlab="x", ylab="f(x)")


# Fitted survival function and Turnbull survival function 
SpliceTB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# Log-log plot with Turnbull survival function and fitted survival function
SpliceLL_TB(x, L=Z, U=U, censored=censored, splicefit=splicefit)


# PP-plot of Turnbull survival function and fitted survival function
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

# PP-plot of Turnbull survival function and 
# fitted survival function with log-scales
SplicePP_TB(L=Z, U=U, censored=censored, splicefit=splicefit, log=TRUE)

# QQ-plot using Turnbull survival function and fitted survival function
SpliceQQ_TB(L=Z, U=U, censored=censored, splicefit=splicefit)

## End(Not run)

Turnbull estimator

Description

Computes the Turnbull estimator for the survival function of interval censored data.

Usage

Turnbull(x, L, R, censored, trunclower = 0, truncupper = Inf,
         conf.type = "plain", conf.int = 0.95)

Arguments

Details

We consider the random interval censoring model where one observes L \le R and where the variable of interest X lies between L and R.

Right censored data should be entered as L=l and R=truncupper, and right censored data should be entered as L=trunclower and R=r.

This function calls survfit.formula from survival.

See Section 4.3.2 in Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Turnbull, B. W. (1974). "Nonparametric Estimation of a Survivorship Function with Doubly Censored Data." Journal of the American Statistical Association, 69, 169–173.

Turnbull, B. W. (1976). "The Empirical Distribution Function with Arbitrarily Grouped, Censored and Truncated Data." Journal of the Royal Statistical Society: Series B (Methodological), 38, 290–295.

See Also

survfit.formula, KaplanMeier

Examples

L <- 1:10
R <- c(1, 2.5, 3, 4, 5.5, 6, 7.5, 8.25, 9, 10.5)
censored <- c(0, 1, 0, 0, 1, 0, 1, 1, 0, 1)

x <- seq(0, 12, 0.1)

# Turnbull estimator
plot(x, Turnbull(x, L, R, censored)$cdf, type="s", ylab="Turnbull estimator")

VaR of splicing fit

Description

Compute Value-at-Risk (VaR_{1-p}=Q(1-p)) of the fitted spliced distribution.

Usage

VaR(p, splicefit)

Arguments

Details

See Reynkens et al. (2017) and Section 4.6 of Albrecher et al. (2017) for details.

Note that VaR(p, splicefit) corresponds to qSplice(p, splicefit, lower.tail = FALSE).

Value

Vector of quantiles VaR_{1-p}=Q(1-p).

Author(s)

Tom Reynkens with R code from Roel Verbelen for the mixed Erlang quantiles.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

Reynkens, T., Verbelen, R., Beirlant, J. and Antonio, K. (2017). "Modelling Censored Losses Using Splicing: a Global Fit Strategy With Mixed Erlang and Extreme Value Distributions". Insurance: Mathematics and Economics, 77, 65–77.

Verbelen, R., Gong, L., Antonio, K., Badescu, A. and Lin, S. (2015). "Fitting Mixtures of Erlangs to Censored and Truncated Data Using the EM Algorithm." Astin Bulletin, 45, 729–758

See Also

qSplice, CTE, SpliceFit, SpliceFitPareto, SpliceFiticPareto, SpliceFitGPD

Examples

## Not run: 

# Pareto random sample
X <- rpareto(1000, shape = 2)

# Splice ME and Pareto
splicefit <- SpliceFitPareto(X, 0.6)



p <- seq(0,1,0.01)

# Plot of quantiles
plot(p, qSplice(p, splicefit), type="l", xlab="p", ylab="Q(p)")

# Plot of VaR
plot(p, VaR(p, splicefit), type="l", xlab="p", ylab=bquote(VaR[1-p]))

## End(Not run)

Weibull quantile plot

Description

Computes the empirical quantiles of the log-transform of a data vector and the theoretical quantiles of the standard Weibull distribution. These quantiles are then plotted in a Weibull QQ-plot with the theoretical quantiles on the x-axis and the empirical quantiles on the y-axis.

Usage

WeibullQQ(data, plot = TRUE, main = "Weibull QQ-plot", ...)

Arguments

Details

The Weibull QQ-plot is given by

(\log(-\log(1-i/(n+1))),\log X_{i,n})

for i=1,...,n, with X_{i,n} the i-th order statistic of the data.

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

WeibullQQ_der, ExpQQ, LognormalQQ, ParetoQQ

Examples

data(norwegianfire)

# Weibull QQ-plot for Norwegian Fire Insurance data for claims in 1976.
WeibullQQ(norwegianfire$size[norwegianfire$year==76])

# Derivative of Weibull QQ-plot for Norwegian Fire Insurance data for claims in 1976.
WeibullQQ_der(norwegianfire$size[norwegianfire$year==76])

Derivative plot of the Weibull QQ-plot

Description

Computes the derivative plot of the Weibull QQ-plot. These values can be plotted as a function of the data or as a function of the tail parameter k.

Usage

WeibullQQ_der(data, k = FALSE, plot = TRUE, 
              main = "Derivative plot of Weibull QQ-plot", ...)

Arguments

Details

The derivative plot of a Weibull QQ-plot is

(k, H_{k,n}/W_{k,n})

or

(\log X_{n-k,n}, H_{k,n}/W_{k,n})

with H_{k,n} the Hill estimates and

W_{k,n}=1/k\sum_{j=1}^k \log(\log((n+1)/j)) - \log(\log((n+1)/(k+1))).

See Section 4.1 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

See Also

WeibullQQ, Hill, MeanExcess, LognormalQQ_der, ParetoQQ_der

Examples

data(norwegianfire)

# Weibull QQ-plot for Norwegian Fire Insurance data for claims in 1976.
WeibullQQ(norwegianfire$size[norwegianfire$year==76])

# Derivative of Weibull QQ-plot for Norwegian Fire Insurance data for claims in 1976.
WeibullQQ_der(norwegianfire$size[norwegianfire$year==76])

EPD estimator for right censored data

Description

Computes the EPD estimates adapted for right censored data.

Usage

cEPD(data, censored, rho = -1, beta = NULL, logk = FALSE, 
     plot = FALSE, add = FALSE, main = "EPD estimates of the EVI", ...)

Arguments

Details

The function EPD uses \tau which is equal to -\beta.

This estimator is only suitable for right censored data.

Value

A list with following components:

Author(s)

Tom Reynkens based on R code from Anastasios Bardoutsos.

References

Beirlant, J., Bardoutsos, A., de Wet, T. and Gijbels, I. (2016). "Bias Reduced Tail Estimation for Censored Pareto Type Distributions." Statistics & Probability Letters, 109, 78–88.

Fraga Alves, M.I. , Gomes, M.I. and de Haan, L. (2003). "A New Class of Semi-parametric Estimators of the Second Order Parameter." Portugaliae Mathematica, 60, 193–214.

See Also

EPD, cProbEPD, cGPDmle

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# EPD estimator adapted for right censoring
cepd <- cEPD(Z, censored=censored, plot=TRUE)

Exponential quantile plot for right censored data

Description

Exponential QQ-plot adapted for right censored data.

Usage

cExpQQ(data, censored, plot = TRUE, main = "Exponential QQ-plot", ...)

Arguments

Details

The exponential QQ-plot adapted for right censoring is given by

( -\log(1-F_{km}(Z_{j,n})), Z_{j,n} )

for j=1,\ldots,n-1, with Z_{i,n} the i-th order statistic of the data and F_{km} the Kaplan-Meier estimator for the CDF. Hence, it has the same empirical quantiles as an ordinary exponential QQ-plot but replaces the theoretical quantiles -\log(1-j/(n+1)) by -\log(1-F_{km}(Z_{j,n})).

This QQ-plot is only suitable for right censored data.

In Beirlant et al. (2007), only a Pareto QQ-plot adapted for right-censored data is proposed. This QQ-plot is constructed using the same ideas, but is not described in the paper.

Value

A list with following components:

Author(s)

Tom Reynkens

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.

See Also

ExpQQ, cLognormalQQ, cParetoQQ, cWeibullQQ, KaplanMeier

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Exponential QQ-plot adapted for right censoring
cExpQQ(Z, censored=censored)

GPD-ML estimator for right censored data

Description

Computes ML estimates of fitting GPD to peaks over a threshold adapted for right censoring.

Usage

cGPDmle(data, censored, start = c(0.1,1), warnings = FALSE, logk = FALSE, 
        plot = FALSE, add = FALSE, main = "POT estimates of the EVI", ...)
      
cPOT(data, censored, start = c(0.1,1), warnings = FALSE, logk = FALSE, 
     plot = FALSE, add = FALSE, main = "POT estimates of the EVI", ...)

Arguments

Details

The GPD-MLE estimator for the EVI adapted for right censored data is equal to the ordinary GPD-MLE estimator for the EVI divided by the proportion of the k largest observations that is non-censored. The estimates for \sigma are the ordinary GPD-MLE estimates for \sigma.

This estimator is only suitable for right censored data.

cPOT is the same function but with a different name for compatibility with POT.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Einmahl, J.H.J., Fils-Villetard, A. and Guillou, A. (2008). "Statistics of Extremes Under Random Censoring." Bernoulli, 14, 207–227.

See Also

GPDmle, cProbGPD, cQuantGPD, cEPD

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# GPD-ML estimator adapted for right censoring
cpot <- cGPDmle(Z, censored=censored, plot=TRUE)

Hill estimator for right censored data

Description

Computes the Hill estimator for positive extreme value indices, adapted for right censoring, as a function of the tail parameter k (Beirlant et al., 2007). Optionally, these estimates are plotted as a function of k.

Usage

cHill(data, censored, logk = FALSE, plot = FALSE, add = FALSE, 
      main = "Hill estimates of the EVI", ...)

Arguments

Details

The Hill estimator adapted for right censored data is equal to the ordinary Hill estimator H_{k,n} divided by the proportion of the k largest observations that is non-censored.

This estimator is only suitable for right censored data, use icHill for interval censored data.

See Section 4.3.2 of Albrecher et al. (2017) for more details.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Albrecher, H., Beirlant, J. and Teugels, J. (2017). Reinsurance: Actuarial and Statistical Aspects, Wiley, Chichester.

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.

See Also

Hill, icHill, cParetoQQ, cProb, cQuant

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Hill estimator adapted for right censoring
chill <- cHill(Z, censored=censored, plot=TRUE)

Log-normal quantile plot for right censored data

Description

Log-normal QQ-plot adapted for right censored data.

Usage

cLognormalQQ(data, censored, plot = TRUE, main = "Log-normal QQ-plot", ...)

Arguments

Details

The log-normal QQ-plot adapted for right censoring is given by

( \Phi^{-1}(F_{km}(Z_{j,n})), \log(Z_{j,n}) )

for j=1,\ldots,n-1, with Z_{i,n} the i-th order statistic of the data, \Phi^{-1} the quantile function of the standard normal distribution and F_{km} the Kaplan-Meier estimator for the CDF. Hence, it has the same empirical quantiles as an ordinary log-normal QQ-plot but replaces the theoretical quantiles \Phi^{-1}(j/(n+1)) by \Phi^{-1}(F_{km}(Z_{j,n})).

This QQ-plot is only suitable for right censored data.

In Beirlant et al. (2007), only a Pareto QQ-plot adapted for right-censored data is proposed. This QQ-plot is constructed using the same ideas, but is not described in the paper.

Value

A list with following components:

Author(s)

Tom Reynkens

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.

See Also

LognormalQQ, cExpQQ, cParetoQQ, cWeibullQQ, KaplanMeier

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Log-normal QQ-plot adapted for right censoring
cLognormalQQ(Z, censored=censored)

MOM estimator for right censored data

Description

Computes the Method of Moment estimates adapted for right censored data.

Usage

cMoment(data, censored, logk = FALSE, plot = FALSE, add = FALSE, 
        main = "Moment estimates of the EVI", ...)

Arguments

Details

The moment estimator adapted for right censored data is equal to the ordinary moment estimator divided by the proportion of the k largest observations that is non-censored.

This estimator is only suitable for right censored data.

Value

A list with following components:

Author(s)

Tom Reynkens

References

Einmahl, J.H.J., Fils-Villetard, A. and Guillou, A. (2008). "Statistics of Extremes Under Random Censoring." Bernoulli, 14, 207–227.

See Also

Moment, cProbMOM, cQuantMOM

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Moment estimator adapted for right censoring
cmom <- cMoment(Z, censored=censored, plot=TRUE)

Pareto quantile plot for right censored data

Description

Pareto QQ-plot adapted for right censored data.

Usage

cParetoQQ(data, censored, plot = TRUE, main = "Pareto QQ-plot", ...)

Arguments

Details

The Pareto QQ-plot adapted for right censoring is given by

( -\log(1-F_{km}(Z_{j,n})), \log Z_{j,n} )

for j=1,\ldots,n-1, with Z_{i,n} the i-th order statistic of the data and F_{km} the Kaplan-Meier estimator for the CDF. Hence, it has the same empirical quantiles as an ordinary Pareto QQ-plot but replaces the theoretical quantiles -\log(1-j/(n+1)) by -\log(1-F_{km}(Z_{j,n})).

This QQ-plot is only suitable for right censored data, use icParetoQQ for interval censored data.

Value

A list with following components:

Author(s)

Tom Reynkens

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.

See Also

ParetoQQ, icParetoQQ, cExpQQ, cLognormalQQ, cWeibullQQ, cHill, KaplanMeier

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Pareto QQ-plot adapted for right censoring
cParetoQQ(Z, censored=censored)

Estimator of small exceedance probabilities and large return periods using censored Hill

Description

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.

Usage

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

Details

The probability is estimated as

\hat{P}(X>q)=(1-km) \times (q/Z_{n-k,n})^{-1/H_{k,n}^c}

with Z_{i,n} the i-th order statistic of the data, H_{k,n}^c the Hill estimator adapted for right censoring and km the Kaplan-Meier estimator for the CDF evaluated in Z_{n-k,n}.

Value

A list with following components:

Author(s)

Tom Reynkens

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.

See Also

cHill, cQuant, Prob, KaplanMeier

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# Hill estimator adapted for right censoring
chill <- cHill(Z, censored=censored, plot=TRUE)

# Small exceedance probability
q <- 10
cProb(Z, censored=censored, gamma1=chill$gamma1, q=q, plot=TRUE)

# Return period
cReturn(Z, censored=censored, gamma1=chill$gamma1, q=q, plot=TRUE)

Estimator of small exceedance probabilities and large return periods using censored EPD

Description

Computes estimates of a small exceedance probability P(X>q) or large return period 1/P(X>q) using the parameters from the EPD fit adapted for right censoring.

Usage

cProbEPD(data, censored, gamma1, kappa1, beta, q, plot = FALSE, add = FALSE,
         main = "Estimates of small exceedance probability", ...)

cReturnEPD(data, censored, gamma1, kappa1, beta, q, plot = FALSE, add = FALSE,
           main = "Estimates of large return period", ...)          

Arguments

Details

The probability is estimated as

\hat{P}(X>q)= (1-km) \times (1-F(q))

with F the CDF of the EPD with estimated parameters \hat{\gamma}_1, \hat{\kappa}_1 and \hat{\tau}=-\hat{\beta} and km the Kaplan-Meier estimator for the CDF evaluated in Z_{n-k,n} (the (k+1)-th largest data point).

Value

A list with following components:

Author(s)

Tom Reynkens.

References

Beirlant, J., Bardoutsos, A., de Wet, T. and Gijbels, I. (2016). "Bias Reduced Tail Estimation for Censored Pareto Type Distributions." Statistics & Probability Letters, 109, 78–88.

See Also

cEPD, ProbEPD, Prob, KaplanMeier

Examples

# Set seed
set.seed(29072016)

# Pareto random sample
X <- rpareto(500, shape=2)

# Censoring variable
Y <- rpareto(500, shape=1)

# Observed sample
Z <- pmin(X, Y)

# Censoring indicator
censored <- (X>Y)

# EPD estimator adapted for right censoring
cepd <- cEPD(Z, censored=censored, plot=TRUE)

# Small exceedance probability
q <- 10
cProbEPD(Z, censored=censored, gamma1=cepd$gamma1,
        kappa1=cepd$kappa1, beta=cepd$beta, q=q, plot=TRUE)

# Return period
cReturnEPD(Z, censored=censored, gamma1=cepd$gamma1,
        kappa1=cepd$kappa1, beta=cepd$beta, q=q, plot=TRUE)        

Estimator of small exceedance probabilities and large return periods using censored generalised Hill

Description

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 generalised Hill estimator adapted for right censoring.

Usage

cProbGH(data, censored, gamma1, q, plot = FALSE, add = FALSE, 
        main = "Estimates of small exceedance probability", ...)

cReturnGH(data, censored, gamma1, q, plot = FALSE, add = FALSE, 
          main = "Estimates of large return period", ...)        

Arguments