pFailure: Probability of falling into a failure region

In ExtremalDep: Extremal Dependence Models

Probability of falling into a failure region

Description

This function computes the empirical estimate of the probability of falling into two types of failure regions.

Usage

pFailure(n, beta, u1, u2, mar1, mar2, type, plot, xlab, ylab, 
         nlevels=10)

Arguments

n	An integer indicating the number of points generated to compute the empirical probability.
beta	A vector representing the Bernstein polynomial coefficients.
u1, u2	Vectors of positive reals representing some high thresholds.
mar1, mar2	Vectors of marginal (GEV) parameters
type	A character string indicating if the failure region includes at least one exceedance ("or"), or both exceednaces ("and"). If "both" then probabilities to fall in both regions are computed.
plot	A logical value; if TRUE contour plots of the probalities are displayed.
xlab, ylab	A character string equivalent representing the graphical parameters as in par.
nlevels	The number of contour levels to be displayed.

Details

The two failure regions are:

Au = {(v1,v2): v1>u1 or v2>u2}

and

Bu = {(v1,v2): v1>u1 and v2>u2}

where (v1,v2) \in R_+^2 and u1, u2>0.

Exceedances samples y11, ..., yn1 and y12, ..., yn2 are generating according to Algorithm 3 of Marcon et al. (2017) and the the estimates of the probability of falling in Au and Bu are given by

P^Au = 1/n * ∑_{i=1}^n I{yi1 >u1* or yi2 > u2*}

and

P^Bu = 1/n * ∑_{i=1}^n I{yi1 >u1* and yi2 > u2*}.

Value

A list containing AND and/or OR, depending on the type argument. Each element is a length(u1) x length(u2) matrix.

Author(s)

Simone Padoan, [e-mail], https://faculty.unibocconi.it/simonepadoan/; Boris Beranger, [e-mail] https://www.borisberanger.com;

References

Marcon, G., Naveau, P. and Padoan, S.A. (2017) A semi-parametric stochastic generator for bivariate extreme events Stat, 6, 184-201.

See Also

dExtDep, rExtDep, fExtDep, fExtDep.np

Examples

# Example wind speed and wind gust

data(WindSpeedGust)
years <- format(ParcayMeslay$time, format="%Y")
attach(ParcayMeslay[which(years %in% c(2004:2013)),])
  
WS_th <- quantile(WS,.9)
DP_th <- quantile(DP,.9)
  
pars.WS <- evd::fpot(WS, WS_th, model="pp")$estimate
pars.DP <- evd::fpot(DP, DP_th, model="pp")$estimate
  
data_uf <- trans2UFrechet(cbind(WS,DP), type="Empirical")
  
rdata <- rowSums(data_uf)
r0 <- quantile(rdata, probs=.90)
extdata <- data_uf[rdata>=r0,]
  
SP_mle <- fExtDep.np(method="Frequentist", data=extdata, k0=10, 
                     type="maxima")

if(interactive()){
  pF <- pFailure(n=50000, beta=SP_mle$Ahat$beta,
                 u1=seq(from=19, to=28, length=200), mar1=pars.WS,
                 u2=seq(from=40, to=60, length=200), mar2=pars.DP, 
                 type="both", plot=TRUE, 
                 xlab="Daily-maximum Wind Speed (m/s)",
                 ylab="Differential of Pressure (mbar)", nlevels=15)
}

ExtremalDep documentation built on March 7, 2023, 3:16 p.m.