Stress-strength reliability estimation for the inverted exponentiated Rayleigh distribution under unified progressive hybrid censoring with application

: In this paper, we studied the estimation of a stress-strength reliability model ( 𝑅 (cid:3404) 𝑃(cid:4666)𝑋 (cid:3408) 𝑌(cid:4667) ) based on inverted exponentiated Rayleigh distribution under the unified progressive hybrid censoring scheme (unified PHCS). The maximum likelihood estimates of the unknown parameters were obtained using the stochastic expectation-maximization algorithm (stochastic EMA). The asymptotic confidence intervals were also created. Under squared error and Linex and generalized entropy loss functions, the Gibbs sampler together with Metropolis-Hastings algorithm was provided to compute the Bayes estimates and the credible intervals. Extensive simulations were performed to see the effectiveness of the proposed estimation methods. Also, parallel to the development of reliability studies, it is necessary to study its application in different sciences such as engineering. Therefore, droplet splashing data under two nozzle pressures were proposed to exemplify the theoretical outcomes.


Introduction
In reliability study under an ideal state, a set of the complete sample can be observed, and experimenters make statistical inferences using the entire data. However, in real life, due to substantial improvement of science and technology, we cannot get the adequate number of failure times of some units which are put on the life testing experiment. Therefore, censoring schemes are put forward for practical significance. The forms of censoring schemes are varied. Type-I and Type-II schemes are the two fundamental censoring schemes. Type-I censoring refers to the removal of units that have not failed after a preset time T. The Type-II censoring is to discard the remaining units when the number of failed units reaches a preset number m. Furthermore, if Type-I and Type-II are combined together (see [1]), it is hybrid censoring scheme. The main limitation of the hybrid censoring scheme and its generalizations is that a small effective sample size may be acquired. Therefore, Balakrishnan et al. [2] proposed Type-I and Type-II unified hybrid (unified Hy) censoring schemes. However, none of these censoring schemes (CSs) allows the removal of the units while in a life-test experiment. The progressive censoring (PrC) scheme allows the removal of units and observing a certain number of failures during the experiment, which makes life testing efficient; see [3]. Later, a generalization of the PrC known as the progressive hybrid was proposed in [4]. Needless to say, there has been a great deal of research done on estimation problems of different distributions under the above censoring schemes. Panahi [5] studied the ML and Bayesian statistical inference problems of unified Hy censored model for the Burr Type III distribution. Lee and Seo [6] considered the estimates for the Gompertz parameters under PrC sample. Similar tests were considered in [7] for Weibull distribution when right CS is available. Starling et al. [8] for improving Weibull distribution using generalized Type-I CS. Also, for accelerated life-time reliability model (ALTM), Wang et al. [9] considered estimation methods under PrC. In [10] and [11] presented statistical inferences of the ALTM under adaptive type-I progressively hybrid and unified Hy censoring schemes respectively. Further, Lone et al. [12] considered a reliability model based on a unified hybrid censoring scheme. The other developments on censoring schemes can be seen in [13,14]. Moreover, the developments and needs in engineering, manufacturing and technology inspire more improved censoring schemes (CS). Recently, Gorny and Cramer [15] introduced a novel CS called the unified PHCS, which is used in reliability studies as follows: An experiment starts with n identical items. The effective number of failures and ; are predetermined, and the Prg , , . . . , with ∑ is also fixed in advance. Let and be the two time points. This experiment will be stopped at * : : , , : : , . So, the joint likelihood function of the three forms of failure times is defined as where , is , : : , , , , : : and , for cases I, II, III and IV, respectively. Figure 1 shows the structure of the unified PHCS.
The problem of estimating the stress-strength reliability model has received significant interest in reliability theory. This reliability system is introduced in [16] and is defined as the probability of the strength component (X) overcoming the stress (Y) imposed on it, namely, . Therefore, if the stress exceeds the strength, system component will fail, and the probability of not failing of the system component is proposed as a measure of system reliability. The reliability model has wide applications in different areas of scientific study, such as engineering, medicine, physics, quality control and clinical trials. In engineering studies, R gives the probability that the strength of a system component exceeds the stress on the system from external factors. In medical studies, if X represents control, and Y represents treatment outcomes, then can be interpreted as the effectiveness of treatment. In addition to usages in engineering and clinical studies, the probability R comes up in various fields of study, such as bioequivalence assessment [17], heritability of a genetic trait [18], etc. Many researchers have considered the problem of estimating R with complete and censored samples. Bhattacharyya and Johnson [19] estimated the multicomponent S-S reliability model. Rao et al. [20] considered independent inverse Rayleigh random variables as stress and strength and estimated the . Ghitany et al. [21] investigated the S-S reliability model using power Lindley distributions. Akgul and Senoglu [22] assumed that stress and strength are random variables following Weibull distributions and estimated R using ranked set sampling. Bai [24] estimated an S-S model of a consecutive k-out-of-n system with proportional hazard rate family. Asadi and Panahi [25] discussed reliability estimation and presented the application of this model in coating processes. Zhang et al. [26] investigated the multicomponent reliability model using Marshall-Olkin Weibull distribution. De La Cruz [27] and Wang et al. [28] considered the reliability estimation under unit-half-normal and generalized exponential distributions, respectively, and Zhuang et al. [29] studied the reliability based on Bayesian deep learning.
This paper discusses the problem of reliability parameters estimation using the inverted exponentiated Rayleigh (IERayleigh) distribution under Unified PHCS. The probability density function PDF f . and CDF F . , of the IERayleigh distribution are given by and The IERayleigh distribution is a widely used life distribution in reliability and life test research, especially in analyzing the data of life-testing experiments in physics, medicine, biological sciences and engineering sciences. The hazard functions ( 2 / / 1 / ) have non-monotonic forms (see Figure 2), which are quite suitable for data sets with a nonmonotonic failure rate in actual production. For example, this model can be used in evaluating data of survival after successful surgery. In this case, the initial risk increases immediately postoperatively for possible reasons, and then the risk decreases steadily as the patient recovers. Due to the non-monotonic form of the hazard rate, this distribution behaves similarly to some well-known models: namely, the log-normal distribution, the inverse Weibull distribution and the generalized inverse exponential distribution.
Since then, IERayleigh has been widely used in many fields. For example, Kayal et al. [30] showed that the IERayleigh distribution can be considered in modeling the endurance of deep groove ball bearings data. Gao and Gui [31] applied this distribution to carbon fibers data. Panahi and Moradi [32] found that the IERayleigh distribution provided a good model for engineering problems, such as the maximum spreading diameter of nano-droplet impact on hydrophobic surface. Fan and Gui [33] showed that this distribution can be applied in the African aluminum industry.
Let the model consist of strength variable with ℎ , and independent stress variable with ℎ , . Then, the S-S reliability model can be written as where, . and . are defined in Eqs (2) and (3), respectively. This model is widely utilized in reliability engineering to determine the probability of the system working properly. Also, the use of different types of data in reliability encourages researchers to study the application of reliability models in different sciences. In this study, we try to expand the applications of the reliability model in the engineering field. Hence, this study aims to look into the reliability model through different estimation methods based on Unified PHCS and consider its application in the splashing data. To the best of our knowledge, no one has applied the IERayleigh distribution to S-S reliability model under unified PHCS. That is why we consider the estimation of the reliability model, in the case of each component having an IERayleigh distribution with common shape parameter β, and both are exposed to unified PHCS. Another importance is the implementation of the obtained methods to the splashing diameter data under two pressures. So, this study has novelty with regard to the considered S-S reliability model under unified PHCS as well as real data application of the obtained estimates.
The S-S reliability model parameters are acquired using MLEs. The stochastic EMA is applied to attain the MLE of the S-S model. Using the normal approximation, approximate confidence intervals for the reliability model's parameters are obtained. Under the squared error and Linex and generalized entropy loss functions, the Bayesian estimates are calculated by Gibbs sampler together with the MeHa algorithm. This algorithm is considered to compute the associated credible intervals. The droplet splashing data under two nozzle pressures are compared by using the considered S-S reliability model. The rest of this paper is organized as follows. The stochastic EMA and ApCIs estimates for the reliability model are derived in Section 2. In Section 3, the Bayesian estimates by the Gibbs sampler together with the MeHa algorithm under three loss functions are proposed. In Section 4, we offer a simulation study, under which the above estimation methods are comparatively analyzed. In Section 5, an illustration of how the proposed model and methods may be utilized in engineering problems is presented. A summary and some conclusions are given in Section 6.  The Unified PHCS sample forY. ⎴ : : , . . . , : :

Nomenclature and Acronyms
where : : The MLEs can be obtained by taking derivatives of , and and making them equal to 0. Using the invariance property of the MLE helps us to obtain MLE of R, denoted by , replacing the parameters in Eq (4) with their estimates to be found by using and . That is, can be written as .
Theorem 2.1. Suppose the failure times come from IERayleigh distribution. The MLEs of and , given , exist and have the following expressions. Proof. See Appendix. It is observed that the MLEs of parameters cannot be solved explicitly. Therefore, we propose to use stochastic EMA, which is an appropriate method for the censored data sets.

Stochastic EMA
The Newton-Raphson algorithm is a frequently used method to derive the MLE of model parameters. In this method, the second derivatives of the log-likelihood function are needed, which may be complicated in some situations. Therefore, the EM algorithm is employed for comparison purposes. The EM algorithm has two steps: the expectation (E) step and the maximization (M) step. The E-step involves computation of the pseudo log-likelihood function. The M-step involves maximization of the pseudo log-likelihood function. The sample collected from unified PHCS may be viewed as incomplete data, and thus we proceed to explore the EM algorithm to obtain the MLEs of model parameters.
Based on E-step of the EM algorithm, we replace the missing information data with the conditional expectations of corresponding random variables. Then, the | can be shown as follows: It is observed that the calculations of the conditional expectations of Eq (8) are more complicated and take a long time. So, by borrowing the idea of Diebolt and Celeux [34], we apply the stochastic EM algorithm (stochastic EMA), in which the E-step has been replaced by the stochastic step. This algorithm is a more suitable approximate idea and less computationally expensive than the EM algorithm. To calculate the stochastic EMA, we consider * * | * * , , , Therefore, the MLEs of the parameters can be obtained by differentiating the following loglikelihood function partially with respect to the parameters:

Asymptotic confidence intervals
In this subsection, we derive the approximate asymptotic distribution of MLE to construct approximate confidence intervals for . Theoretically speaking, the variance-covariance of , , can be derived from the Fisher matrix , , , whose elements are shown as , , , , , | ; , 1,2,3 The detailed derivations are omitted to maintain brevity. It is worth mentioning that the expected values in the Fisher information matrix cannot be obtained explicitly since the distribution of the MLEs under the unified PHCS cannot be obtained explicitly. Therefore, the observed information matrix , , is applied in the asymptotic normality of the MLE. Based on regularity conditions, the estimator asymptotically is normal with mean and variance .

,
where and . Thus, the 100 1 % asymptotic confidence interval of R can be created by where is / 100 percentile of normal distribution ( 0,1 ).

Bayesian estimation
In contrast to frequentist methods, the Bayesian approaches take advantage of available data and incorporate prior information of model parameters, thereby attracting much attention in statistical inference. So, in this section, the approximate BEs of R model are attained when all the parameters , and have independent gamma distributions with parameters , ; 1,2,3, respectively, as prior distributions. The gamma distribution is very flexible, and the Jeffreys prior, which is one of the commonly used non-informative priors, is a special case of it. Further, in the case of noninformative priors, very small non-negative values of the hyper-parameters, i.e., 0.0001, are used, which are almost like Jeffreys priors, but they are proper, inversely.
Using the likelihood function (5) Based on the Linex loss function, the Bayesian estimator for R is given by Under the GeE loss function, the Bayesian estimator for R is expressed as follows: The integrals (11)-(13) are very hard to solved analytically, so we will utilize the MCMC method. Important sampling (ImS) and Metropolis-Hastings (MeHa) are two approaches of the MCMC. It is easy to calculate the BEs based on the ImS method, which is important in practice. However, if the form of the PDF becomes complicated, the conditional posterior distribution of the parameters cannot be easily attained, so the MeHa algorithm is a more appropriate choice for calculating the BEs. This algorithm is an efficient method to calculate the BEs. It can be easily conducted and help to reduce the operational complexity of high dimensional distributions. So, the BEs and the corresponding CrIs of R based on unified PHCS are computed by using the Gibbs sampler together with MeHa in [35,36] the algorithm. To apply this method, the marginal posterior PDFs of , and can be re-expressed as  It is clearly seen that Eq (16) does not show standard forms, and therefore it is not possible to sample directly by standard process. In that case, if the posterior density function is roughly symmetric, a normal distribution can be used to approximate it. So, we will use the Gibbs sampler together with MeHa technique to calculate the BE of R. The conditional posterior distributions , and are generated using Algorithm 1. We borrow the idea of Ren and Gui [37] and use the trace plots for an MCMC chain to assess the convergence of the iteration process. If the MCMC chain is mixing adequately, the trace plot will fluctuate randomly around the true value as the lag value increases. The trace plot for n = 30, m = 15, 0.7, 1.5 is displayed in Figure 3. From Figure 3, it can be seen that the estimation results are reasonable.

Algorithm 1: Gibbs sampler together with MeHa method
Step 1: Set 1 and start by using the initial values of ( , , ).
Step 8: Compute the Bayesian estimates under SqE, Linex and GeE as: where NB is the burn-in period.

Simulation study
This section is devoted to the comparative study of the proposed estimates under different unified PHCS. Different estimation methods of reliability parameters are applied using frequency and Bayesian estimations. For comparison and analysis, 10,000 cycles of simulations are repeated for the whole process. The performance of the estimates is compared based on the following criteria: 1, 2 for the stress-strength variables are assumed to be the same for both components to account for the simultaneous testing time. We can see that the increasing of T1 and T2 for other fixed values causes smaller MSE values in all cases since it provides more testing time to let more failure be observed. From the censoring scheme, it is observed that the CS III always has the worst performance. In parallel to increasing sample sizes, all estimates provide closer results to each other. The Bayesian method performs better than the frequency method obviously, and its estimates have smaller MSEs. The MSEs of the estimations based on the informative prior are smaller than the non-informative cases. The performances of the ML and Bayesian estimators are observed close to each other under non-informative prior.
Moreover, it is noted from Tables 3 and 4 that the average lengths of the ApCIs and CrIs are decreasing with larger sample sizes. The CrIs have also smaller average length than ApCIs. The coverage probabilities (CP) of interval estimates under different censoring schemes show no apparent difference. Also, the CPs of the ApCIs are closer to the significance level (0.95) than the CrIs.
From the above analysis of the results, we present the conclusion that the results of the BEs, especially under the Linex loss function, perform better than ML estimates in the proposed model environment.
Electronic Research Archive V olume 31, Issue 7, 4011-4033.   Table 4. Average lengths and CPs of 95% ACI and HPD credible intervals for R when 1.5, 2.8. Table 5. Estimations of R for the real data example under different CSs.

Real data analysis
One of the engineering applications of reliability modeling is to consider the diameter of the splashing droplets under two different atomizing pressures. The data sets are taken from Planche et al. [38], which represents the splashing degree under different atomizing pressures (MPa). We consider splashing degree for atomizing pressures of 1.0 MPa (dataset 1) and 1.5 MPa (dataset 2). The splashing droplets under 1.5 and 1.0 MPa pressures are considered as stress data X and strength data Y, respectively. If the S-S reliability is greater than 0.50, we can say that the 1.5 pressure will increase the splashing phenomenon. Moreover, if the system reliability is less than 0.50, we will think of the reverse result of the aforementioned scenario. In order to find the fitting performance of the IERayleigh distribution, we adopt the Kolmogorov-Smirnov statistic (KoSm) with the p-value. The KoSm statistics and corresponding p-values (given in brackets) are computed based on the MLEs of the model parameters. The KoSm values (p-values) for dataset I and dataset II are 0.077004 (0.8287) and 0.077926 (0.7162), respectively. Furthermore, PP and histogram plots for observations under dataset I and dataset II samples have been drawn in Figures 4 and 5, which present the same law as the numerical results. We assumed that this distribution has the same shape parameter to illustrate the proposed method. For testing that , we perform a likelihood ratio test, and the respective value is 1.786. This value cannot get the critical point . 1 3.84. Thus, both stress and strength variables can be modeled by using the IERayleigh with equal parameter.   It is observed that since all the estimates of values are greater than 0.50, the 1.5 MPa pressure should be used for decreasing the splashing phenomenon on the considered scenario. We assume different unified PHCS censoring plans as The point estimation results along with 95% ApCI and CrI are reported in Table 5 for this real data study. We observed close estimates to each other for both estimation methods. Moreover, both approximate Bayes estimates of R are similar. The CrI lengths (CrILs) are obtained smaller than ApCI lengths (ApCILs). It is observed that since all the estimates of R values are greater than 0.50, the 1.5 MPa increases the splashing phenomenon.
We also consider the diagnostic plot for showing the convergence of the simulated Markov chains. We consider the trace plot which is a plot of the iteration number i against the value of R i at each iteration. We present the trace plot of the stress-strength reliability in MCMC chains with 14,000 iterations and first 5000 time being burn-in time in Figure 6. The trace plot shows that the Markov chain fluctuates around estimation with a similar variation. Further, as we see from Figure 7, the density plot seems to have asymmetrical and unimodal shape. It has a unique peak, and it determines the mode of the distribution. We present the diagnostic plots for SCII. The plots of the other censoring schemes are similar.

Conclusions
The issue of a stress-strength reliability model for IERayleigh distribution under unified PHCS was examined in this study. The unified PHCS is a flexible progressive censoring that can be used to reduce cost and time. It provides many different options for reliability studies and helps to overcome many problems in reliability engineering. We used classical and Bayesian methods for the estimation procedures, and we observed consistent results. Since MLE of R cannot be obtained in a closed form, the stochastic EMA is proposed to estimate the unknown parameters. The asymptotic intervals are created based on the Fisher matrix. The Bayesian estimates are calculated using the Gibbs sampler together with the MeHa algorithm with three loss functions and gamma priors. In addition, the constriction of the HPD credible intervals is provided based on MCMC samples. An extensive simulation study is implemented to validate the accuracy of the produced estimators. We found that, for point and interval estimations, the Bayes estimators using MCMC method outperformed MLEs. Therefore, if previous knowledge about the data is available, one may consider employing the Bayesian approach using the Gibbs sampler together with the MeHa algorithm for all practical purposes when analyzing data; otherwise, one may turn to MLEs or Bayesian estimation based on noninformative prior. Finally, we consider a real engineering example to show the potential application of our reliability model problem. The results show that increasing the pressure reduces the splashing phenomenon and thus improves the surface quality. Finally, the considered problem in this study can be extended in different ways such as considering the correlation between the stress and strength random variables and multicomponent stress-strength reliability model.