On the Stress–Strength Reliability of Transmuted GEV Random Variables with Applications to Financial Assets Selection

Abstract

In reliability contexts, probabilities of the type $R=P(X<Y)$, where X and Y are random variables, have shown to be useful tools to compare the performance of these stochastic entities. By considering that both X and Y follow a transmuted generalized extreme-value (TGEV) distribution, new analytical relationships were derived for R in terms of special functions. The results hereby obtained are more flexible when compared to similar results found in the literature. To highlight the applicability and correctness of our results, we conducted a Monte-Carlo simulation study and investigated the use of the reliability measure $P(X<Y)$ to select among financial assets whose returns were characterized by the random variables X and Y. Our results highlight that R is an interesting alternative to modern portfolio theory, which usually relies on the contrast of involved random variables by a simple comparison of their means and standard deviations.

1. Introduction

Consider a component with a strength Y and subjected to a stress X. The component fails if the stress X exceeds the component strength Y; otherwise, it works properly. For independent components, the stress–strength reliability (SSR) R, also referred to as stress–strength probability, is given by:

$$\begin{array}{c}\hfill R=P(X<Y)={\int }_{-\infty }^{+\infty }{F}_X\left(x\right)f_Y\left(x\right)dx,\end{array}$$

(1)

where $F_X$ and $f_Y$ denote, respectively, the cumulative distribution function (CDF) of X and the probability density function (PDF) of Y.

Although R was initially applied in the context of engineering, the interest in such metric spreads to several areas, such as household financial fragility [1], stock marketmodeling [2], asset selection [3], among others. We refer the reader to [4] for further details on stress–strength models.

The choice of an appropriate distribution to model both X and Y directly influences the calculation and estimation of R. In Finance, we have strong evidence that asset returns are better modeled by either $\alpha$-stable processes (heavy-tailed alternative to Brownian motion [5]) or by heavy-tailed time series models [6,7]. The Extreme-Value Theory (EVT) made available a body of knowledge around heavy-tailed distributions, like the definition of the extreme-value distributions [8] (and its variations), which can be used as a proxy of various fat-tailed distributions. Several studies have successfully applied EVT to model financial data [3,9,10,11,12], where it has been shown that EVT-based models may provide adequate risk management strategies. Therefore, generalizations of extreme-value distributions may enhance the quality of the models, and this is the general contribution of the present work.

In this paper, we are interested in further exploring R calculations in the context of extreme-value distributions. The estimation of R, when X and Y are independent random variables following extreme-value distributions, has been extensively studied. For example, Ref. [13] derived the expression of R for the extreme distributions Gumbel, Fréchet, and Weibull, Ref. [14] considered a Bayesian analysis of the Fréchet stress–strength model, Ref. [15] discussed Bayesian estimation of $P(Y<X)$ for the Weibull distribution with arbitrary parameters and [16] improved the estimation for R by not using transformations in the data and eliminating the constraints on the parameters in the case of the Weibull models. Closed-form expressions for R when X and Y follow generalized extreme-value (GEV) distributions were obtained in [3], who also proposed an estimation procedure for R by not using transformations in the data and with as few parameter restrictions as possible.

Several generalizations of the extreme-value distribution have been proposed, but in the present paper, the so-called transmuted generalized extreme-value (TGEV) distribution shall be considered. The TGEV distribution, initially proposed by [17], has since been extensively studied and applied in various modeling scenarios. Significant contributions to its application and understanding have been made by [18,19]. Essentially, the TGEV distribution is a modification of the generalized extreme-value distribution (GEV), whose CDF is given by:

$$G\left(x\right)=\left\{\begin{array}{cc}exp\left(-{(1+\gamma \frac{x-\mu }{\sigma })}^{-1/\gamma }\right),& 1+\gamma \left(\frac{x-\mu }{\sigma }\right)>0\mathrm{and}\gamma \ne 0,\\ exp\left(-exp\left(-\frac{x-\mu }{\sigma }\right)\right),& x\in \mathbb{R}\mathrm{and}\gamma =0,\end{array}\right.$$

(2)

where $\gamma >0$ is the shape parameter, $\mu \in \mathbb{R}$ is the location parameter and $\sigma >0$ is the scale parameter. Then, the TGEV distribution is obtained as follows: given the GEV distribution $G\left(x\right)$, the transmuted distribution function F is given by:

$$\begin{array}{c}\hfill F\left(x\right)=(1+\lambda )G\left(x\right)-\lambda {\left[G\left(x\right)\right]}^2,\left|\lambda \right|<1.\end{array}$$

(3)

Properties such as moments, quartiles, tail behavior, and order statistics, among others, were studied in [19]. They also showed its applicability in modeling log-returns of stock prices.

In [19], the TGEV parameters $(\mu ,\sigma ,\gamma ,\lambda )$ were estimated by a maximum likelihood approach. In contrast, this work proposes a two-step estimation procedure. First, a GEV model is estimated to yield $(\mu ,\sigma ,\gamma )$ parameters. Then, a refinement step is taken by estimating the $\lambda$ parameter in an attempt to improve the first step fit and to reduce the overall computational effort to estimate the TGEV parameters.

In this paper, we consider the problem of estimating the stress–strength parameter R when X and Y are independent TGEV random variables. In addition, our framework does not require transformations in the data and, to the best of our knowledge, allows for the fewest parameter restrictions.

Our main contributions are

• to derive R analytically in terms of special functions;
• to derive closed-form expressions for multicomponent system reliability;
• to propose an estimation procedure for R and validate such procedure via a simulation study and
• to apply the theoretical results in asset selection problems in finance.

The paper is organized as follows: in Section 2, we define the extreme-value $\mathbb{H}$-function and the H-function, and we explicitly present the CDF and PDF of the TGEV distribution. Section 3 deals with the derivation of R when X and Y are independent TGEV random variables. The maximum likelihood estimation for R is presented in Section 4. In Section 5, we discuss a simulation study and a stock price modeling application for asset selection. The last section presents the conclusions.

2. Preliminaries

In this section, we give some definitions and results which will be used subsequently.

2.1. Special Functions

Recently, the extreme-value $\mathbb{H}$-function was introduced in [20]. This function is defined as follows:

$$\mathbb{H}(a_1,a_2,a_3,a_4,a_5,a_6):={\int }_0^{\infty }{y}^{a_6}exp(-a_{1}y-{(a_2{y}^{a_3}+a_4)}^{a_5})dy,$$

(4)

where $\Re \left(a_1\right),\Re \left(a_2\right),\Re \left(a_4\right)\in {\mathbb{R}}_{+}$, and $a_3,a_5\in \mathbb{C}$. It is important to note that both $\Re \left(a_1\right)$ and $\Re \left(a_2\right)$ cannot be equal to zero simultaneously. Moreover, $\Re \left(a_6\right)>-1$ when $a_1\ne 0$, or when $a_1=0$ and $\mathrm{sign}\left(a_3\right)=\mathrm{sign}\left(a_5\right)$. Conversely, $\Re \left(a_6\right)<-1$ when $a_1=0$ and $\mathrm{sign}\left(a_3\right)\ne \mathrm{sign}\left(a_5\right)$. In this context, $\mathbb{R}$, $\mathbb{C}$, and ℜ, respectively, denote the set of real numbers, complex numbers, and the real part of a complex number.

Another important special function is the H-function, which can be defined by:

$$H_{p,q}^{m,n}\left[z|\begin{array}{c}(a_1,A_1),\cdots ,(a_p,A_p)\\ (b_1,B_1),\cdots ,(b_q,B_q)\end{array}\right]=\frac{1}{2\pi i}{\int }_L\frac{{\prod }_{k=1}^m\Gamma (b_j+B_{j}s){\prod }_{j=1}^n\Gamma (1-a_j-A_{j}s)}{{\prod }_{k=m+1}^q\Gamma (1-b_j-B_{j}s){\prod }_{j=n+1}^p\Gamma (a_j+A_{j}s)}{z}^{-s}ds,$$

where $0\le m\le q$, $0\le n\le p$ (not both m and n zeros simultaneously), $A_j>0$ ($j=1,\cdots ,p$), $B_k>0$ ($k=1,\cdots ,q$), $a_j$ and $b_k$ are complex numbers such that no poles of $\Gamma (b_k+B_{k}s)$ ($k=1,\cdots ,m)$ coincide with poles of $\Gamma (1-a_j-A_{j}s)$ ($j=1,\cdots ,n)$. L is a suitable contour $w-i\infty$ to $w+i\infty$, $w\in \mathbb{R}$, separating the poles of the two types mentioned above. For more details, see [21]. As special cases, we have:

$${\int }_0^{\infty }exp\{-ay-b{y}^c\}dy=\frac{1}{b^{1/c}c}{H}_{1,1}^{1,1}\left[\frac{a}{b^{1/c}}|\begin{array}{c}(\frac{c-1}{c},\frac{1}{c})\\ (0,1)\end{array}\right]$$

(5)

for $a>0$, $b>0$ and $c>0$ and

$${\int }_0^{\infty }{y}^{d-1}exp\{-ay-b{y}^c\}dy=\frac{1}{a^d}{H}_{1,1}^{1,1}\left[b{a}^{-c}|\begin{array}{c}(1-d,c)\\ (0,1)\end{array}\right],$$

(6)

where $a>0$, $b>0$, $c>0$ and $d>0$.

2.2. Transmuted GEV Distribution

The CDF and the PDF of the TGEV distribution are given, respectively, by:

$$F(x;\mu ,\sigma ,\gamma ,\lambda )=\left\{\begin{array}{cc}\hfill exp(-w^{-1/\gamma })\left[(1+\lambda )-\lambda exp(-w^{-1/\gamma })\right],& \gamma \ne 0;\hfill \\ \hfill exp\left(-exp\left(-\left(\frac{w-1}{\gamma }\right)\right)\right)\left[(1+\lambda )-\lambda exp\left(-\left(\frac{w-1}{\gamma }\right)\right)\right],& \gamma =0,\hfill \end{array}\right.$$

and

$$f(x;\mu ,\sigma ,\gamma ,\lambda )=\left\{\begin{array}{cc}\hfill \left[\frac{{\left(w\right)}^{-1-1/\gamma }exp(-w^{-1/\gamma })}{\sigma }\right]\left[(1+\lambda )-2\lambda exp(-w^{-1/\gamma })\right],& \gamma \ne 0;\hfill \\ \hfill \frac{exp(-\left(\frac{w-1}{\gamma }\right))exp\left(-exp(-\left(\frac{w-1}{\gamma }\right)\right)}{\sigma }\left[(1+\lambda )-2\lambda exp\left(-exp\left(-\left(\frac{w-1}{\gamma }\right)\right)\right)\right],& \gamma =0,\hfill \end{array}\right.$$

(7)

where $w=1+\left(\gamma \right(x-\mu )/\sigma )$ and $supp\left(F\right)=supp\left(f\right)=\{x;w>0\}$.

Let $X\sim TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and $Y\sim TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$ be independent random variables. The stress–strength probability is given by:

$$\begin{array}{c}\hfill R=P(X<Y)={\int }_{-\infty }^{+\infty }{F}_X(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)f_Y(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du.\end{array}$$

In Section 5, we will apply the density f to modeling log-return stock prices. Furthermore, the support of f depends on the parameters. Thus, the maximum likelihood estimation is not as straightforward as in the usual cases. Figure 1 shows the behavior of the PDF of TGEV random variables for some choices of parameters. Proper parameter selection can be explored to represent the extremal models as particular cases of TGEV distribution, as shown in

Table 1. Particular cases of TGEV model.

Distribution	CDF
GEV	$F(x;0,1,\gamma ,0)$
Fréchet	$F\left(\right(x-1)/\gamma ;0,1,\gamma ,0)$, $\gamma >0$
Reversed Weibull	$F(-(1+x)/\gamma ;0,1,\gamma ,0)$, $\gamma <0$
Gumbel	$F(x;0,1,0,0)$

.

3. Main Results

In this section, the reliability of two independent TGEV random variables is derived in terms of $\mathbb{H}$-functions. In addition, with suitable parameter restrictions, simpler expressions in terms of the H-function are also obtained. First, we consider the case of two independent TGEV with $\mathrm{sign}\left({\gamma }_x\right)=\mathrm{sign}\left({\gamma }_y\right)\ne 0$.

Theorem 1.

Let X and Y be independent random variables, respectively, with distribution $TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and $TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$, ${\mu }_j\in \mathbb{R}$, ${\sigma }_j\in {\mathbb{R}}_{+}$, ${\gamma }_j\in \mathbb{R}$$({\gamma }_j\ne 0)$, ${\lambda }_j\in [-1,1]$, $j\in \{x,y\}$. Then

• When ${\gamma }_j>0$, $j\in \{x,y\}$:

$$\begin{array}{ccc}\hfill R=P(X<Y)& =& (1+{\lambda }_x)(1+{\lambda }_y)\mathbb{H}\left(1,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right),-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & -2{\lambda }_y(1+{\lambda }_x)\mathbb{H}\left(2,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right),-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & -{\lambda }_x(1+{\lambda }_y)\mathbb{H}\left(1,2^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,2^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right],-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & +2{\lambda }_x{\lambda }_y\mathbb{H}\left(2,2^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,2^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right],-\frac{1}{{\gamma }_x},0\right),\hfill \end{array}$$

(8)

• provided that ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}\ge {\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$. When ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}<{\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$:

$$\begin{array}{ccc}\hfill R=P(X<Y)& =& 1-(1+{\lambda }_y)(1+{\lambda }_x)\mathbb{H}\left(1,\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right),-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & +2{\lambda }_x(1+{\lambda }_y)\mathbb{H}\left(2,\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right),-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & +{\lambda }_y(1+{\lambda }_x)\mathbb{H}\left(1,2^{-{\gamma }_y}\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,2^{-{\gamma }_y}\left[1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & -2{\lambda }_y{\lambda }_x\mathbb{H}\left(2,2^{-{\gamma }_y}\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,2^{-{\gamma }_y}\left[1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_y},0\right).\hfill \end{array}$$

(9)

• When ${\gamma }_j<0$, $j\in \{x,y\}$:

$$\begin{array}{ccc}\hfill R=P(X<Y)& =& (1+{\lambda }_x)(1+{\lambda }_y)\mathbb{H}\left(1,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right),-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & -2{\lambda }_y(1+{\lambda }_x)\mathbb{H}\left(2,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right),-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & -{\lambda }_x(1+{\lambda }_y)\mathbb{H}\left(1,2^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,2^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right],-\frac{1}{{\gamma }_x},0\right)\hfill \\ \hfill & & +2{\lambda }_x{\lambda }_y\mathbb{H}\left(2,2^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,2^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right],-\frac{1}{{\gamma }_x},0\right),\hfill \end{array}$$

(10)

provided that ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}<{\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$. When ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}\ge {\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$:

$$\begin{array}{ccc}\hfill R=P(X<Y)& =& 1-(1+{\lambda }_y)(1+{\lambda }_x)\mathbb{H}\left(1,\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right),-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & +2{\lambda }_x(1+{\lambda }_y)\mathbb{H}\left(2,\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right),-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & +{\lambda }_y(1+{\lambda }_x)\mathbb{H}\left(1,2^{-{\gamma }_y}\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,2^{-{\gamma }_y}\left[1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_y},0\right)\hfill \\ \hfill & & -2{\lambda }_y{\lambda }_x\mathbb{H}\left(2,2^{-{\gamma }_y}\frac{{\gamma }_y{\sigma }_x}{{\gamma }_x{\sigma }_y},-{\gamma }_x,2^{-{\gamma }_y}\left[1+\frac{{\gamma }_y}{{\sigma }_y}\left({\mu }_x-{\mu }_y-\frac{{\sigma }_x}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_y},0\right).\hfill \end{array}$$

(11)

Proof. 

Set $S=supp{F}_X\cap supp{f}_Y$. Then

$$S=\left\{\begin{array}{cc}(M,\infty ),& \mathrm{if}sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=1,\\ (-\infty ,m),& \mathrm{if}sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=-1,\\ \mathbb{R},& \mathrm{if}sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=0,\end{array}\right.$$

(12)

where $M=max\{{\mu }_x-{\sigma }_x/{\gamma }_x,{\mu }_y-{\sigma }_y/{\gamma }_y\}$ and $m=min\{{\mu }_x-{\sigma }_x/{\gamma }_x,{\mu }_y-{\sigma }_y/{\gamma }_y\}$.

Note that

$$\begin{array}{ccc}\hfill R=P(X<Y)& =& {\int }_{-\infty }^{\infty }{F}_X(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)f_Y(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du\hfill \\ \hfill & =& {\int }_S\left(exp(-w_x^{-1/{\gamma }_x})\left[(1+{\lambda }_x)-{\lambda }_{x}exp(-w_x^{-1/{\gamma }_x})\right]\right.\hfill \\ \hfill & & \times \left.\left[\frac{{\left(w_y\right)}^{-1-1/{\gamma }_y}exp(-w_y^{-1/{\gamma }_y})}{{\sigma }_y}\right]\left[(1+{\lambda }_y)-2{\lambda }_{y}exp(-w_y^{-1/{\gamma }_y})\right]\right)du,\hfill \end{array}$$

(13)

where $w_i=1+\frac{{\gamma }_i}{{\sigma }_i}(u-{\mu }_i)$, $i\in \{x,y\}$. We have four cases to consider:

• ${\gamma }_x>0$ and ${\gamma }_y>0$

  (a)

  ${\mu }_y-{\sigma }_y/{\gamma }_y\ge {\mu }_x-{\sigma }_x/{\gamma }_x$;

  (b)

  ${\mu }_y-{\sigma }_y/{\gamma }_y<{\mu }_x-{\sigma }_x/{\gamma }_x$;

• ${\gamma }_x<0$ and ${\gamma }_y<0$

  (a)

  ${\mu }_y-{\sigma }_y/{\gamma }_y<{\mu }_x-{\sigma }_x/{\gamma }_x$;

  (b)

  ${\mu }_y-{\sigma }_y/{\gamma }_y\ge {\mu }_x-{\sigma }_x/{\gamma }_x$.

Let us consider case 1(a). Substituting $v=w_y^{-1/{\gamma }_y}$, it follows from (13) that

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y){\int }_0^{\infty }exp\left(-v-{\left[\frac{{\gamma }_x{\sigma }_y}{{\sigma }_x{\gamma }_y}{v}^{-{\gamma }_y}+1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right]}^{-1/{\gamma }_x}\right)dv\hfill \\ \hfill & & -2(1+{\lambda }_x){\lambda }_y{\int }_0^{\infty }exp\left(-2v-{\left[\frac{{\gamma }_x{\sigma }_y}{{\sigma }_x{\gamma }_y}{v}^{-{\gamma }_y}+1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right]}^{-1/{\gamma }_x}\right)dv\hfill \\ \hfill & & -{\lambda }_x(1+{\lambda }_y){\int }_0^{\infty }exp\left(-v-2{\left[\frac{{\gamma }_x{\sigma }_y}{{\sigma }_x{\gamma }_y}{v}^{-{\gamma }_y}+1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right]}^{-1/{\gamma }_x}\right)dv\hfill \\ \hfill & & +2{\lambda }_x{\lambda }_y{\int }_0^{\infty }exp\left(-2v-2{\left[\frac{{\gamma }_x{\sigma }_y}{{\sigma }_x{\gamma }_y}{v}^{-{\gamma }_y}+1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right]}^{-1/{\gamma }_x}\right)dv.\hfill \end{array}$$

(14)

Therefore, (8) follows from (4) and (14). For case 1(b), it suffices to notice that $P(X<Y)=1-P(Y<X)$ and apply the result in (8) with interchanged sub-indices. For cases 2(a) and 2(b), the same rationale can be applied, just noticing that in such cases, x mostly takes negative values.    □

Remark 1.

Note that if we take ${\lambda }_x={\lambda }_y=0$, X and Y are random variables with GEV distributions, then our Theorem 1 generalizes the Theorem 3.1 in [3].

Remark 2.

In a practical scenario, the estimates $({\widehat{\mu }}_x,{\widehat{\sigma }}_x,{\widehat{\gamma }}_x,{\widehat{\lambda }}_x,{\widehat{\mu }}_y,{\widehat{\sigma }}_y,{\widehat{\gamma }}_y,{\widehat{\lambda }}_y)$ should be obtained. Then, if $sign\left({\widehat{\gamma }}_x\right)=sign\left({\widehat{\gamma }}_y\right)\ne 0$, the conditions ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}\ge {\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$ or ${\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}<{\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}$ must be verified and the corresponding R expression should be used.

Remark 3.

It follows from (5) that if ${\mu }_x-{\sigma }_x/{\gamma }_x={\mu }_y-{\sigma }_y/{\gamma }_y$ and $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)\ne 0$, then (8) can be written in terms of H-function as:

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y)\frac{{\gamma }_x}{{\gamma }_y}{\left(\frac{{\sigma }_y{\gamma }_x}{{\sigma }_x{\gamma }_y}\right)}^{1/{\gamma }_y}{H}_{1,1}^{1,1}\left[{\left(\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}\right)}^{1/{\gamma }_y}|\begin{array}{c}(1-\frac{{\gamma }_x}{{\gamma }_y},\frac{{\gamma }_x}{{\gamma }_y})\\ (0,1)\end{array}\right]\hfill \\ \hfill & -& 2{\lambda }_y(1+{\lambda }_x)\frac{{\gamma }_x}{{\gamma }_y}{\left(\frac{{\sigma }_y{\gamma }_x}{{\sigma }_x{\gamma }_y}\right)}^{1/{\gamma }_y}{H}_{1,1}^{1,1}\left[2{\left(\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}\right)}^{1/{\gamma }_y}|\begin{array}{c}(1-\frac{{\gamma }_x}{{\gamma }_y},\frac{{\gamma }_x}{{\gamma }_y})\\ (0,1)\end{array}\right]\hfill \\ \hfill & -& {\lambda }_x(1+{\lambda }_y)\frac{{\gamma }_x}{{\gamma }_y}{2}^{-{\gamma }_x/{\gamma }_y}{\left(\frac{{\sigma }_y{\gamma }_x}{{\sigma }_x{\gamma }_y}\right)}^{1/{\gamma }_y}{H}_{1,1}^{1,1}\left[2^{-{\gamma }_x/{\gamma }_y}{\left(\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}\right)}^{1/{\gamma }_y}|\begin{array}{c}(1-\frac{{\gamma }_x}{{\gamma }_y},\frac{{\gamma }_x}{{\gamma }_y})\\ (0,1)\end{array}\right]\hfill \\ \hfill & +& {\lambda }_x{\lambda }_y\frac{{\gamma }_x}{{\gamma }_y}{2}^{1-{\gamma }_x/{\gamma }_y}{\left(\frac{{\sigma }_y{\gamma }_x}{{\sigma }_x{\gamma }_y}\right)}^{1/{\gamma }_y}{H}_{1,1}^{1,1}\left[2^{1-{\gamma }_x/{\gamma }_y}{\left(\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}\right)}^{1/{\gamma }_y}|\begin{array}{c}(1-\frac{{\gamma }_x}{{\gamma }_y},\frac{{\gamma }_x}{{\gamma }_y})\\ (0,1)\end{array}\right].\hfill \end{array}$$

(15)

In particular, by using a special case of the H-function as seen in [21], if ${\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}={\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}$, $\frac{{\sigma }_y}{{\sigma }_x}<2^{-{\gamma }_x}$ and ${\gamma }_x={\gamma }_y\ne 0$, we have:

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y){\left[1+{\left(\frac{{\sigma }_y}{{\sigma }_x}\right)}^{-1/{\gamma }_x}\right]}^{-1}-2(1+{\lambda }_x){\lambda }_y{\left[2+{\left(\frac{{\sigma }_y}{{\sigma }_x}\right)}^{-1/{\gamma }_x}\right]}^{-1}\hfill \\ \hfill & & -{\lambda }_x(1+{\lambda }_y){\left[1+2{\left(\frac{{\sigma }_y}{{\sigma }_x}\right)}^{-1/{\gamma }_x}\right]}^{-1}+{\lambda }_x{\lambda }_y{\left[1+{\left(\frac{{\sigma }_y}{{\sigma }_x}\right)}^{-1/{\gamma }_x}\right]}^{-1}.\hfill \end{array}$$

(16)

Lastly, we consider the cases of two independent TGEV distributions with ${\gamma }_x={\gamma }_y=0$.

Theorem 2.

Let $X\sim TGEV({\mu }_x,{\sigma }_x,0,{\lambda }_x)$ and $Y\sim TGEV({\mu }_y,{\sigma }_y,0,{\lambda }_y)$ be independent random variables with ${\mu }_j\in \mathbb{R}$, ${\sigma }_j>0$, ${\lambda }_j\in [-1,1]$, $j\in \{x,y\}$. Then

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y)\mathbb{H}\left(1,exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,0\right)\hfill \\ \hfill & -& 2(1+{\lambda }_x){\lambda }_y\mathbb{H}\left(2,exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,0\right)\hfill \\ \hfill & -& {\lambda }_x(1+{\lambda }_y)exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)\mathbb{H}\left(1,exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,\frac{{\sigma }_y}{{\sigma }_x}\right)\hfill \\ \hfill & +& 2{\lambda }_x{\lambda }_{y}exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)\mathbb{H}\left(2,exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,\frac{{\sigma }_y}{{\sigma }_x}\right).\hfill \end{array}$$

(17)

In particular, if we take ${\sigma }_x={\sigma }_y$, we obtain the explicit form

$$R=\frac{(1+{\lambda }_x)(1+{\lambda }_y)}{1+exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)}-\frac{2(1+{\lambda }_x){\lambda }_y}{2+exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)}-\frac{{\lambda }_x(1+{\lambda }_y)exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)}{{\left(1+exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)\right)}^2}+\frac{2{\lambda }_x{\lambda }_{y}exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)}{{\left(2+exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)\right)}^2}.$$

Proof. 

Denote $F_X$ and $f_Y$, respectively, the CDF and PDF function of X and Y. Then

$$\begin{array}{ccc}\hfill R& =& {\int }_{-\infty }^{+\infty }{F}_X(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)f_Y(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du\hfill \\ \hfill & =& {\int }_{-\infty }^{+\infty }exp\left(-exp\left(-\frac{u-{\mu }_x}{{\sigma }_x}\right)\right)\left[(1+{\lambda }_x)-{\lambda }_{x}exp\left(-\frac{u-{\mu }_x}{{\sigma }_x}\right)\right]\hfill \\ \hfill & \times & \frac{exp\left(-\frac{u-{\mu }_y}{{\sigma }_y}\right)exp\left(-exp\left(-\frac{u-{\mu }_y}{{\sigma }_y}\right)\right)}{{\sigma }_y}\left[(1+{\lambda }_y)-2{\lambda }_{y}exp\left(-exp\left(-\frac{u-{\mu }_y}{{\sigma }_y}\right)\right)\right]du.\hfill \end{array}$$

(18)

Substituting $v=exp\left(-\frac{u-{\mu }_y}{{\sigma }_y}\right)$, we can rewrite (18) as

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y){\int }_0^{\infty }exp\left(-v-exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)v^{{\sigma }_y/{\sigma }_x}\right)dv\hfill \\ \hfill & -& 2(1+{\lambda }_x){\lambda }_y{\int }_0^{\infty }exp\left(-2v-exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)v^{{\sigma }_y/{\sigma }_x}\right)dv\hfill \\ \hfill & -& {\lambda }_x(1+{\lambda }_y)exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right){\int }_0^{\infty }{v}^{{\sigma }_y/{\sigma }_x}exp\left(-v-exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)v^{{\sigma }_y/{\sigma }_x}\right)dy\hfill \\ \hfill & +& 2{\lambda }_x{\lambda }_{y}exp\left(-\frac{{\mu }_y-{\mu }_x}{{\sigma }_x}\right){\int }_0^{\infty }{v}^{{\sigma }_y/{\sigma }_x}exp\left(-2v-exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)v^{{\sigma }_y/{\sigma }_x}\right)dv.\hfill \end{array}$$

(19)

Hence, (17) follows from (4) and (19).    □

Remark 4.

It follows from (5) and (6) that (17) can be rewritten in terms of H-function as

$$\begin{array}{ccc}\hfill R& =& (1+{\lambda }_x)(1+{\lambda }_y)\frac{{\sigma }_x}{{\sigma }_y}exp\left(\frac{{\mu }_y-{\mu }_x}{{\sigma }_y}\right)H_{1,1}^{1,1}\left[exp\left(\frac{{\mu }_y-{\mu }_x}{{\sigma }_y}\right)|\begin{array}{c}(\frac{{\sigma }_y-{\sigma }_x}{{\sigma }_y},\frac{{\sigma }_x}{{\sigma }_y})\\ (0,1)\end{array}\right]\hfill \\ \hfill & -& 2(1+{\lambda }_x){\lambda }_y\frac{{\sigma }_x}{{\sigma }_y}exp\left(\frac{{\mu }_y-{\mu }_x}{{\sigma }_y}\right)H_{1,1}^{1,1}\left[2exp\left(\frac{{\mu }_y-{\mu }_x}{{\sigma }_y}\right)|\begin{array}{c}(\frac{{\sigma }_y-{\sigma }_x}{{\sigma }_y},\frac{{\sigma }_x}{{\sigma }_y})\\ (0,1)\end{array}\right]\hfill \\ \hfill & -& {\lambda }_x(1+{\lambda }_y)exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)H_{1,1}^{1,1}\left[exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)|\begin{array}{c}(-\frac{{\sigma }_y}{{\sigma }_x},\frac{{\sigma }_y}{{\sigma }_x})\\ (0,1)\end{array}\right]\hfill \\ \hfill & +& {\lambda }_x{\lambda }_{y}exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)2^{-\frac{{\sigma }_y}{{\sigma }_x}}{H}_{1,1}^{1,1}\left[2^{-\frac{{\sigma }_y}{{\sigma }_x}}exp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right)|\begin{array}{c}(-\frac{{\sigma }_y}{{\sigma }_x},\frac{{\sigma }_y}{{\sigma }_x})\\ (0,1)\end{array}\right].\hfill \end{array}$$

(20)

Multicomponent System Reliability

Let $X_1,\cdots ,X_n$ be independent and identically distributed random variables with distribution $TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and Y be an independent random variable with distribution $TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$. Set $M_n=max\{X_1,\cdots ,X_n\}$. Then,

$$P(M_n\le u)=F^n(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$$

and we have

$$P(X_1<Y,\cdots ,X_n<Y)=P(M_n\le Y)={\int }_{-\infty }^{\infty }{F}^n(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)f(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du=:I_n.$$

(21)

In a broader context, consider independent random variables $Y,X_1,\cdots ,X_k$ with

$$Y\sim TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)and{X}_j\sim TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x),j=1,\cdots ,k.$$

The reliability in a multicomponent stress–strength model is given by

$$\begin{array}{ccc}\hfill R_{s,k}& =& P\left(\mathrm{at} \mathrm{least} \mathrm{s}\mathrm{of}(X_1,\cdots ,X_k)\mathrm{exceed}Y\right)\hfill \\ \hfill & =& \sum _{j=s}^k\left(\begin{array}{c}k\\ j\end{array}\right){\int }_{-\infty }^{\infty }{\left(1-F(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)\right)}^j{\left(F(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)\right)}^{k-j}f(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du.\hfill \end{array}$$

Using a binomial expansion, we obtain

$$R_{s,k}=\sum _{j=s}^k\sum _{r=0}^j\left(\begin{array}{c}k\\ j\end{array}\right)\left(\begin{array}{c}j\\ r\end{array}\right){(-1)}^{j-r}{\int }_{-\infty }^{\infty }{\left(F(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)\right)}^{k-r}f(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)du.$$

(22)

Note that the integral terms in (22) is the same as (21) when $n=k-r$. Therefore,

$$R_{s,k}=\sum _{j=s}^k\sum _{r=0}^j\left(\begin{array}{c}k\\ j\end{array}\right)\left(\begin{array}{c}j\\ r\end{array}\right){(-1)}^{j-r}{I}_{k-r}.$$

Closed expressions for (21) are presented below.

Theorem 3.

Let $X_1,\cdots ,X_n$ be independent and identically distributed random variables with distribution $TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and Y be an independent random variable with distribution $TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$. Then

• When $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=1$:

$$I_n=\sum _{l=0}^n\left(\begin{array}{c}n\hfill \\ l\hfill \end{array}\right){(-1)}^{n-l}{\left(1+{\lambda }_x\right)}^l{\lambda }_x^{n-l}\left[(1+{\lambda }_y){\tilde{I}}_{n,l}-2{\lambda }_y{\widehat{I}}_{n,l}\right],$$

(23)

where

$${\tilde{I}}_{n,l}=\mathbb{H}\left(1,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}{(2n-l)}^{-{\gamma }_x},-{\gamma }_y,\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y}\right)\right]{(2n-l)}^{-{\gamma }_x},-\frac{1}{{\gamma }_x},0\right)$$

and

$${\widehat{I}}_{n,l}=\mathbb{H}\left(2,{(2n-l)}^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,{(2n-l)}^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_x},0\right),$$

provided that ${\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}\le {\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}$.

• When $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=-1$, (23) holds provided that ${\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}\ge {\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}$.
• When ${\gamma }_x={\gamma }_y=0$:

$$I_n=\sum _{l=0}^n\left(\begin{array}{c}n\hfill \\ l\hfill \end{array}\right){(-1)}^{n-l}{\left(1+{\lambda }_x\right)}^l{\lambda }_x^{n-l}exp\left(\frac{(n-l)({\mu }_x-{\mu }_y)}{{\sigma }_x}\right)\left[(1+{\lambda }_y){\tilde{J}}_{n,l}-2{\lambda }_y{\widehat{J}}_{n,l}\right],$$

(24)

where

$${\tilde{J}}_{n,l}=\mathbb{H}\left(1,nexp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,\frac{(n-l){\sigma }_y}{{\sigma }_x}\right)$$

and

$${\widehat{J}}_{n,l}=\mathbb{H}\left(2,nexp\left(\frac{{\mu }_x-{\mu }_y}{{\sigma }_x}\right),\frac{{\sigma }_y}{{\sigma }_x},0,1,\frac{(n-l){\sigma }_y}{{\sigma }_x}\right).$$

Proof. 

For simplicity of notations, denote

$$\begin{array}{cc}& F_x\left(u\right)=F\left(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x\right),\hfill \\ & F_y\left(u\right)=F\left(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y\right),\hfill \\ & G_x\left(u\right)=G\left(u;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x\right),\hfill \end{array}$$

and

$$G_y\left(u\right)=G\left(u;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y\right).$$

It follows from (3) and (21) that

$$\begin{array}{cc}\hfill I_n& ={\int }_{-\infty }^{+\infty }{F}_x^n\left(u\right)f_y\left(u\right)du\hfill \\ & ={\int }_{-\infty }^{+\infty }{\left[\left(1+{\lambda }_x\right)G_x\left(u\right)-{\lambda }_x{G}_x^2\left(u\right)\right]}^n{f}_y\left(u\right)du.\hfill \end{array}$$

By binomial expansion

$$I_n=\sum _{l=0}^n\left(\begin{array}{c}n\hfill \\ l\hfill \end{array}\right){(-1)}^{n-l}{\left(1+{\lambda }_x\right)}^l{\lambda }_x^{n-l}{\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{f}_y\left(u\right)du.$$

(25)

Observe that $f_y\left(u\right)=F_y^{\prime }\left(u\right)$, which implies

$$\begin{array}{ccc}\hfill {\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{f}_y\left(u\right)du& =& {\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}\left[\left(1+{\lambda }_y\right)g_y\left(u\right)-2{\lambda }_y{G}_y\left(u\right)g_y\left(u\right)\right]du\hfill \\ \hfill & =& \left(1+{\lambda }_y\right){\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{g}_y\left(u\right)du\hfill \\ \hfill & & -2{\lambda }_y{\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{G}_y\left(u\right)g_y\left(u\right)du.\hfill \end{array}$$

(26)

If $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=1$, it follows from (3.15) in [3] that

$$\begin{array}{ccc}\hfill & & {\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{g}_y\left(u\right)du\hfill \\ \hfill & & =\mathbb{H}\left(1,\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x}{(2n-l)}^{-{\gamma }_x},-{\gamma }_y,\left[1+\frac{{\gamma }_x}{{\sigma }_x}({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_y})\right]{(2n-l)}^{-{\gamma }_x},-\frac{1}{{\gamma }_x},0\right),\hfill \end{array}$$

(27)

provided that ${\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}\le {\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}$. If $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=-1$, (3.17) in [3] implies (27) since ${\mu }_x-\frac{{\sigma }_x}{{\gamma }_x}\ge {\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}$.

Observe that the integration range can be simplified using the results for the intersection of the supports of $G_x$ and $g_y$, such that:

$$S=supp{G}_x\cap supp{g}_y=\left\{\begin{array}{cc}\left(M,+\infty \right),\hfill & {\gamma }_x>0\mathrm{and}{\gamma }_y>0,\hfill \\ (-\infty ,m),\hfill & {\gamma }_x<0\mathrm{and}{\gamma }_y<0,\hfill \\ \mathbb{R},\hfill & {\gamma }_x={\gamma }_y=0,\hfill \end{array}\right.$$

where $M=max\{{\mu }_x-\frac{{\sigma }_x}{{\gamma }_y},{\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}\}$ and $m=min\{{\mu }_x-\frac{{\sigma }_x}{{\gamma }_y},{\mu }_y-\frac{{\sigma }_y}{{\gamma }_y}\}$. Then, if ${\gamma }_x>0$ and ${\gamma }_y>0$ (case ${\gamma }_y<0$ and ${\gamma }_y<0$ is analogous), we have that $\tilde{I}:={\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{G}_y\left(u\right)g_y\left(u\right)du$ is given by

$$\begin{array}{ccc}\hfill \tilde{I}& =& {\int }_M^{+\infty }exp\left\{-(2n-l){[1+\frac{{\gamma }_x}{{\sigma }_x}\left(u-{\mu }_x\right)]}^{-1/{\gamma }_x}-2{[1+\frac{{\gamma }_y}{{\sigma }_y}\left(u-{\mu }_y\right)]}^{-1/{\gamma }_y}\right\}\hfill \\ & & \times {[1+\frac{{\gamma }_y}{{\sigma }_y}\left(u-{\mu }_y\right)]}^{-1/{\gamma }_y-1}\frac{du}{{\sigma }_y}.\hfill \end{array}$$

Substituting $v={\left[1+\frac{{\gamma }_y}{{\sigma }_y}\left(u-{\mu }_y\right)\right]}^{-1/{\gamma }_y}$, we obtain

$$\begin{array}{ccc}\hfill & & -2{\lambda }_y{\int }_{-\infty }^{+\infty }{G}_x{\left(u\right)}^{2n-l}{G}_y\left(u\right)g_y\left(u\right)du\hfill \\ \hfill & =& -2{\lambda }_y\mathbb{H}\left(2,{(2n-l)}^{-{\gamma }_x}\frac{{\gamma }_x{\sigma }_y}{{\gamma }_y{\sigma }_x},-{\gamma }_y,{(2n-l)}^{-{\gamma }_x}\left[1+\frac{{\gamma }_x}{{\sigma }_x}\left({\mu }_y-{\mu }_x-\frac{{\sigma }_y}{{\gamma }_x}\right)\right],-\frac{1}{{\gamma }_x},0\right).\hfill \end{array}$$

(28)

Hence, (23) follows from (25), (26), (27) and (28). On the other hand, when ${\gamma }_x={\gamma }_y=0$, the proof follows the same rationale as in the case of the proof of Theorem 2, just considering the binomial expansion in the process. This proof is omitted for simplicity.    □

4. Estimation

This section deals with parameter estimation for $R=P(X<Y)$ given two independent TGEV random variables. The literature presents maximum likelihood estimators (MLEs) for R considering explicit forms of R obtained after severe parameter restrictions on extreme-value distributions (such as [14,15,22]). Those approaches require the estimation of the parameters to be done jointly in the two samples and require a series of transformations to be properly applied for TGEV components. For the TGEV distribution, we have two cases to consider: $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)\ne 0$ and $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)=0$. The first case requires ${\mu }_y-{\sigma }_y/{\gamma }_y\ge {\mu }_x-{\sigma }_x/{\gamma }_x$ or ${\mu }_y-{\sigma }_y/{\gamma }_y<{\mu }_x-{\sigma }_x/{\gamma }_x$ (Theorem 1). On the other hand, if ${\gamma }_x={\gamma }_y=0$, we release any restrictions on the parameters for the expressions, as a single formula can be used to obtain R in terms of $\mathbb{H}$ functions (Theorem 2).

4.1. MLE for R

Let $X\sim TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and $Y\sim TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$ independent random variables with $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)\ne 0$. Theorem 1 indicates that $R=R\left(\mathit{\theta }\right)$, where we denote $\mathit{\theta }=({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x,{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y).$ Thus, let $\mathbf{x}=(X_1,\cdots ,X_n)$ be a random sample of $TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and consider an independent random sample $\mathbf{y}=(Y_1,\cdots ,Y_m)$ of $TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$, with $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)\ne 0$. Let $\widehat{\mathit{\theta }}=({\widehat{\mu }}_x,{\widehat{\sigma }}_x,{\widehat{\gamma }}_x,{\widehat{\lambda }}_x,{\widehat{\mu }}_y,{\widehat{\sigma }}_y,{\widehat{\gamma }}_y,{\widehat{\lambda }}_y)$ be the estimates of $\mathit{\theta }$. Since Theorem 1 describes R in terms of integrals (hence continuous and measurable functions), we can estimate R simply as $\widehat{R}=R\left(\widehat{\mathit{\theta }}\right)$ due to the invariance property of MLE.

4.2. Parameters Estimation of TGEV Samples

Consider the PDF $f(·;\mu ,\sigma ,\gamma ,\lambda )$ defined in (7). Take $\mathbf{x}=(X_1,\cdots ,X_n)$ and $\mathbf{y}=(Y_1,\cdots ,Y_m)$ independent random samples of sizes n and m, respectively. The likelihood function is given by:

$$L(\mathit{\theta };\mathbf{x},\mathbf{y})=\prod _{j=1}^{n}f(X_j;{\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)\prod _{i=j}^{m}f(Y_j;{\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y).$$

(29)

When ${\gamma }_x={\gamma }_y=0$, the support of f does not depend on unknown parameters and the Theorem 2 does not require parameter restrictions. The log-likelihood function is given by

$$\begin{array}{ccc}\hfill l(\mathit{\theta };\mathbf{x},\mathbf{y})& =& -n{\sigma }_x-m{\sigma }_y-\sum _{j=1}^n\frac{x_j-{\mu }_x}{{\sigma }_x}-\sum _{k=1}^m\frac{y_k-{\mu }_y}{{\sigma }_y}\hfill \\ \hfill & & -\sum _{j=1}^{n}exp\left(-\frac{x_j-{\mu }_x}{{\sigma }_x}\right)+\sum _{j=1}^{n}log\left[(1+{\lambda }_x)-2{\lambda }_{x}exp\left(-exp\left(-\frac{x_j-{\mu }_x}{{\sigma }_x}\right)\right)\right]\hfill \\ \hfill & & -\sum _{k=1}^{m}exp\left(-\frac{y_k-{\mu }_y}{{\sigma }_y}\right)+\sum _{k=1}^{m}log\left[(1+{\lambda }_y)-2{\lambda }_{y}exp\left(-exp\left(-\frac{y_k-{\mu }_y}{{\sigma }_y}\right)\right)\right].\hfill \end{array}$$

(30)

Then, the MLE can be obtained by the log-likelihood function (30), equating its gradient to zero and finding its critical points.

When $sign\left({\gamma }_x\right)=sign\left({\gamma }_y\right)\ne 0$, the support of f depends on the unknown parameter $(\mu ,\sigma ,\gamma )$. Then, we are not able to obtain the MLE explicitly, so an additional numeric procedure is required to perform the likelihood maximization. This is similar to what happens with the GEV distribution (see [6] for a more detailed discussion). The likelihood function becomes:

$$\begin{array}{ccc}\hfill L(\mathit{\theta };\mathbf{x},\mathbf{y})& =& {\sigma }_x^{-n}{\sigma }_y^{-m}exp\left\{-\sum _{j=1}^n{w}_{x,j}^{-1/{\gamma }_x}\right\}exp\left\{-\sum _{k=1}^m{w}_{y,k}^{-1/{\gamma }_y}\right\}\hfill \\ \hfill & \times & \prod _{j=1}^n\frac{\left[(1+{\lambda }_x)-2{\lambda }_{x}exp\left(-w_{x,j}^{-1/{\gamma }_x}\right)\right]}{w_{x,j}^{1+1/{\gamma }_x}}{\mathbb{1}}_{(0,\infty )}\left(w_{x,j}\right)\hfill \\ \hfill & \times & \prod _{k=1}^m\frac{\left[(1+{\lambda }_y)-2{\lambda }_{y}exp\left(-w_{y,k}^{-1/{\gamma }_y}\right)\right]}{w_{y,k}^{1+1/{\gamma }_y}}{\mathbb{1}}_{(0,\infty )}\left(w_{y,k}\right),\hfill \end{array}$$

(31)

where $w_{x,j}=1+{\gamma }_x(x_j-{\mu }_x)/{\sigma }_x$ and $w_{y,k}=1+{\gamma }_y(y_k-{\mu }_y)/{\sigma }_y$. Note that ${\prod }_{j=1}^n{\mathbb{1}}_{(0,\infty )}\left(w_{x,j}\right)>0$ if and only if $w_{x,j}\in (0,\infty )$ for all $j=1,\cdots ,n$. A similar restriction should be observed for $w_{y,k}$. Numerical procedures must be applied to overcome the unavailability of an explicit MLE expression.

4.3. A Two-Step Estimation and Confidence Intervals

We introduce an alternative method for estimating $\mathit{\theta }$ through a two-step process outlined below:

Example 1.

Given the samples $\mathbf{x}$ and $\mathbf{y}$,

• Step 1 We estimate $({\widehat{\mu }}_x,{\widehat{\sigma }}_x,{\widehat{\gamma }}_x)$ and $({\widehat{\mu }}_y,{\widehat{\sigma }}_y,{\widehat{\gamma }}_y)$ using MLE for the GEV (This estimation can be carried out utilizing the extRemes package within the R software version 4.3.3 [23].) distribution;
• Step 2 The parameters $({\lambda }_x,{\lambda }_y)$ are estimated by determining

  $$({\widehat{\lambda }}_x,{\widehat{\lambda }}_y)=arg\underset{{[-1,1]}^2}{max}L({\lambda }_x,{\lambda }_y;\mathbf{x},\mathbf{y}),$$
  where $L({\lambda }_x,{\lambda }_y;\mathbf{x},\mathbf{y})$ is derived from (29) using the estimated parameters from Step 1 as initial guesses.

To choose between Theorems 1 and 2 to obtain $R\left(\widehat{\mathit{\theta }}\right)$, we need to verify if:

(a)

${\widehat{\gamma }}_x\approx {\widehat{\gamma }}_y\approx 0$;

or,

(b)

${\widehat{\mu }}_y-{\widehat{\sigma }}_y/{\widehat{\gamma }}_y\ge {\widehat{\mu }}_x-{\widehat{\sigma }}_x/{\widehat{\gamma }}_x$ or ${\widehat{\mu }}_y-{\widehat{\sigma }}_y/{\widehat{\gamma }}_y<{\widehat{\mu }}_x-{\widehat{\sigma }}_x/{\widehat{\gamma }}_x$.

Despite this additional verification, the computational time required for theExample 1 is expected to be less than that required for directly maximizing (29) and (30).

Example 2 describes the approach used in Section 5 to obtain confidence intervals (CIs) for the estimates of R.

Example 2.

Let $(\mathbf{x},\mathbf{y})$ be a sample of size n and M be a positive integer denoting the number of bootstrap repetitions.

• Step 1 Generate bootstrap samples ${(\mathbf{x},\mathbf{y})}_i$.
• Step 2 Compute the estimates $\widehat{{\mathit{\theta }}_i}={({\widehat{\mu }}_x,{\widehat{\sigma }}_x,{\widehat{\gamma }}_x,{\widehat{\lambda }}_x,{\widehat{\mu }}_y,{\widehat{\sigma }}_y,{\widehat{\gamma }}_y,{\widehat{\lambda }}_y)}_i$ based on ${(\mathbf{x},\mathbf{y})}_i$. In this case, the parameters of each bootstrap sample are individually estimated using Example 1.
• Step 3 Obtain ${\widehat{R}}_i=R\left({\widehat{\mathit{\theta }}}_i\right)$ using Theorem 1 or 2.
• Step 4 Repeat Steps $\mathbf{1}$ to $\mathbf{3}$ M times.
• Step 5 The approximate $100(1-\alpha )\%$ confidence interval of $\widehat{R}$ is given by $[{\widehat{R}}_M(\alpha /2),{\widehat{R}}_M(1-\alpha /2)]$, where ${\widehat{R}}_M\left(\alpha \right)\approx {\widehat{G}}^{-1}\left(\alpha \right)$ and $\widehat{G}$ is the cumulative distribution function of $\widehat{R}$.

For the problem of asset selection using stress–strength reliability, only a single time series of observed returns is available for each asset. Then, the maximum likelihood estimation approach above is of utmost importance. To illustrate the suitability of the analytical closed-form expressions hereby derived, a simulation study is carried out in the next section. In such a case, several samples of size n can be drawn from each random variable, which is then used to estimate the value of R and can be repeated several times.

5. Applications

In this section, we provide a study involving Monte-Carlo simulations that analyze the performance of estimator $\widehat{R}=R\left(\widehat{\mathit{\theta }}\right)$. Additionally, we apply the stress–strength reliability model discussed in the preceding sections to actual real-world data.

5.1. Simulation Study

To evaluate the performance of the estimator $\widehat{R}=R\left(\widehat{\mathit{\theta }}\right)$, we fix several values of the parameters ${\mu }_x$,${\sigma }_x$,${\gamma }_x$, ${\lambda }_x$, ${\mu }_y$,${\sigma }_y$,${\gamma }_y$, ${\lambda }_y$, and then we generate $N\in \{100;1000;10,000\}$ Monte-Carlo samples, each of which of size $n=100$, of the random variables $X\sim TGEV({\mu }_x,{\sigma }_x,{\gamma }_x,{\lambda }_x)$ and $Y\sim TGEV({\mu }_y,{\sigma }_y,{\gamma }_y,{\lambda }_y)$. We analyze the estimates $\widehat{R}$, bias, and root mean squared error (RMSE).

As described by [19], random samples of TGEV distribution can be generated by the inversion method using the quantiles

$$F^{-1}\left(U\right)=\left\{\begin{array}{cc}\frac{\sigma }{\gamma }\left\{-1+{\left[-log\left(\frac{1+\lambda -\sqrt{{(1+\lambda )}^2-4\lambda U}}{2\lambda }\right)\right]}^{-\gamma }\right\},& \gamma \ne 0,\\ \mu +\sigma \left\{-log\left[-log\left(\frac{1+\lambda -\sqrt{{(1+\lambda )}^2-4\lambda U}}{2\lambda }\right)\right]\right\}& \gamma =0,\end{array}\right.$$

where U is a uniform random variable in $[0,1]$.

For the simulation, for each line in the

Table 2. Negative-shape mean. bias and root mean squared error (RMSE) of ${\widehat{R}}_{MC}$$(N=100\mathrm{and}n=100)$.

${\mathit{\mu }}_{\mathit{x}}$	${\mathit{\sigma }}_{\mathit{x}}$	${\mathit{\gamma }}_{\mathit{x}}$	${\mathit{\lambda }}_{\mathit{x}}$	${\mathit{\mu }}_{\mathit{y}}$	${\mathit{\sigma }}_{\mathit{y}}$	${\mathit{\gamma }}_{\mathit{y}}$	${\mathit{\lambda }}_{\mathit{y}}$	R	${\widehat{\mathit{R}}}_{\mathit{MC}}$	Bias	RMSE
0	0.7	−0.1	0.3	0	0.5	−0.2	0.1	0.50945	0.50151	−0.00793	0.00030
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.1	0.71367	0.70374	−0.00993	0.00034
0	0.7	−0.1	0.3	0	0.7	−0.2	0.1	0.52511	0.51976	−0.00534	0.00031
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.1	0.70569	0.70066	−0.00503	0.00026
0	0.7	−0.1	0.3	0	0.5	−0.4	0.1	0.49137	0.48478	−0.00659	0.00029
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.1	0.69977	0.69183	−0.00794	0.00030
0	0.7	−0.1	0.3	0	0.7	−0.4	0.1	0.50809	0.50191	−0.00618	0.00030
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.1	0.69163	0.68617	−0.00546	0.00033
0	0.7	−0.1	0.3	0	0.5	−0.2	0.2	0.49602	0.48709	−0.00893	0.00036
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.2	0.70372	0.69562	−0.00810	0.00034
0	0.7	−0.1	0.3	0	0.7	−0.2	0.2	0.50870	0.50454	−0.00416	0.00029
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.2	0.69299	0.68870	−0.00429	0.00031
0	0.7	−0.1	0.3	0	0.5	−0.4	0.2	0.47863	0.47081	−0.00782	0.00030
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.2	0.68979	0.68537	−0.00442	0.00027
0	0.7	−0.1	0.3	0	0.7	−0.4	0.2	0.49230	0.49000	−0.00230	0.00033
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.2	0.67876	0.67348	−0.00529	0.00030

,

Table 3. Positive-shape mean. bias and root mean squared error (RMSE) of ${\widehat{R}}_{MC}$$(N=100\mathrm{and}n=100)$.

${\mathit{\mu }}_{\mathit{x}}$	${\mathit{\sigma }}_{\mathit{x}}$	${\mathit{\gamma }}_{\mathit{x}}$	${\mathit{\lambda }}_{\mathit{x}}$	${\mathit{\mu }}_{\mathit{y}}$	${\mathit{\sigma }}_{\mathit{y}}$	${\mathit{\gamma }}_{\mathit{y}}$	${\mathit{\lambda }}_{\mathit{y}}$	R	${\widehat{\mathit{R}}}_{\mathit{MC}}$	Bias	RMSE
0	0.7	0.1	0.3	0	0.5	0.2	0.1	0.53183	0.52584	−0.00600	0.00030
0	0.7	0.1	0.3	0.5	0.5	0.2	0.1	0.71971	0.71210	−0.00761	0.00029
0	0.7	0.1	0.3	0	0.7	0.2	0.1	0.54076	0.53666	−0.00410	0.00030
0	0.7	0.1	0.3	0.5	0.7	0.2	0.1	0.71272	0.70832	−0.00440	0.00030
0	0.7	0.1	0.3	0	0.5	0.4	0.1	0.54736	0.54049	−0.00686	0.00029
0	0.7	0.1	0.3	0.5	0.5	0.4	0.1	0.73060	0.72135	−0.00926	0.00031
0	0.7	0.1	0.3	0	0.7	0.4	0.1	0.55478	0.55053	−0.00425	0.00023
0	0.7	0.1	0.3	0.5	0.7	0.4	0.1	0.72470	0.71852	−0.00618	0.00027
0	0.7	0.1	0.3	0	0.5	0.2	0.2	0.51767	0.51685	−0.00082	0.00028
0	0.7	0.1	0.3	0.5	0.5	0.2	0.2	0.71043	0.70462	−0.00581	0.00030
0	0.7	0.1	0.3	0	0.7	0.2	0.2	0.52393	0.52110	−0.00284	0.00028
0	0.7	0.1	0.3	0.5	0.7	0.2	0.2	0.70105	0.69495	−0.00610	0.00029
0	0.7	0.1	0.3	0	0.5	0.4	0.2	0.53280	0.52937	−0.00343	0.00031
0	0.7	0.1	0.3	0.5	0.5	0.4	0.2	0.72130	0.71108	−0.01021	0.00031
0	0.7	0.1	0.3	0	0.7	0.4	0.2	0.53775	0.53576	−0.00199	0.00027
0	0.7	0.1	0.3	0.5	0.7	0.4	0.2	0.71325	0.71132	−0.00193	0.00024

,

Table 4. Negative-shape mean, bias, and root mean squared error (RMSE) of ${\widehat{R}}_{MC}$$(N=1000\mathrm{and}n=100)$.

${\mathit{\mu }}_{\mathit{x}}$	${\mathit{\sigma }}_{\mathit{x}}$	${\mathit{\gamma }}_{\mathit{x}}$	${\mathit{\lambda }}_{\mathit{x}}$	${\mathit{\mu }}_{\mathit{y}}$	${\mathit{\sigma }}_{\mathit{y}}$	${\mathit{\gamma }}_{\mathit{y}}$	${\mathit{\lambda }}_{\mathit{y}}$	R	${\widehat{\mathit{R}}}_{\mathit{MC}}$	Bias	RMSE
0	0.7	−0.1	0.3	0	0.5	−0.2	0.1	0.50945	0.50135	−0.00809	0.00031
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.1	0.71367	0.70749	−0.00618	0.00027
0	0.7	−0.1	0.3	0	0.7	−0.2	0.1	0.52511	0.51830	−0.00680	0.00032
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.1	0.70569	0.69953	−0.00616	0.00028
0	0.7	−0.1	0.3	0	0.5	−0.4	0.1	0.49137	0.48216	−0.00921	0.00033
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.1	0.69977	0.69122	−0.00855	0.00032
0	0.7	−0.1	0.3	0	0.7	−0.4	0.1	0.50809	0.50197	−0.00612	0.00028
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.1	0.69163	0.68582	−0.00581	0.00031
0	0.7	−0.1	0.3	0	0.5	−0.2	0.2	0.49602	0.48986	−0.00616	0.00030
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.2	0.70372	0.69768	−0.00604	0.00028
0	0.7	−0.1	0.3	0	0.7	−0.2	0.2	0.50870	0.50619	−0.00250	0.00027
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.2	0.69299	0.68945	−0.00354	0.00027
0	0.7	−0.1	0.3	0	0.5	−0.4	0.2	0.47863	0.47359	−0.00504	0.00030
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.2	0.68979	0.68290	−0.00689	0.00031
0	0.7	−0.1	0.3	0	0.7	−0.4	0.2	0.49230	0.48817	−0.00413	0.00031
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.2	0.67876	0.67517	−0.00359	0.00030

,

Table 5. Positive-shape mean, bias, and root mean squared error (RMSE) of ${\widehat{R}}_{MC}$$(N=1000\mathrm{and}n=100)$.

${\mathit{\mu }}_{\mathit{x}}$	${\mathit{\sigma }}_{\mathit{x}}$	${\mathit{\gamma }}_{\mathit{x}}$	${\mathit{\lambda }}_{\mathit{x}}$	${\mathit{\mu }}_{\mathit{y}}$	${\mathit{\sigma }}_{\mathit{y}}$	${\mathit{\gamma }}_{\mathit{y}}$	${\mathit{\lambda }}_{\mathit{y}}$	R	${\widehat{\mathit{R}}}_{\mathit{MC}}$	Bias	RMSE
0	0.7	0.1	0.3	0	0.5	0.2	0.1	0.53183	0.52561	−0.00623	0.00031
0	0.7	0.1	0.3	0.5	0.5	0.2	0.1	0.71971	0.71151	−0.00820	0.00031
0	0.7	0.1	0.3	0	0.7	0.2	0.1	0.54076	0.53576	−0.00500	0.00031
0	0.7	0.1	0.3	0.5	0.7	0.2	0.1	0.71272	0.70567	−0.00705	0.00029
0	0.7	0.1	0.3	0	0.5	0.4	0.1	0.54736	0.54108	−0.00628	0.00031
0	0.7	0.1	0.3	0.5	0.5	0.4	0.1	0.73060	0.72414	−0.00646	0.00027
0	0.7	0.1	0.3	0	0.7	0.4	0.1	0.55478	0.54712	−0.00765	0.00030
0	0.7	0.1	0.3	0.5	0.7	0.4	0.1	0.72470	0.71812	−0.00658	0.00032
0	0.7	0.1	0.3	0	0.5	0.2	0.2	0.51767	0.51462	−0.00305	0.00029
0	0.7	0.1	0.3	0.5	0.5	0.2	0.2	0.71043	0.70328	−0.00715	0.00030
0	0.7	0.1	0.3	0	0.7	0.2	0.2	0.52393	0.52199	−0.00194	0.00032
0	0.7	0.1	0.3	0.5	0.7	0.2	0.2	0.70105	0.69579	−0.00526	0.00031
0	0.7	0.1	0.3	0	0.5	0.4	0.2	0.53280	0.52864	−0.00416	0.00031
0	0.7	0.1	0.3	0.5	0.5	0.4	0.2	0.72130	0.71402	−0.00728	0.00029
0	0.7	0.1	0.3	0	0.7	0.4	0.2	0.53775	0.53627	−0.00148	0.00028
0	0.7	0.1	0.3	0.5	0.7	0.4	0.2	0.71325	0.70860	−0.00464	0.00030

and

Table 6. Negative-shape mean, bias, and root mean squared error (RMSE) of ${\widehat{R}}_{MC}$ (N = 10,000 and n = 100).

${\mathit{\mu }}_{\mathit{x}}$	${\mathit{\sigma }}_{\mathit{x}}$	${\mathit{\gamma }}_{\mathit{x}}$	${\mathit{\lambda }}_{\mathit{x}}$	${\mathit{\mu }}_{\mathit{y}}$	${\mathit{\sigma }}_{\mathit{y}}$	${\mathit{\gamma }}_{\mathit{y}}$	${\mathit{\lambda }}_{\mathit{y}}$	R	${\widehat{\mathit{R}}}_{\mathit{MC}}$	Bias	RMSE
0	0.7	−0.1	0.3	0	0.5	−0.2	0.1	0.5094	0.5021	−0.00737	0.00030
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.1	0.7137	0.7058	−0.00783	0.00030
0	0.7	−0.1	0.3	0	0.7	−0.2	0.1	0.5251	0.5198	−0.00529	0.00030
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.1	0.7057	0.6994	−0.00629	0.00030
0	0.7	−0.1	0.3	0	0.5	−0.4	0.1	0.4914	0.4841	−0.00731	0.00029
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.1	0.6998	0.6921	−0.00764	0.00029
0	0.7	−0.1	0.3	0	0.7	−0.4	0.1	0.5081	0.5022	−0.00594	0.00030
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.1	0.6916	0.6854	−0.00621	0.00029
0	0.7	−0.1	0.3	0	0.5	−0.2	0.2	0.4960	0.4908	−0.00521	0.00029
0	0.7	−0.1	0.3	0.5	0.5	−0.2	0.2	0.7037	0.6986	−0.00512	0.00028
0	0.7	−0.1	0.3	0	0.7	−0.2	0.2	0.5087	0.5059	−0.00282	0.00029
0	0.7	−0.1	0.3	0.5	0.7	−0.2	0.2	0.6930	0.6892	−0.00376	0.00028
0	0.7	−0.1	0.3	0	0.5	−0.4	0.2	0.4786	0.4728	−0.00588	0.00029
0	0.7	−0.1	0.3	0.5	0.5	−0.4	0.2	0.6898	0.6837	−0.00610	0.00029
0	0.7	−0.1	0.3	0	0.7	−0.4	0.2	0.4923	0.4885	−0.00379	0.00029
0	0.7	−0.1	0.3	0.5	0.7	−0.4	0.2	0.6788	0.6750	−0.00381	0.00028

the following procedure was carried out:

(1)

for each Monte-Carlo sample, the estimate $\widehat{R}=R\left(\widehat{\mathit{\theta }}\right)$ is computed;

(2)

${\widehat{R}}_{MC}$ is evaluated by taking the sample mean of the Monte-Carlo samples $\widehat{R}$;

(3)

the bias is computed as the difference between the theoretical R value and ${\widehat{R}}_{MC}$. The same applies to the root mean squared error, which also considers the true value as the analytically obtained one.

The TGEV distribution with negative-shape parameters is treated in Table 2, Table 4 and Table 6 (for $N=100$, 1000 and 10,000, respectively), while Table 3 and Table 5 deal with positive-shape parameters. In both cases, the estimator shows good behavior with minimal bias and low root mean squared error. Furthermore, it is clear that increasing the number of replications N leads to the same conclusions.

5.2. Real Data Set Application

Asset selection is addressed to evaluate the proposed framework. To guide the selection of financial assets when managing a portfolio, we adopt metrics of the type $P(X<Y)$.

We start by modeling stock price log-returns as TGEV distributions, and afterward, we compare log-returns from tickers (companies) of different economic sectors and traded on BOVESPA (São Paulo Stock Exchange): BBAS3.SA (banking: Banco do Brasil S.A.), ITUB4.SA (banking: Itaú Unibanco Holding S.A.), VALE3.SA (mining: Vale S.A.) and VIIA3.SA (retail: Via Varejo S.A). From now on, we will omit the “.SA” suffix present on the tickers under analysis. The time series for each ticker represents the daily closing prices in Brazilian currency (R$, BRL) covering the period from 1 January 2022 to 30 April 2023. The analyzed data comprises a total of 331 daily prices.

Figure 2 presents the stock prices for each ticker, highlighting their distinct value scales and volatility. Subsequently, we aim to compare the returns using the expression $P(X<Y)$.

It is important to point out that these data sets were analyzed previously in the literature [3], and here we show that TGEV distribution fits the log-returns better than GEV, according to information criteria. The daily closing prices were imported directly through the software R by the command:

• ticker = "BBAS3.SA"
• quantmod::getSymbols(ticker, src = "yahoo", auto.assign = FALSE,
• from = ’2022-01-01’, to = ’2023-04-30’, return.class = ’xts’).

We assume that the returns are independent. To meet independence requirements, we measure log-returns correlations using Pearson, Kendal, and Spearman methods comparing pairs of stocks that are not correlated. The pairs selected for comparison– VALE3-BBAS3, VALE3-ITUB4, and VALE3-VIIA3–all exhibited correlation measures equal to or below $0.25$, as detailed in the Appendix of [3]. Figure 3 displays the autocorrelation function of log-returns, indicating an absence of discernible temporal correlations among the returns.

Descriptive statistics for the four data sets are presented in

Table 7. Summary statistics for the stock prices log-returns.

Data Set	Min.	1st Qu	Median	Mean	3rd Qu.	Max.	Std. dv.	Skewness	Kurtosis
BBAS3	−0.1057	−0.0097	0.0019	0.0012	0.0136	0.0736	0.0204	−0.3452	5.7413
ITUB4	−0.0492	−0.0105	0.0004	0.0006	0.0109	0.0794	0.0172	0.3809	4.4864
VALE3	−0.0689	−0.0140	0.0001	−0.0002	0.0128	0.0989	0.0231	0.4092	4.5967
VIIA3	−0.1075	−0.0344	−0.0059	−0.0030	0.0231	0.1504	0.0447	0.6144	3.6044

, showing certain symmetry of log-returns around zero and greater variability for VIIA3. The sample size was $n=330$ and each sample informs the daily closing stock price log-return. Across all data sets, the existence of extreme values is a consistent characteristic, which is in accordance with the nature of financial data.

Quintino et al. [3] showed that the GEV distribution adequately fits the data. Our interest lies in determining if the addition of the $\lambda$ parameter provided by the TGEV distribution will improve the model fit. To accomplish this, we employed the two-step estimation method, described in theExample 1. Considering the different number of parameters of the GEV and TGEV models, to conduct a comparative analysis between these models, we utilized the information criteria Akaike Information Criterion (AIC), Bayesian Information Criterion (BIC), and Efficient Determination Criterion (EDC). Parameter estimates for the stock prices log-returns are presented in

Table 8. Parameter estimates for the stock prices log-returns: BBAS4, ITUB4, VIIA3 and VALE3.

Data Set	$\widehat{\mathit{\lambda }}$	$\widehat{\mathit{\mu }}$	$\widehat{\mathit{\sigma }}$	$\widehat{\mathit{\gamma }}$
BBAS3	0.0103	−0.0063	0.0219	−0.2535
ITUB4	−0.0373	−0.0064	0.0165	−0.1545
VALE3	−0.0088	−0.0095	0.0222	−0.1631
VIIA3	−0.0058	−0.0217	0.0396	−0.1170

, while

Table 9. Information criteria and Kolmogorov–Smirnov (KS) p-values for GEV and TGEV models.

Data Set	Distribution	AIC	BIC	EDC	KS p-Value
BBAS3	TGEV	−1615.84	−1654.23	−1621.90	0.0155
	GEV	−1613.83	−1642.63	−1618.38	0.0133
ITUB4	TGEV	−1749.33	−1787.72	−1755.39	0.3303
	GEV	−1747.28	−1776.07	−1751.82	0.4422
VALE3	TGEV	−1556.20	−1594.59	−1562.26	0.2388
	GEV	−1554.18	−1582.97	−1558.73	0.2313
VIIA3	TGEV	−1143.34	−1181.73	−1149.40	0.6851
	GEV	−1141.40	−1170.19	−1145.94	0.7091

shows that all criteria indicate that there was an improvement in the fit when using the TGEV distribution, compared to GEV.

The adequacy of the fitted TGEV distributions can be assessed through graphical evaluation methods. This includes plotting the theoretical PDF over the histogram (Figure 4), comparing the theoretical CDF against the empirical CDF (ECDF) (Figure 5), and examining the Normal Quantile-Quantile plots of the residuals (Figure 6). Although the Kolmogorov–Smirnov test rejects the TGEV adjustment for BBAS3 data, a visual examination of the histogram and ECDF might not discredit the suitability of the TGEV distribution. Furthermore, the Kolmogorov–Smirnov test tends to be overly sensitive, particularly for medium to large sample sizes, leading to its responsiveness even to minor deviations, which might account for this discrepancy.

Reliability measures, denoted as $R=P(X<Y)$, play a pivotal role in an investor’s decision-making process. To simplify, when X and Y symbolize profit from log-returns and $R<1/2$, the investor tends to favor selecting the financial asset corresponding to X. Conversely, if $R>1/2$, the investor leans toward the opposite choice. However, when $R=1/2$, the decision becomes inconclusive. In this sense,

Table 10. Stress–strength probability estimates and Bootstrap confidence interval (CI) for log-returns following TGEV distribution.

X	Y	$\widehat{\mathit{R}}$	95% CI
VALE3	BBAS3	0.53	(0.40; 0.59)
VALE3	ITUB4	0.52	(0.40; 0.59)
VALE3	VIIA3	0.45	(0.38; 0.52)

presents the estimates of $P(X<Y)$ and the $95\%$ Bootstrap confidence intervals, obtained by $\widehat{R}$ and Example 2.

Utilizing the GEV distribution, reliability estimates ${\widehat{R}}^{GEV}$ of 0.54, 0.54, and 0.43 for the VALE3-BBAS3, VALE3-ITUB4, and VALE3-VIIA3 pairs were obtained in [3], respectively. These values closely resembled those outlined in Table 10 for TGEV distribution. Regarding confidence intervals, there was a reduction in the interval size for the last pair, while the sizes remained consistent for the others.

Point estimates can also be compared with an empirical estimator that does not depend on the estimation of parameters or the choice of a probabilistic model. Let one consider the estimator:

$${\widehat{R}}_{NP}=\frac{1}{n}\sum _{j=1}^n{\mathbb{1}}_{\{x_j\le y_j\}},$$

where ${\mathbb{1}}_A$ denotes the indicator function on the set A and n is the sample size. The estimates obtained are, respectively, 0.55, 0.55, and 0.43 for the pairs presented in Table 10, showing the proximity of parametric and nonparametric estimates.

In Table 10, all the confidence intervals crossed the $0.5$ edge; however, reliability measurement for the pair VALE3xVIIA3 brings some evidence that VALE3 should be the asset to be selected.

6. Conclusions

In this paper, we studied the stress–strength reliability $R=P(X<Y)$ when both X and Y follow independent TGEV distributions. Thus, exact expressions for R have been obtained in terms of the extreme-value $\mathbb{H}$-function with minimal parameter restrictions. With additional restrictions, it was shown that R can be calculated in terms of H-functions.

The present work evaluated the advantages yielded by adding a $\lambda$ parameter to the GEV distribution and modelling data sets with the TGEV distribution. The added parameter brought a more complex analytical derivation of $R=P(X<Y)$ and an expected increase in the computational effort to estimate it. To avoid the computational burden of an added parameter, we proposed a two-step estimation where we first fit a GEV model and then estimate the TGEV parameter $\lambda$. Notwithstanding the complexities of an added parameter, information criteria demonstrated the superiority of TGEV models when compared to GEV ones. This advantage is also perceived when estimating probabilities $R=P(X<Y)$ by obtaining better estimates.

Monte-Carlo simulations attested to the performance of the analytical closed-form expressions hereby derived. By applying our methodology to real-world financial data, we could orient a stock selection procedure by calculating $P(X<Y)$ when both X and Y represent stock returns. In summary, when X and Y represent the return of the stock prices and $R<1/2$, the investor should choose the variable X. If $R>1/2$, the opposite occurs. The case $R=1/2$ is inconclusive.

The framework we explored in this work can be a starting point to study probabilities $R=P(X<Y)$ for recently proposed extreme-value distributions like bimodal Gumbell, bimodal Weibull, bimodal GEV, and extreme-value bivariate models.