Product of Spacing Estimation of Stress–Strength Reliability for Alpha Power Exponential Progressively Type-II Censored Data

Abstract

The present study focuses on estimating the stress–strength parameter when the parent distribution is the alpha power exponential model and the available data are progressively Type-II censored. As a starting point, the usual maximum likelihood approach is applied to obtain point and interval estimates of the model parameters, as well as the stress–strength parameter. Another competing strategy employed in this paper is the maximum product of spacing method, which may be thought of as a rival to the maximum likelihood method. The product of spacing approach is used to obtain point and interval estimates for the various parameters. The asymptotic properties of both methods are used to obtain interval estimates of the model parameter and the stress–strength parameter, and the variance of the stress–strength parameter is approximated using the well-known delta method. Two parametric bootstrap confidence intervals are provided based on the suggested classical estimation procedures. A simulation study is also used to assess the performance of various point and interval estimations. For illustrative purposes, the proposed methods are applied to two real data sets, one for the kidney patients’ recurrence times to infection and the other for breaking the strength of jute fibers.

1. Introduction

Since the work of Birnbaum [1], the issue of estimating the stress–strength parameter $R=P(Y<X)$ and its associated inference when X and Y are two independent random variables has attracted a lot of interest. If X reflects a component’s strength under the stress of Y, then R can be thought of as an indicator of the performance of the system. When the stress on the system surpasses its strength, the system loses control. Since R reflects a relationship between a system’s stress and strength, it is frequently referred to as the system’s stress–strength parameter. The estimation problems of the stress–strength parameter for numerous known distributions have been discussed by several authors in the statistical literature. Based on the complete sample data, Chung [2], Kundu and Gupta [3], Raqab et al. [4] and Sharma [5] addressed the estimation of R for the exponential, Weibull, three-parameter generalized exponential, and generalized inverse Lindley distributions, respectively. The estimation of R based on record and censored data has also been taken into consideration by many authors. For example, in dealing with lower record values, Baklizi [6] estimated R for the generalized exponential distribution. Valiollahi et al. [7] and Asgharzadeh et al. [8] studied the estimation of R for the Weibull distribution, given progressively Type-II and hybrid censored data, respectively. See also the work of Yadav et al. [9], Rostamian and Nematollahi [10], Ghanbari [11], and Asadi and Panahi [12].

In a lot of lifespan tests, the experimenter may remove the cases they are looking at before the experiment ends. The removal procedure could be unintentional or preplanned by the experimenter. To further explain the situation, it should be noted that during clinical trials, participants may decline to complete for personal reasons or the trial may end before all patients whose lifetimes last longer years have been observed. The experimenter would then encounter censored data. Therefore, a censoring scheme has to be specifically defined. Many censoring schemes are available in the literature. One of the most important types of censoring samples is the progressively Type-II censored sample, which operates as follows: Let $r=(r_1,r_2,\dots ,r_n)$ be a provided vector defined by the researcher that shows the line of removal of some experimental elements before ending the test. In this censored procedure, $n(n<N)$ stands for the number of observed failures, while N represents the number of experimental elements. Upon the initial failure, $r_1$ of the surviving units is arbitrarily eliminated from the remaining units in the experiment. When the second failure happens, $r_2$ of the remaining units is randomly eliminated, and the experiment proceeds until the $n^{th}$ failure, then all of the remaining units are eliminated. One can see that the complete sample can be derived as a special case of the progressively Type-II censored sample by setting $(r_1=\cdots =r_n=0)$. Also, the Type-II censored sample can be given by putting $(r_1=\cdots =r_{n-1}=0)$ and $r_n=N-n$. Readers are highly recommended to take a look at the book authored by Balakrishnan and Aggarwala [13] for additional information on the theory, methodology, and applications of progressively Type-II censored data. In addition, one may refer to the book authored by Balakrishnan and Cramer [14] for more forms of censoring schemes, inference of parameters, and its applications in reliability.

Adding a parameter to well-known distributions is common in the literature to increase the flexibility of the existing ones. Mahdavi and Kundu [15] recently proposed a new variant of the conventional exponential distribution by adding a new shape parameter, naming the new distribution the alpha power exponential (APE) distribution. In terms of real data analysis, they showed that the APE distribution has more flexibility than some well-known distributions like Weibull and gamma distributions in modeling various types of data. The APE distribution with shape parameter $\alpha$ and scale parameter $\sigma$, denoted by $APE(\alpha ,\theta )$ has the following probability density function (PDF):

$$f(x;\alpha ,\sigma )=\frac{\sigma log\left(\alpha \right)e^{-\sigma x}{\alpha }^{1-e^{-\sigma x}}}{\alpha -1},x>0,\sigma ,\alpha >0,\alpha \ne 1,$$

(1)

and the corresponding cumulative distribution function (CDF) is given by

$$F(x;\alpha ,\sigma )=\frac{{\alpha }^{1-e^{-\sigma x}}-1}{\alpha -1}.$$

(2)

Despite the importance of the APE distribution, little work has been performed on the estimation of the distribution parameters and some of its indices. For example, Nassar et al. [16] studied some classical parameter estimation methods. Alotaibi et al. [17] considered some statistical inferences using adaptive progressive censoring. Nassar et al. [18] studied the competing risks model using progressively Type-II censored samples. However, no one has paid attention to the estimation of the stress–strength parameter under censoring in the case of the APE distribution. Furthermore, the majority of studies that estimated the stress–strength parameter used the maximum likelihood (ML) method as the sole approach among the classical ones, whereas the maximum product of spacing (MPS) method can be considered an effective alternative, particularly for small sample sizes as pointed out by Anatolyev and Kosenok [19]. As a result, the main focus of this study is to analyze the estimation issues of the stress–strength parameter of the APE distribution, with a common scale parameter, in the case of progressively Type-II censored data. To accomplish this, two classical estimation methods are considered, namely ML and MPS. Besides the point estimates, the approximate confidence intervals (ACIs) are investigated based on the asymptotic properties of the ML estimates (MLEs) and MPS estimates (MPSEs). Moreover, two parametric bootstrap confidence intervals are obtained based on both the MLEs and MPSEs.

The remainder of the paper is divided into the following sections: Section 2 presents the MLEs and the related ACIs of the unknown parameters and the stress–strength parameter. Section 3 discusses the MPSEs and ACIs constructed using the MPSEs of the various parameters. Section 4 displays the bootstrap confidence intervals based on MLEs and MPSEs. Section 5 describes the simulation research that was conducted to compare the performance of the proposed estimators. Section 6 demonstrates the study’s applicability by examining two real-world data sets. Section 7 summarizes the findings of the study.

2. Maximum Likelihood Estimation

In this section, both point and interval estimations of R are investigated by employing the ML method. Let X and Y be independent APE random variables, where $X\sim APE({\alpha }_1,\sigma )$ and $Y\sim APE({\alpha }_2,\sigma )$, respectively, then according to Mahdavi and Kundu [15], we can express the parameter $R=P(Y<X)$ as follows:

$$R=\frac{({\alpha }_1{\alpha }_2-1)log\left({\alpha }_2\right)}{({\alpha }_1-1)({\alpha }_2-1)log\left({\alpha }_1{\alpha }_2\right)}-\frac{1}{{\alpha }_1-1}.$$

(3)

One can see from (3) that the parameter R is a function of ${\alpha }_1$ and ${\alpha }_2$. Then, the point estimator of R can be obtained based on the invariance property of the MLEs after obtaining the MLEs of ${\alpha }_1$ and ${\alpha }_2$. Let $(X_{1:n:N},\dots ,X_{n:n:N})$ be a progressively Type-II censored sample randomly taken from a population that follows $APE({\alpha }_1,\sigma )$ distribution with the censoring scheme $(r_1,\dots ,r_n)$, and $(Y_{1:m:M},\dots ,Y_{m:m:M})$ be a progressively Type-II censored sample randomly selected from the $APE({\alpha }_2,\sigma )$ population with the censoring scheme $(s_1,\dots ,s_m)$. Then, we can write the likelihood function of the observed data as

$$\begin{array}{ccc}\hfill L({\alpha }_1,{\alpha }_2,\sigma )& =& A_1{A}_2\left\{\prod _{i=1}^{n}f\left(x_i\right){[1-F\left(x_i\right)]}^{r_i}\right\}\left\{\prod _{j=1}^{m}f\left(y_j\right){[1-F\left(y_j\right)]}^{s_j}\right\},\hfill \end{array}$$

(4)

where $A_1$ and $A_2$ are two constants that do not depend on the unknown parameters, $x_i=x_{i:n:N}$ and $y_i=y_{j:m:M}$ for simplicity. From (1), (2), and (4), we can write the likelihood function of ${\alpha }_1,{\alpha }_2$, and $\sigma$, ignoring the constant term, in the following form:

$$\begin{array}{ccc}\hfill L({\alpha }_1,{\alpha }_2,\sigma )& =& {\left(\frac{{\alpha }_1}{{\alpha }_1-1}\right)}^N{\left(\frac{{\alpha }_2}{{\alpha }_2-1}\right)}^M{\sigma }^{n+m}{log}^n\left({\alpha }_1\right){log}^m\left({\alpha }_2\right)\hfill \\ & \times & exp\left[-\sigma \left(\sum _{i=1}^n{x}_i+\sum _{j=1}^m{y}_j\right)-log\left({\alpha }_1\right)\sum _{i=1}^n{e}^{-\sigma x_i}-log\left({\alpha }_2\right)\sum _{j=1}^m{e}^{-\sigma y_j}\right]\hfill \\ & \times & \prod _{i=1}^n{\left(1-{\alpha }_1^{-e^{-\sigma x_i}}\right)}^{r_i}\prod _{j=1}^m{\left(1-{\alpha }_2^{-e^{-\sigma y_j}}\right)}^{s_j}.\hfill \end{array}$$

(5)

The natural logarithm of (5) can be expressed as shown below:

$$\begin{array}{ccc}\hfill l({\alpha }_1,{\alpha }_2,\sigma )& =& Nlog\left({\alpha }_1\right)+Mlog\left({\alpha }_2\right)-Nlog({\alpha }_1-1)-Mlog({\alpha }_2-1)+(n+m)log\left(\sigma \right)\hfill \\ & +& nlog[log\left({\alpha }_1\right)]+mlog[log\left({\alpha }_2\right)]-\sigma \left(\sum _{i=1}^n{x}_i+\sum _{j=1}^m{y}_j\right)-log\left({\alpha }_1\right)\sum _{i=1}^n{e}^{-\sigma x_i}\hfill \\ & -& log\left({\alpha }_2\right)\sum _{j=1}^m{e}^{-\sigma y_j}+\sum _{i=1}^n{r}_{i}log\left(1-{\alpha }_1^{-e^{-\sigma x_i}}\right)+\sum _{j=1}^m{s}_{j}log\left(1-{\alpha }_2^{-e^{-\sigma y_j}}\right).\hfill \end{array}$$

(6)

The MLEs of ${\alpha }_1,{\alpha }_2$, and $\sigma$, denoted by ${\widehat{\alpha }}_1,{\widehat{\alpha }}_2$, and $\widehat{\sigma }$, can be acquired by maximizing (6) with respect to the unknown parameters, alternatively by solving the following normal equations:

$$\begin{array}{ccc}\hfill \frac{\partial l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1}& =& \frac{1}{{\alpha }_1}\left[N+\frac{n}{log\left({\alpha }_1\right)}\right]-\frac{N}{{\alpha }_1-1}-\frac{1}{{\alpha }_1}\sum _{i=1}^n{e}^{-\sigma x_i}+\frac{1}{{\alpha }_1}\sum _{i=1}^n\frac{r_i{e}^{-\sigma x_i}}{{\alpha }_1^{e^{-\sigma x_i}}-1}=0,\hfill \end{array}$$

(7)

$$\begin{array}{ccc}\hfill \frac{\partial l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2}& =& \frac{1}{{\alpha }_2}\left[M+\frac{m}{log\left({\alpha }_2\right)}\right]-\frac{M}{{\alpha }_2-1}-\frac{1}{{\alpha }_2}\sum _{j=1}^m{e}^{-\sigma y_j}+\frac{1}{{\alpha }_2}\sum _{j=1}^m\frac{s_j{e}^{-\sigma y_j}}{{\alpha }_2^{e^{-\sigma y_j}}-1}=0\hfill \end{array}$$

(8)

and

$$\begin{array}{ccc}\hfill \frac{\partial l({\alpha }_1,{\alpha }_2,\sigma )}{\partial \sigma }& =& \frac{n+m}{\sigma }-\left(\sum _{i=1}^n{x}_i+\sum _{j=1}^m{y}_j\right)+log\left({\alpha }_1\right)\sum _{i=1}^n{x}_i{e}^{-\sigma x_i}\left(1-\frac{r_i}{{\alpha }_1^{e^{-\sigma x_i}}-1}\right)\hfill \\ & +& log\left({\alpha }_2\right)\sum _{j=1}^m{y}_j{e}^{-\sigma y_j}\left(1-\frac{s_j}{{\alpha }_2^{e^{-\sigma y_j}}-1}\right)=0.\hfill \end{array}$$

(9)

It is clear that the normal equations in (7)–(9) require a numerical solution in order to yield the MLEs ${\widehat{\alpha }}_1,{\widehat{\alpha }}_2$, and $\widehat{\sigma }$. This is because the normal equations do not yield explicit solutions for ${\alpha }_1,{\alpha }_2$, and $\sigma$. After obtaining the MLEs, one can easily obtain the MLE of R from (3) by substituting the unknown parameters with their respective MLEs, as shown below, owing to the invariance property of the MLEs:

$$\widehat{R}=\frac{({\widehat{\alpha }}_1{\widehat{\alpha }}_2-1)log\left({\widehat{\alpha }}_2\right)}{({\widehat{\alpha }}_1-1)({\widehat{\alpha }}_2-1)log\left({\widehat{\alpha }}_1{\widehat{\alpha }}_2\right)}-\frac{1}{{\widehat{\alpha }}_1-1}.$$

The normal approximation of the MLEs is used in this part to derive the ACI of the unknown parameters. As a result of the asymptotic properties of MLEs, we can investigate the interval estimation for R. It is known, under some mild regularity conditions, that the asymptotic distribution of the MLEs $({\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma })$ is a trivariate normal distribution with mean $({\alpha }_1,{\alpha }_2,\sigma )$ and a variance–covariance matrix $\mathit{J}({\alpha }_1,{\alpha }_2,\sigma )$. In practice, we can estimate $\mathit{J}({\alpha }_1,{\alpha }_2,\sigma )$ by $\mathit{J}({\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma })$ due to the difficulties of obtaining the expected values of the second-order derivatives of the log-likelihood function. In this case, the estimated variance–covariance matrix $\mathit{J}({\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma })$ is obtained based on the empirical Fisher information matrix as

$$\begin{array}{ccc}\hfill \mathit{J}({\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma })& =& {\left(\begin{array}{ccc}-\frac{{\partial }^{2}l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1^2}& 0& -\frac{{\partial }^{2}l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1\sigma }\\ & -\frac{{\partial }^{2}l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2^2}& -\frac{{\partial }^{2}l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2\partial \sigma }\\ & & -\frac{{\partial }^{2}l({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\sigma }^2}\end{array}\right)}_{({\alpha }_1,{\alpha }_2,\sigma )=({\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma })}^{-1}\hfill \\ & =& \left(\begin{array}{ccc}{\widehat{V}}_{11}& 0& {\widehat{V}}_{13}\\ & {\widehat{V}}_{22}& {\widehat{V}}_{23}\\ & & {\widehat{V}}_{33}\end{array}\right),\hfill \end{array}$$

(10)

where

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1^2}& =& \frac{N(2{\alpha }_1-1)}{{\left[{\alpha }_1({\alpha }_1-1)\right]}^2}-\frac{n[1+log\left({\alpha }_1\right)]}{{[{\alpha }_{1}log\left({\alpha }_1\right)]}^2}+\frac{1}{{\alpha }_1^2}\sum _{i=1}^n{e}^{-\sigma x_i}+\sum _{i=1}^n\frac{r_i{e}^{-\sigma x_i}{\varphi }_{1i}}{{\left({\alpha }_1-{\alpha }_1^{1+e^{-\sigma x_i}}\right)}^2},\hfill \end{array}$$

(11)

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2^2}& =& \frac{M(2{\alpha }_2-1)}{{\left[{\alpha }_2({\alpha }_2-1)\right]}^2}-\frac{m[1+log\left({\alpha }_2\right)]}{{[{\alpha }_{2}log\left({\alpha }_2\right)]}^2}+\frac{1}{{\alpha }_2^2}\sum _{j=1}^m{e}^{-\sigma y_j}+\sum _{j=1}^m\frac{s_j{e}^{-\sigma y_j}{\varphi }_{2j}}{{\left({\alpha }_2-{\alpha }_2^{1+e^{-\sigma y_j}}\right)}^2},\hfill \end{array}$$

(12)

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\sigma }^2}& =& -\frac{n+m}{{\sigma }^2}-log\left({\alpha }_1\right)\sum _{i=1}^n{x}_i^2{e}^{-\sigma x_i}(1+r_i{\phi }_{1i})-log\left({\alpha }_2\right)\sum _{j=1}^m{y}_j^2{e}^{-\sigma y_j}(1+s_j{\phi }_{2j}),\hfill \end{array}$$

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1\partial \sigma }& =& \frac{1}{{\alpha }_1}\sum _{i=1}^n{x}_i{e}^{-\sigma x_i}(1+r_i{\phi }_{1i})\hfill \end{array}$$

and

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2\partial \sigma }& =& \frac{1}{{\alpha }_2}\sum _{j=1}^m{y}_j{e}^{-\sigma y_j}(1+s_j{\phi }_{2j}),\hfill \end{array}$$

where

$${\varphi }_{kt}=1-(1+e^{-\sigma z_t}){\alpha }_k^{e^{-\sigma z_t}}\mathrm{and}{\phi }_{kt}=\frac{1+{\alpha }_k^{e^{-\sigma z_t}}[log\left({\alpha }_k\right)e^{-\sigma z_t}-1]}{{\left(1-{\alpha }_k^{e^{-\sigma z_t}}\right)}^2},$$

with $k=1,2$, $(z,t)=(x,i)$ for $k=1$ and $(z,t)=(y,j)$ for $k=2$, $i=1,\dots ,n,j=1,\dots ,m$. Therefore, the $100(1-\tau )$ ACIs for ${\alpha }_1,{\alpha }_2$, and $\sigma$ can be formulated, respectively, as

$${\widehat{\alpha }}_1\pm z_{\tau /2}\sqrt{{\widehat{V}}_{11}},{\widehat{\alpha }}_2\pm z_{\tau /2}\sqrt{{\widehat{V}}_{22}}\mathrm{and}\widehat{\theta }\pm z_{\tau /2}\sqrt{{\widehat{V}}_{33}},$$

(13)

where $z_{\tau /2}$ is the upper ${(\tau /2)}^{th}$ percentile point of a standard normal distribution. Let $n=m$; using the delta method, the distribution of $\widehat{R}$ as n tends to infinity is a normal distribution with mean R and variance $V_R$, where $V_R=\acute{\mathbf{\Lambda }}\mathit{J}({\alpha }_1,{\alpha }_2,\sigma )\mathbf{\Lambda }$, where $\acute{\mathbf{\Lambda }}=\left(R_{{\alpha }_1},R_{{\alpha }_2},0\right)$ is the vector of the first-order partial derivativesof R with respect to the unknown parameters; see, for more detail, the work of Davison [20], and Nadarajah and Jia [21]. In practice, the variance $V_R$ can be estimated by replacing the unknown parameters with their MLEs as follows:

$$\begin{array}{ccc}\hfill {\widehat{V}}_R& =& \left({\widehat{R}}_{{\alpha }_1},{\widehat{R}}_{{\alpha }_2},0\right)\left(\begin{array}{ccc}{\widehat{V}}_{11}& 0& {\widehat{V}}_{13}\\ & {\widehat{V}}_{22}& {\widehat{V}}_{23}\\ & & {\widehat{V}}_{33}\end{array}\right)\left(\begin{array}{c}{\widehat{R}}_{{\alpha }_1}\\ {\widehat{R}}_{{\alpha }_2}\\ 0\end{array}\right),\hfill \end{array}$$

(14)

where

$${\widehat{R}}_{{\alpha }_1}=\frac{1}{{({\widehat{\alpha }}_1-1)}^2}+\omega \left\{{\widehat{\alpha }}_2-({\widehat{\alpha }}_1{\widehat{\alpha }}_2-1)\left[\frac{1}{{\widehat{\alpha }}_1-1}+\frac{1}{{\widehat{\alpha }}_{1}log\left({\widehat{\alpha }}_1{\widehat{\alpha }}_2\right)}\right]\right\}$$

(15)

and

$${\widehat{R}}_{{\alpha }_2}=\frac{\omega ({\widehat{\alpha }}_1{\widehat{\alpha }}_2-1)}{{\widehat{\alpha }}_2}\left[\frac{1}{log\left({\widehat{\alpha }}_2\right)}+\frac{{\widehat{\alpha }}_1{\widehat{\alpha }}_2}{{\widehat{\alpha }}_1{\widehat{\alpha }}_2-1}-\frac{{\widehat{\alpha }}_2}{{\widehat{\alpha }}_2-1}-\frac{1}{log\left({\widehat{\alpha }}_1{\widehat{\alpha }}_2\right)}\right],$$

(16)

where

$$\omega =\frac{log\left({\widehat{\alpha }}_2\right)}{({\widehat{\alpha }}_1-1)({\widehat{\alpha }}_2-1)log\left({\widehat{\alpha }}_1{\widehat{\alpha }}_2\right)}.$$

From (14), the approximate estimated variance of $\widehat{R}$ can be as follows:

$${\widehat{V}}_R={\widehat{R}}_{{\alpha }_1}^2{\widehat{V}}_{11}+{\widehat{R}}_{{\alpha }_2}^2{\widehat{V}}_{22}.$$

This outcome can be used to build the ACI of R. With a $100\%(1-\tau )$ confidence level, the ACI of R can be obtained as

$$\widehat{R}\pm z_{\tau /2}\sqrt{{\widehat{V}}_R}.$$

3. Maximum Product of Spacing Estimation

Cheng and Amin [22] proposed the MPS approach as an alternative to the ML method. The MPS method has all of the large sample properties of the ML method and keeps most of the ML method’s properties under more general conditions; for more information, see Cheng and Iles [23] and Cheng and Traylor [24]. The MPSEs are calculated by selecting parameter values that maximize the product of the distances between the distribution function values at neighboring ordered points. Numerous authors recently used the MPS approach to estimate various lifetime distributions, such as Basu et al. [25,26] and Nassar et al. [27]. In this part, the MPS procedure is used to derive point and interval estimators of R as an alternative to the ML approach. As in the case of the ML method, the point and interval estimations of the parameter R are acquired based on some properties of MPSEs, including the invariance property and normal approximation. Let $(X_{1:n:N},\dots ,X_{n:n:N})$ be a progressively Type-II censored sample drawn at random from a population with an $APE({\alpha }_1,\sigma )$ distribution with the censoring scheme $(r_1,\dots ,r_n)$, and let $(Y_{1:m:M},\dots ,Y_{m:m:M})$ be a progressively Type-II censored sample randomly chosen from the $APE({\alpha }_2,\sigma )$ population with the censoring scheme $(s_1,\dots ,s_m)$. To apply the MPS method, we first need to define the following quantities based on the CDF given by (2):

$$\begin{array}{ccc}\hfill {▵}_{1i}& =& F\left(x_{i:n:N}\right|{\alpha }_1,\sigma )-F\left(x_{i-1:n:N}\right|{\alpha }_1,\sigma )\hfill \\ & =& \frac{{\alpha }_1}{{\alpha }_1-1}\left({\alpha }_1^{-e^{-\sigma x_i}}-{\alpha }_1^{-e^{-\sigma x_{i-1}}}\right)\hfill \end{array}$$

(17)

and

$$\begin{array}{ccc}\hfill {▵}_{2j}& =& F\left(y_{j:m:M}\right|{\alpha }_2,\sigma )-F\left(y_{j-1:m:M}\right|{\alpha }_2,\sigma )\hfill \\ & =& \frac{{\alpha }_2}{{\alpha }_2-1}\left({\alpha }_2^{-e^{-\sigma y_j}}-{\alpha }_2^{-e^{-\sigma y_{j-1}}}\right),\hfill \end{array}$$

(18)

where $x_i=x_{i:n:N}$ and $y_i=y_{j:m:M}$. Based on (17) and (18), we can write the product of the spacing function to be maximized as shown below:

$$\begin{array}{ccc}\hfill P({\alpha }_1,{\alpha }_2,\sigma )& =& {\left(\frac{{\alpha }_1}{{\alpha }_1-1}\right)}^{N+1}{\left(\frac{{\alpha }_2}{{\alpha }_2-1}\right)}^{M+1}\prod _{i=1}^{n+1}{D}_{1i}\prod _{j=1}^{m+1}{D}_{2j}\prod _{i=1}^n{\left(1-{\alpha }_1^{-e^{-\sigma x_i}}\right)}^{r_i}\hfill \\ & \times & \prod _{j=1}^m{\left(1-{\alpha }_2^{-e^{-\sigma y_j}}\right)}^{s_j}.\hfill \end{array}$$

(19)

where $D_{1i}=({\alpha }_1-1){▵}_{1i}/{\alpha }_1$ and $D_{2j}=({\alpha }_2-1){▵}_{2j}/{\alpha }_2$. From the product of the spacing function in (19), we can write the natural logarithm of $P({\alpha }_1,{\alpha }_2,\sigma )$ as follows:

$$\begin{array}{ccc}\hfill p({\alpha }_1,{\alpha }_2,\sigma )& =& N^{*}log\left({\alpha }_1\right)+M^{*}log\left({\alpha }_2\right)-N^{*}log({\alpha }_1-1)-M^{*}log({\alpha }_2-1)+\sum _{i=1}^{n+1}log\left(D_{1i}\right)\hfill \\ & +& \sum _{j=1}^{m+1}log\left(D_{2j}\right)+\sum _{i=1}^n{r}_{i}log\left(1-{\alpha }_1^{-e^{-\sigma x_i}}\right)+\sum _{j=1}^m{s}_{j}log\left(1-{\alpha }_2^{-e^{-\sigma y_j}}\right),\hfill \end{array}$$

(20)

where $N^{*}=N+1$ and $M^{*}=M+1$. The MPSEs of ${\alpha }_1,{\alpha }_2$, and $\sigma$, denoted by ${\tilde{\alpha }}_1,{\tilde{\alpha }}_2$, and $\tilde{\sigma }$, can be acquired by solving the following normal equations:

$$\begin{array}{ccc}\hfill \frac{\partial p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1}& =& \frac{N^{*}}{{\alpha }_1}-\frac{N^{*}}{{\alpha }_1-1}-\sum _{i=1}^{n+1}\frac{{\Psi }_{1i}}{D_{1i}}+\frac{1}{{\alpha }_1}\sum _{i=1}^n\frac{r_i{e}^{-\sigma x_i}}{{\alpha }_1^{e^{-\sigma x_i}}-1}=0,\hfill \end{array}$$

(21)

$$\begin{array}{ccc}\hfill \frac{\partial p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2}& =& \frac{M^{*}}{{\alpha }_2}-\frac{M^{*}}{{\alpha }_2-1}-\sum _{j=1}^{m+1}\frac{{\Psi }_{2j}}{D_{2j}}+\frac{1}{{\alpha }_2}\sum _{j=1}^m\frac{s_j{e}^{-\sigma y_j}}{{\alpha }_2^{e^{-\sigma y_j}}-1}=0\hfill \end{array}$$

(22)

and

$$\begin{array}{ccc}\hfill \frac{\partial p({\alpha }_1,{\alpha }_2,\sigma )}{\partial \sigma }& =& \sum _{i=1}^{n+1}\frac{{\xi }_{1i}}{D_{1i}}+\sum _{j=1}^{m+1}\frac{{\xi }_{2j}}{D_{2j}}+log\left({\alpha }_1\right)\sum _{i=1}^n{x}_i{e}^{-\sigma x_i}\left(1-\frac{r_i}{{\alpha }_1^{e^{-\sigma x_i}}-1}\right)\hfill \\ & +& log\left({\alpha }_2\right)\sum _{j=1}^m{y}_j{e}^{-\sigma y_j}\left(1-\frac{s_j}{{\alpha }_2^{e^{-\sigma y_j}}-1}\right)=0.\hfill \end{array}$$

(23)

where ${\Psi }_{kt}=\frac{e^{-\sigma z_t}}{{\alpha }_k^{1+e^{-\sigma z_t}}}-\frac{e^{-\sigma z_{t-1}}}{{\alpha }_k^{1+e^{-\sigma z_{t-1}}}}$ and ${\xi }_{kt}=\frac{z_t{e}^{-\sigma z_t}log\left({\alpha }_k\right)}{{\alpha }_k^{e^{-\sigma z_t}}}-\frac{z_{t-1}{e}^{-\sigma z_{t-1}}log\left({\alpha }_k\right)}{{\alpha }_k^{e^{-\sigma z_{t-1}}}}$, and the values of $k,z$, and t as mentioned before in the previous section. As (21)–(23) show, the MPSEs of ${\alpha }_1,{\alpha }_2$, and $\sigma$ do not have closed forms; therefore, they can be solved numerically as in the case of MLEs to obtain ${\tilde{\alpha }}_1,{\tilde{\alpha }}_2$, and $\tilde{\sigma }$. According to Ranneby [28] and Coolen and Newby [29], the MPSEs possess a similar invariance property as MLEs. Utilizing this property, the MPSE of R can be acquired as follows:

$$\tilde{R}=\frac{({\tilde{\alpha }}_1{\tilde{\alpha }}_2-1)log\left({\tilde{\alpha }}_2\right)}{({\tilde{\alpha }}_1-1)({\tilde{\alpha }}_2-1)log\left({\tilde{\alpha }}_1{\tilde{\alpha }}_2\right)}-\frac{1}{{\tilde{\alpha }}_1-1}.$$

One can simply derive the ACIs of the model parameters as well as the parameter R based on the asymptotic normality of the MPSEs. Studies on the asymptotic normality of the MSPEs can be found in the works of Ranneby [28], Cheng and Stephens [30] and Ghosh and Jammalamadaka [31]. According to this property, the MPSEs $({\tilde{\alpha }}_1,{\tilde{\alpha }}_2,\tilde{\sigma })$ are asymptotically normally distributed, i.e., $({\tilde{\alpha }}_1,{\tilde{\alpha }}_2,\tilde{\sigma })\sim N[({\alpha }_1,{\alpha }_2,\sigma ),\mathit{I}({\alpha }_1,{\alpha }_2,\sigma )]$. In practice, we use $\mathit{I}({\tilde{\alpha }}_1,{\tilde{\alpha }}_2,\tilde{\sigma })$ to estimate the variance–covariance matrix $\mathit{I}({\alpha }_1,{\alpha }_2,\sigma )$, where

$$\begin{array}{ccc}\hfill \mathit{I}({\tilde{\alpha }}_1,{\tilde{\alpha }}_2,\tilde{\sigma })& =& {\left(\begin{array}{ccc}-\frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1^2}& 0& -\frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1\sigma }\\ & -\frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2^2}& -\frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2\partial \sigma }\\ & & -\frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\sigma }^2}\end{array}\right)}_{({\alpha }_1,{\alpha }_2,\sigma )=({\tilde{\alpha }}_1,{\tilde{\alpha }}_2,\tilde{\sigma })}^{-1}\hfill \\ & =& \left(\begin{array}{ccc}{\tilde{V}}_{11}& 0& {\tilde{V}}_{13}\\ & {\tilde{V}}_{22}& {\tilde{V}}_{23}\\ & & {\tilde{V}}_{33}\end{array}\right),\hfill \end{array}$$

(24)

where

$$\begin{array}{ccc}\hfill \frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1^2}& =& \frac{N(2{\alpha }_1-1)}{{\left[{\alpha }_1({\alpha }_1-1)\right]}^2}-\sum _{i=1}^{n+1}\frac{D_{1i}{\acute{\Psi }}_{1i}-{\Psi }_{1i}^2}{D_{1i}^2}+\sum _{i=1}^n\frac{r_i{e}^{-\sigma x_i}{\varphi }_{1i}}{{\left({\alpha }_1-{\alpha }_1^{1+e^{-\sigma x_i}}\right)}^2},\hfill \end{array}$$

$$\begin{array}{ccc}\hfill \frac{{\partial }^{2}p({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2^2}& =& \frac{M(2{\alpha }_2-1)}{{\left[{\alpha }_2({\alpha }_2-1)\right]}^2}-\sum _{j=1}^{m+1}\frac{D_{2j}{\acute{\Psi }}_{2j}-{\Psi }_{2j}^2}{D_{2j}^2}+\sum _{j=1}^m\frac{s_j{e}^{-\sigma y_j}{\varphi }_{2j}}{{\left({\alpha }_2-{\alpha }_2^{1+e^{-\sigma y_j}}\right)}^2},\hfill \end{array}$$

$$\begin{array}{ccc}\hfill \frac{\partial l^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\sigma }^2}& =& \sum _{i=1}^{n+1}\frac{D_{1i}{\acute{\xi }}_{1i}-{\xi }_{1i}^2}{D_{1i}^2}+\sum _{j=1}^{m+1}\frac{D_{2j}{\acute{\xi }}_{2j}-{\xi }_{2j}^2}{D_{2j}^2}-log\left({\alpha }_1\right)\sum _{i=1}^n{x}_i^2{e}^{-\sigma x_i}(1+r_i{\phi }_{1i})\hfill \\ & -& log\left({\alpha }_2\right)\sum _{j=1}^m{y}_j^2{e}^{-\sigma y_j}(1+s_j{\phi }_{2j}),\hfill \end{array}$$

$$\begin{array}{ccc}\hfill \frac{\partial p^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_1\partial \sigma }& =& \sum _{i=1}^{n+1}\frac{D_{1i}{\nu }_{1i}-{\Psi }_{1i}{\xi }_{1i}}{D_{1i}^2}+\frac{1}{{\alpha }_1}\sum _{i=1}^n{x}_i{e}^{-\sigma x_i}(1+r_i{\phi }_{1i})\hfill \end{array}$$

and

$$\begin{array}{ccc}\hfill \frac{\partial p^2({\alpha }_1,{\alpha }_2,\sigma )}{\partial {\alpha }_2\partial \sigma }& =& \sum _{j=1}^{m+1}\frac{D_{2j}{\nu }_{2j}-{\Psi }_{2j}{\xi }_{2j}}{D_{2j}^2}+\frac{1}{{\alpha }_2}\sum _{j=1}^m{y}_j{e}^{-\sigma y_j}(1+s_j{\phi }_{2j}),\hfill \end{array}$$

where

$${\acute{\Psi }}_{kt}=\frac{(1+e^{-\sigma z_{t-1}})}{e^{\sigma z_{t-1}}{\alpha }_k^{2+e^{-\sigma z_{t-1}}}}-\frac{(1+e^{-\sigma z_t})}{e^{\sigma z_t}{\alpha }_k^{2+e^{-\sigma z_t}}},$$

$${\acute{\xi }}_{kt}=\frac{z_t^{2}log\left({\alpha }_k\right)(e^{-\sigma z_t}log\left({\alpha }_k\right)-1)}{e^{\sigma z_t}{\alpha }_k^{e^{-\sigma z_t}}}-\frac{z_{t-1}^{2}log\left({\alpha }_k\right)(e^{-\sigma z_{t-1}}log\left({\alpha }_k\right)-1)}{e^{\sigma z_{t-1}}{\alpha }_k^{e^{-\sigma z_{t-1}}}}$$

and

$${\nu }_{kt}=\frac{z_t(e^{-\sigma z_t}log\left({\alpha }_k\right)-1)}{e^{\sigma z_t}{\alpha }_k^{1+e^{-\sigma z_t}}}-\frac{z_{t-1}(e^{-\sigma z_{t-1}}log\left({\alpha }_k\right)-1)}{e^{\sigma z_{t-1}}{\alpha }_k^{1+e^{-\sigma z_{t-1}}}}.$$

Thus, the $100(1-\tau )$ ACIs for ${\alpha }_1,{\alpha }_2$, and $\sigma$ can be obtained, respectively, as follows:

$${\tilde{\alpha }}_1\pm z_{\tau /2}\sqrt{{\tilde{V}}_{11}},{\tilde{\alpha }}_2\pm z_{\tau /2}\sqrt{{\tilde{V}}_{22}}\mathrm{and}\tilde{\theta }\pm z_{\tau /2}\sqrt{{\tilde{V}}_{33}}.$$

(25)

On the other hand, based on the delta method, one can approximate the estimate of the variance for $\tilde{R}$ using the MPSEs from (14) as follows:

$${\tilde{V}}_R={\tilde{R}}_{{\alpha }_1}^2{\tilde{V}}_{11}+{\tilde{R}}_{{\alpha }_2}^2{\tilde{V}}_{22},$$

where ${\tilde{R}}_{{\alpha }_1}$ and ${\tilde{R}}_{{\alpha }_2}$ as given by (15) and (16), respectively, but evaluated at the MPSEs. Therefore, the $100\%(1-\tau )$ ACI of R takes the form

$$\tilde{R}\pm z_{\tau /2}\sqrt{{\tilde{V}}_R},$$

4. Bootstrap Confidence Intervals

In this section, we look at two parametric bootstrap confidence intervals. The bootstrap confidence intervals are used because the ACIs can be inaccurate, especially when the sample size is small. The first one is the percentile bootstrap confidence interval (PBCI), which was developed by Efron [32]. The second is the student-t bootstrap confidence interval (SBCI) by Hall [33]. It should be noted that we obtain these two bootstrap confidence intervals using the MLEs and MPSEs. To compute the PBCIs and SBCIs, we can utilize the next procedures. It is essential to note that the following processes are used to obtain the PBCIs and SBCIs based on the MLEs. On the other side, the same processes can be used to obtain the needed intervals based on the MPSEs.

• (A) PBCIs

  (1)

  Use the original data $n,m,x_i,y_j,r_i$, and $s_j$ to calculate ${\widehat{\alpha }}_1,{\widehat{\alpha }}_2,\widehat{\sigma }$, and $\widehat{R}$.

  (2)

  Based on the same $r_i$ and $s_j$ and the estimates ${\widehat{\alpha }}_1,{\widehat{\alpha }}_2$, and $\widehat{\sigma }$, generate two progressively Type-II censored samples.

  (3)

  Use the simulated bootstrap samples in step (2) to compute ${\widehat{\alpha }}_1^{*},{\widehat{\alpha }}_2^{*},{\widehat{\sigma }}^{*}$, and ${\widehat{R}}^{*}$.

  (4)

  Repeat steps 2 and 3 B times to compute ${\widehat{\alpha }}_1^{*\left(i\right)}$, ${\widehat{\alpha }}_2^{*\left(i\right)}$, ${\widehat{\sigma }}^{*\left(i\right)}$, and ${\widehat{R}}^{*\left(i\right)}$, $i=1,\dots ,B$.

  (5)

  Arrange the estimates in (t) to obtain (${\widehat{\alpha }}_1^{*\left[1\right]},\dots ,{\widehat{\alpha }}_1^{*\left[B\right]}$), (${\widehat{\alpha }}_2^{*\left[1\right]},\dots ,{\widehat{\alpha }}_2^{*\left[B\right]}$), (${\widehat{\sigma }}^{*\left[1\right]},\dots ,{\widehat{\sigma }}^{*\left[B\right]}$), and (${\widehat{R}}^{*\left[1\right]},\dots ,{\widehat{R}}^{*\left[B\right]}$).

  (6)

  For any parameter, say $\theta$, the $100(1-\tau )\%$ PBCI is obtained as follows:

  $$\left\{{\widehat{\theta }}^{*[B\tau /2]},{\widehat{\theta }}^{*\left[B\right(1-\tau /2\left)\right]}\right\}.$$

• (B) SBCIs

  (1–4)

  As displayed in the PBCIs.

  (5)

  Obtain the statistics $T_1^{*\left(i\right)}=\frac{{\widehat{\alpha }}_1^{*\left(i\right)}-{\widehat{\alpha }}_1}{\sqrt{{\widehat{V}}_{11}\left({\widehat{\alpha }}_1^{*\left(i\right)}\right)}}$, $T_2^{*\left(i\right)}=\frac{{\widehat{\alpha }}_2^{*\left(i\right)}-{\widehat{\alpha }}_2}{\sqrt{{\widehat{V}}_{22}\left({\widehat{\alpha }}_2^{*\left(i\right)}\right)}}$, $T_3^{*\left(i\right)}=\frac{{\widehat{\sigma }}^{*\left(i\right)}-\widehat{\sigma }}{\sqrt{{\widehat{V}}_{33}\left({\widehat{\sigma }}^{*\left(i\right)}\right)}}$, and $T_4^{*\left(i\right)}=\frac{{\widehat{R}}^{*\left(i\right)}-\widehat{R}}{\sqrt{{\widehat{V}}_R\left({\widehat{R}}^{*\left(i\right)}\right)}}$, $i=1,\dots ,B$.

  (6)

  Arrange the values in step 5 to get ($T_1^{*\left[1\right]},\dots ,T_1^{*\left[B\right]}$), ($T_2^{*\left[1\right]},\dots ,T_2^{*\left[B\right]}$), ($T_3^{*\left[1\right]},\dots ,T_3^{*\left[B\right]}$), and ($T_4^{*\left[1\right]},\dots ,T_4^{*\left[B\right]}$).

  (7)

  The $100(1-\tau )$ SBCIs of ${\alpha }_1,{\alpha }_2,\sigma$ and R are given, respectively, as

  $$\begin{array}{ccc}& & \left\{{\widehat{\alpha }}_1+T_1^{*[B\tau /2]}\sqrt{{\widehat{V}}_{11}\left({\widehat{\alpha }}_1\right)},{\widehat{\alpha }}_1+T_1^{*\left[B\right(1-\tau /2\left)\right]}\sqrt{{\widehat{V}}_{11}\left({\widehat{\alpha }}_1\right)}\right\},\hfill \\ & & \left\{{\widehat{\alpha }}_2+T_2^{*[B\tau /2]}\sqrt{{\widehat{V}}_{22}\left({\widehat{\alpha }}_2\right)},{\widehat{\alpha }}_2+T_2^{*\left[B\right(1-\tau /2\left)\right]}\sqrt{{\widehat{V}}_{22}\left({\widehat{\alpha }}_2\right)}\right\},\hfill \\ & & \left\{\widehat{\sigma }+T_3^{*[B\tau /2]}\sqrt{{\widehat{V}}_{33}\left(\widehat{\sigma }\right)},\widehat{\sigma }+T_3^{*\left[B\right(1-\tau /2\left)\right]}\sqrt{{\widehat{V}}_{33}\left(\widehat{\sigma }\right)}\right\}\hfill \\ \hfill \mathrm{and}\\ & & \left\{\widehat{R}+T_4^{*[B\tau /2]}\sqrt{{\widehat{V}}_R\left(\widehat{R}\right)},\widehat{R}+T_4^{*\left[B\right(1-\tau /2\left)\right]}\sqrt{{\widehat{V}}_R\left(\widehat{R}\right)}\right\}.\hfill \end{array}$$

5. Simulation Study

In this section, we perform a simulation study to investigate the performance of parameter estimates and reliability estimates for the progressively censoring stress–strength model. The random variables X and Y in this model are distributed as the APE distribution with a common scale parameter but different shape parameters. The ML and MPS estimation methods are used to obtain the estimates of parameters and reliability. For performance comparison between MLE and MPS methods, we choose different sample sizes of X and Y denoted by N and M, respectively. The shape parameters are set as ${\alpha }_1=2.0$ for X and ${\alpha }_2=0.8$ for Y, and the scale parameter is $\sigma =1.5$. Substituting the shape parameters and the scale parameter into Equation (3), we have the reliability $R=0.58$. Let n and m be the failure numbers of the observations for X and Y. For each observation time, $r_i,i=1,2,\dots ,n$ and $s_j,j=1,2,\dots ,m$ are the progressive censoring schemes for X and Y, respectively. Considering that the $r_i$ for X and $s_i$ for Y are different, in the following, we show two general progressive censoring schemes for $r_i$ and $s_i$, where $p=8\%$ units are censored.

Scheme 1: 

$R_1=r,R_2=R_3=\dots =R_n=0$.

Scheme 2: 

$R_{i_1}=R_{i_2}=\dots =R_{i_n}=1$, and the others $R_i=0$ for $i=\{1,2,\dots ,n\}/\{i_1,i_2,\dots ,i_n\}$.

In the schemes above, $n=N\times (1-p),{\sum }_{i=1}^n{r}_i=N\times p$ and $i_1,i_2,\dots ,i_n$ is the subset randomly selected from the set $\{1,2,\dots ,n\}$. The same censoring schemes are used for Y. We choose the sample size of X as $N=\{30,50,100,150\}$ and the sample size of Y as $M=\{30,40,100,200\}$. We consider the sample pairs denoted by $(N,M)$, covering the cases such as $N<M$, $N>M$ and $N=M$. In the progressively censoring stress–strength reliability model, we conduct the simulation study given the sample pair $(N,M)$ and the specific progressive censoring schemes for X and Y.

We consider the absolute bias (ABias), absolute relative error (ARE), and mean square error (MSE) for point estimates, and the interval length (IL) and coverage probability (CP) for interval estimates at the confidence level 95% to show the estimation performance based on different point and interval estimation methods. The ABias is defined as $|\widehat{\theta }-\theta |$, and the ARE is defined as $\frac{|\widehat{\theta }-\theta |}{\theta }$ given the true parameter $\theta$ and the parameter estimate $\widehat{\theta }$. The estimation performance of model parameters $({\alpha }_1,{\alpha }_2,\sigma )$ and reliability $R=P(Y<X)$ using the MLE and MPS methods, and ACI, PBCI, and SBCI methods are shown in Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5 under Scheme 1 and Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10 under Scheme 2.

From Figure 1 under Scheme 1 and Figure 6 under Scheme 2, we see that the two estimation methods perform well for the APE stress–strength model under the progressively censoring scheme. The ABias, MSE, and ARE are decreasing with the increasing N and M for the specific schemes for X and Y. This indicates that the performance becomes better when $N+M$ is increased. Comparing the estimates under the same sample size pair but different progressively censoring schemes, we observe that the point estimates are slightly affected by the choices of schemes. This indicates the estimation methods are stable to analyze the progressively censoring data in the stress–strength model. For the model parameters, the MLE has lower ABias, MSE, and ARE compared with MPSE. We also observe that the MPSE of the reliability parameter R has the smaller ABias, MSE, and ARE compared with the MLE from Figure 2 and Figure 7. Comparing the ILs of the interval estimates in Figure 3 and Figure 4 under Scheme 1 and Figure 8 and Figure 9 under Scheme 2, we can see that the ILs of the parameters and the reliability based on the MPS method are lower than those based on the MLE for the specific progressively censoring schemes. The interval estimates have minor differences with the choices of schemes of X and Y. The ILs are decreasing with the increasing N and M. For the reliability parameter R, it is seen that the interval estimates obtained based on the MPSEs have the smaller ILs. From Figure 5 under Scheme 1 and Figure 10 under Scheme 2, we see that the CPs of confidence intervals for model parameters and reliability are around 95%, and the CPs of the model parameters using the MLEs are greater than those using the MPSEs.

6. Real Data Applications

In this part, two applications to real data sets are considered to demonstrate the importance and applicability of the proposed approaches.

6.1. Recurrence Times to Infection for Kidney Patients

The first real data consist of two groups of recurrence times from McGilchrist and Aisbett [34], which are the times to infection at the point of insertion of the catheter for each kidney patient using portable dialysis equipment, short for KID data. In the original data set, the event types are 1 for infection occurring and 0 for censoring. See McGilchrist and Aisbett [34] for a more detailed description. To use the introduced APE stress–strength reliability model under the progressively censoring scheme, we first transform the data set by dividing by 300 to fit the APE distribution with stable distribution parameters when the maximum likelihood method and the maximum product of spacing method are used, that is, $X/300$ and $Y/300$ for each patient. The transformed data set is shown in

Table 1. Transformed KID data.

Number	X	${\mathit{r}}_{\mathit{i}}$	Y	${\mathit{s}}_{\mathit{i}}$	Number	X	${\mathit{r}}_{\mathit{i}}$	Y	${\mathit{s}}_{\mathit{i}}$
1	0.0267	0	0.0533	0	20	0.0500	0	0.3600	1
2	0.0767	0	0.0433	1	21	0.5067	0	1.8733	0
3	0.0733	0	0.0933	0	22	1.3400	0	0.0800	1
4	1.4900	0	1.0600	0	23	0.0433	0	0.2200	0
5	0.1000	0	0.0400	0	24	0.1300	0	0.1533	1
6	0.0800	0	0.8167	0	25	0.0400	0	0.1333	0
7	0.0233	0	0.0300	0	26	0.3767	1	0.6700	0
8	1.7033	0	0.1000	0	27	0.4400	0	0.5200	0
9	0.1767	0	0.6533	0	28	0.1133	0	0.1000	0
10	0.0500	0	0.5133	0	29	0.0067	0	0.0833	0
11	0.0233	0	1.1100	0	30	0.4333	0	0.0867	0
12	0.4700	0	0.0267	1	31	0.0900	0	0.1933	0
13	0.3200	0	0.1267	0	32	0.0167	1	0.1433	0
14	0.4967	1	0.2333	1	33	0.5067	0	0.1000	0
15	1.7867	0	0.0833	1	34	0.6333	0	0.0167	1
16	0.0567	0	0.0133	1	35	0.3967	0	0.0267	0
17	0.6167	0	0.5900	0	36	0.1800	1	0.0533	1
18	0.9733	0	0.3800	0	37	0.0200	1	0.2600	0
19	0.0733	1	0.5300	1	38	0.2100	0	0.0267	1

where $r_i$ and $s_i$ denote the progressively censoring schemes for X and Y, respectively. We then compare two such recurrence times for each patient to show the time difference by estimating the reliability, such as the case of $Y<X$ for the transformed data set. We mention that the same transformation in the two groups of recurrence times would not affect the reliability estimates. Specifically, given an increasing function $g(·)$ on X and Y, we have $R=P(Y<X)=P\left(g\right(Y)<g(X\left)\right)$.

To show the goodness-of-fit of APE distribution, the Kolmogorov–Smirnov (KS) test is used for the observations of X and Y in the KID data set. Parameter estimates (shape and scale parameters) based on the two mentioned methods, K-S value and p-value of the fitted distribution for stress–strength random variables (RV) using the transformed KID data are presented in

Table 2. Estimates, K-S value and p-value of X and Y in the transformed KID data.

Method	RV	Shape	Scale	K-S	p-Value
MLE	X	0.2081	3.6256	0.1464	0.5160
	Y	0.2962	3.6256	0.1493	0.6059
MPS	X	0.1036	2.5161	0.1277	0.6846
	Y	0.1921	2.5161	0.1185	0.8503

. We also present the quantile–quantile (Q-Q) plot and the estimated CDF of the APE distribution, given the progressively censored data of X and Y on the left side of Figure 11 based on the ML method and Figure 12 using the MPS method. For performance comparison with APE distribution, we have the estimated CDFs using APE distributions and the empirical CDF using the generalization of the produce-limit method by Michael and Schucany [35] for these progressively censoring KID data on the right side of Figure 11 and Figure 12. For each variable in Figure 11 and Figure 12, we see that most observations are close to the straight line in the Q-Q plots, and the estimates of CDFs perform well for the two methods using APE distribution in the stress–strength model.

Based on the ML and MPS methods, we obtain the point and interval estimates of model parameters and the reliability R for these KID data shown in

Table 3. Point and interval estimates for KID data.

Method	Para	Estm	ACI	PBCI	SBCI
MLE	${\alpha }_1$	0.2081	(0, 0.7097)	(0.0048, 2.0938)	(0, 1.8176)
	${\alpha }_2$	0.2962	(0, 1.043)	(0.004, 3.5501)	(0, 2.8458)
	$\sigma$	3.6256	(0.8434, 6.4078)	(1.3919, 8.3672)	(0.1653, 7.1406)
	R	0.4724	(0.3855, 0.5788)	(0.3202, 0.5697)	(0.3508, 0.6004)
MPS	${\alpha }_1$	0.1036	(0, 0.5072)	(0.0036, 0.6737)	(0.0025, 0.6726)
	${\alpha }_2$	0.1921	(0, 0.9703)	(0.006, 1.6753)	(0, 1.5647)
	$\sigma$	2.5161	(0, 5.7229)	(1.0222, 4.9372)	(1.2234, 5.1384)
	R	0.4544	(0.3643, 0.5668)	(0.3068, 0.535)	(0.3382, 0.5665)

. We observe that both the reliability estimates are less than 0.5 using the ML and MPS methods. But they are near 0.5. This implies that a minor difference exists between the two groups of recurrence times, and the second group of recurrence times Y is slightly larger than the first group of times X. For further illustration, we simulate KID data using the MLEs and MPSEs of parameters in Table 3, respectively. Based on the simulated KID data, the parameter and reliability estimates are obtained in

Table 4. Parameter and reliability estimation for simulated KID data using MLEs.

Method	Para	Estm	SE	ACI	PBCI	SBCI
MLE	${\alpha }_1$	0.6234	0.4095	(0, 1.8753)	(0.0727, 2.5088)	(0, 2.3308)
	${\alpha }_2$	1.7012	0.8346	(0, 4.8596)	(0.4902, 10.1192)	(0, 8.6733)
	$\sigma$	3.3854	0.8927	(1.5711, 5.1997)	(2.3328, 5.5307)	(1.9076, 5.1055)
	R	0.4170	0.0641	(0.337, 0.5161)	(0.2728, 0.5064)	(0.3017, 0.5353)
MPS	${\alpha }_1$	0.2875	0.1958	(0, 1.0799)	(0.0035, 0.9391)	(0.1063, 1.0419)
	${\alpha }_2$	0.8581	0.4699	(0, 2.9505)	(0.01, 4.281)	(0, 4.2514)
	$\sigma$	2.6094	0.7945	(0.4889, 4.73)	(0.7402, 4.1098)	(1.2233, 4.5929)
	R	0.4112	0.0615	(0.3326, 0.5084)	(0.265, 0.5136)	(0.294, 0.5427)

and

Table 5. Parameter and reliability estimation for simulated KID data using MPSEs.

Method	Para	Estm	SE	ACI	PBCI	SBCI
MLE	${\alpha }_1$	0.3333	0.9462	(0, 1.0616)	(0.0101, 2.2839)	(0, 2.0505)
	${\alpha }_2$	0.7203	1.3071	(0, 2.2618)	(0.0371, 6.1638)	(0, 5.1807)
	$\sigma$	1.9772	0.8694	(0.6885, 3.2659)	(0.8767, 3.7931)	(0.4772, 3.3936)
	R	0.4371	0.0567	(0.3547, 0.5388)	(0.2863, 0.5283)	(0.3177, 0.5597)
MPS	${\alpha }_1$	0.0828	0.4884	(0, 0.5641)	(0.003, 0.5729)	(0, 0.5633)
	${\alpha }_2$	0.1923	0.6476	(0, 1.224)	(0.0077, 1.5637)	(0, 1.4653)
	$\sigma$	1.2133	0.6819	(0, 3.4036)	(0.5192, 2.411)	(0.5969, 2.4887)
	R	0.4387	0.0534	(0.3577, 0.5381)	(0.2966, 0.527)	(0.3224, 0.5527)

, where SE indicates the standard error. We see that the estimated parameters using the simulated KID data are close to the estimated parameters using the real KID data by comparing Table 3 and Table 4, and Table 3 and Table 5.

6.2. Breaking Strengths of Jute Fibre Data

We apply the progressively censored stress–strength model to analyze another real data set. This data set is the breaking strengths of jute fiber at gauge lengths 10 mm and 20 mm, short for JFG data; see, for more detail, the work of Xia et al. [36]. Similarly, with the transformation on the KID data set, we show the transformed JFG data in

Table 6. Transformed breaking strengths for X and Y for JFG data.

X	1.3875	1.4093	0.6477	1.5563	0.2461	1.2753	0.7669	0.3030	0.2179	0.1003
	1.3430	0.3663	0.5149	1.4545	0.5825	0.2023	0.7528	0.3268	0.2828	1.4015
	0.5258	0.7065	0.8442	0.0879	1.1810	0.4243	0.6078	1.0132	1.0611	0.3545
Y	0.1429	0.8380	0.5693	1.1711	0.9132	0.2277	0.3757	1.3763	1.3253	0.0912
	1.1572	1.5134	1.1886	0.3330	0.1994	1.4147	1.5303	0.3743	0.2919	0.7014
	1.0949	0.2340	0.7516	1.1632	0.2397	0.0960	0.4003	0.0735	0.4891	0.1671

by dividing by 500. Here, the transformed breaking strengths at gauge length 10 mm are the observations of X under the progressively censoring scheme $r_i=0$ for X and at gauge length 20 mm for Y under the progressively censoring scheme $s_i=0$ for Y, where $i=1,2,\dots ,30$. Therefore, the censoring schemes are not shown in Table 6.

For these JFG data, the estimates and the KS testing values are presented in

Table 7. Estimates, K-S value and p-value of X and Y in the transformed JFG data.

Method	RV	Shape	Scale	K-S	p-Value
MLE	X	16.9289	2.2980	0.0849	0.9692
	Y	7.5969	2.2980	0.1593	0.3903
MPS	X	10.0976	2.0821	0.0961	0.9200
	Y	5.0576	2.0821	0.1428	0.5266

to show the goodness-of-fit of each APE distribution for the stress–strength random variables. The Q-Q plot and the estimated CDFs are also shown in Figure 13 based on the ML method and in Figure 14 using the MPS method. We observe that the CDF estimation performance of APE distribution is good for each variable, and each estimation method is compared with the empirical CDF. This also implies that the estimates of distribution parameters are reasonable using the introduced progressively censoring stress–strength reliability model with APE distribution.

In the stress–strength model, the reliability $R=P(Y<X)$ is estimated given the transformed observations of X and Y given in Table 6. To compare the difference between the breaking strengths at gauge lengths 10 mm and 20 mm, the reliability is defined to measure the difference in a parametric model framework on the censored data. Nadeb et al. [37] used an exponentiated Fréchet distributed stress–strength reliability model to analyze these JFG data to show the difference. Here, we choose the APE distribution to fit this data set. The point and interval estimates of distribution parameters and the reliability R are shown in

Table 8. Point and interval estimates for JFG data.

Method	Para	Estm	ACI	PBCI	SBCI
MLE	${\alpha }_1$	16.9289	(0, 42.034)	(3.1314, 36.8634)	(0, 31.9271)
	${\alpha }_2$	7.5969	(0, 17.9733)	(1.3987, 32.4363)	(0, 28.102)
	$\sigma$	2.2980	(1.8713, 2.7246)	(1.739, 2.828)	(1.7014, 2.7904)
	R	0.5558	(0.4612, 0.6697)	(0.4198, 0.6731)	(0.4257, 0.6789)
MPS	${\alpha }_1$	10.0976	(0, 24.0585)	(0.9452, 29.2365)	(1.4248, 29.7162)
	${\alpha }_2$	5.0576	(0, 11.635)	(0.594, 19.8588)	(0.7678, 20.0326)
	$\sigma$	2.0821	(1.6764, 2.4878)	(1.2571, 2.4602)	(1.4731, 2.6762)
	R	0.5510	(0.457, 0.6643)	(0.4158, 0.6755)	(0.4156, 0.6753)

based on the ML and MPS methods. From Table 8, we observe that the difference of point estimates between the two methods is small, but the MPS method is better than the ML method in terms of the interval lengths. Both the R estimates are larger than 0.5. The reliability estimates are also higher than 0.5 in the work of Nadeb et al. [37]. Similarly, simulated JFG data are generated using the MLEs and MPSEs of the parameters in Table 8, respectively. The parameters and reliability estimates are obtained in

Table 9. Parameter and reliability estimation for simulated JFG data using MLEs.

Method	Para	Estm	SE	ACI	PBCI	SBCI
MLE	${\alpha }_1$	13.4211	10.4318	(0, 36.4709)	(2.0164, 33.6128)	(0, 30.8833)
	${\alpha }_2$	7.0213	7.8096	(0, 24.6551)	(0.7734, 30.4276)	(0, 26.1548)
	$\sigma$	2.1080	0.3062	(1.8108, 2.6475)	(1.4328, 2.5851)	(1.4572, 2.6095)
	R	0.5475	0.0684	(0.466, 0.676)	(0.3974, 0.6645)	(0.4149, 0.6819)
MPS	${\alpha }_1$	6.7101	4.8405	(0, 21.154)	(0.0019, 30.795)	(0, 30.7761)
	${\alpha }_2$	3.8624	3.7854	(0, 9.349)	(0.0049, 14.2785)	(0, 13.335)
	$\sigma$	1.8197	0.2900	(1.6097, 2.4028)	(0.2274, 2.2947)	(0.4093, 2.4766)
	R	0.5436	0.0644	(0.4636, 0.6727)	(0.3696, 0.699)	(0.3843, 0.7138)

and

Table 10. Parameter and reliability estimation for simulated JFG data using MPSEs.

Method	Para	Estm	SE	ACI	PBCI	SBCI
MLE	${\alpha }_1$	16.5095	10.7951	(0, 29.5292)	(4.3397, 38.6413)	(0, 30.0847)
	${\alpha }_2$	12.4487	11.0307	(0, 21.1622)	(1.4595, 30.7399)	(0, 33.1656)
	$\sigma$	2.2606	0.2678	(1.6076, 2.3475)	(1.8619, 2.8978)	(1.1694, 2.601)
	R	0.5290	0.0712	(0.43, 0.6342)	(0.4508, 0.6707)	(0.3833, 0.6351)
MPS	${\alpha }_1$	8.5118	5.1865	(0, 13.8538)	(0.0089, 15.9341)	(0.3599, 16.2851)
	${\alpha }_2$	7.1433	6.9200	(0, 11.2661)	(0.0041, 14.3265)	(0.4654, 14.7878)
	$\sigma$	1.9773	0.2628	(1.3569, 2.0465)	(0.2937, 2.188)	(0.4476, 2.3419)
	R	0.5263	0.0695	(0.4233, 0.6266)	(0.3733, 0.6143)	(0.3692, 0.6102)

based on the simulated JFG data.

7. Conclusions

In this paper, we suggest several point and interval estimators for the stress–strength parameter when the alpha power exponential is the underlying distribution for both stress and strength populations and the available data are progressively Type-II censored. Two classical estimation methods are offered for the point and interval estimation of the model parameters and the stress–strength parameter. The first method is the maximum likelihood method, and the second one is the maximum product of spacing method. The point estimates of the various unknown parameters are acquired using both methods, and the associated confidence intervals are obtained using the asymptotic properties of the proposed estimation methods. Moreover, two parametric confidence intervals are studied based on the two estimation approaches. In order to evaluate the efficiency and applicability of the proposed procedures, a simulation study is conducted, and two real-life data sets are examined. Based on the numerical findings, it is seen that the maximum product of spacing point estimator of the stress–strength reliability performs better than the usual maximum likelihood estimator in terms of the minimum absolute bias, absolute relative error, and mean square error. In terms of minimum interval length, the same pattern is observed when comparing the approximate confidence intervals based on both methods. Finally, it is noted that the bootstrap confidence intervals obtained based on the maximum likelihood estimates have higher coverage probabilities than those obtained using the maximum product of spacing estimates for the model parameters. Finally, it is important to mention that the two selected data sets were used for the practical investigation, which does not necessarily imply the same connection with other data sets of this type.