A comparative inference on reliability estimation for a multi-component stress-strength model under power Lomax distribution with applications

1.   Introduction

The inference of the stress-strength parameter $R = P$ $(Y < X)$ is an interesting topic in reliability analysis. The stress $Y$ and the strength $X$ are considered to be random variables. Let $X$ represents the barrier's fire resistance, and $Y$ represents the severity of the fire to which the barrier is exposed. If a unit's strength exceeds the stress applied to it, it works in the simplest stress-strength paradigm. Several authors have studied the reliability estimate of a single component stress-strength version using various lifetime distributions for the stress-strength random variate. The reliability and estimation of $X$ and $Y$ have been investigated under various distributional assumptions. There is a large literature on estimating $R$ for various stress-strength distributions under various situations. For more examples see: Ghitany et al. ^[1], Chen and Cheng ^[2], Rezaei et al. ^[3] and Sharma ^[4], where the inference of reliability was performed using complete data as examples. The estimation of $R$ under various filtering schemes was established by Genç ^[5], Krishna et al. ^[6] and Babayi and Khorram ^[7]. The inference for $R$ based on upper or lower record data was introduced by Nadar and Kızılaslan ^[8], Tripathi et al.^[9] and Asgharzadeh et al. ^[10]. Under ranked set sampling data, Akgül and Şenoğlu ^[11], Akgül et al. ^[12,13] and Safariyan et al. ^[14] investigated the $R$ estimator. Al-Babtain et al. ^[15] introduced R estimation of stress-strength model for power-modified Lindley distribution. Sabry et al. ^[16] obtained stress-strength model and reliability estimation for extension of the exponential distribution. Yousef and Almetwally ^[17] obtained multi stress-strength R based on progressive first failure for Kumaraswamy model. Rezaei et al. ^[18] discussed estimation of $P$ $[Y < X]$ for generalized Pareto distribution. Kundu and Gupta ^[19] estimated $P$ $[Y < X]$ for Weibull distributions. Jose ^[20] discussed Estimation of stress-strength reliability using discrete phase type distribution. Almetwally et al. ^[21] discussed optimal plan of multi-stress-strength reliability Bayesian and non-Bayesian methods for the alpha power exponential model using progressive first failure. Kotz et al. ^[22], added a significant literature on the topic up to the year $2003$, it can be consulted for more information. In all the previous studies there was a single component in the model, while in our proposed paper we are considering a multi-component system.

A multi-component system refers to a system with more than one component. This system, which is made up of $\; k$ independent and identical strength components, operates if $s$ $(1\leq s\leq k)$ or more of the components operate at the same time. The system is subjected to stress $Y$ in its operating environment. The component's strengths, or the minimal stress required to create failure, are random variables with a distribution function that is independent and identically distributed. The $s$-out-of-$k$: $G$ system corresponds to this model. The $s$-out-of-$k$ system is used in a variety of industrial and military systems.

In a multi-component system with $k$ components, each component has independent and identically distributed (iid) random strengths $X_{1}, X_{2}, ..., X_{k}$ and each component is stressed randomly $Y$. The system would survive if and only if the strengths were greater than the stresses by at least s out of $k, (1\leq s\leq k)$. Let $Y, X_{1}, X_{2}, ..., X_{k}$ be independent random variables, with $G(y)$ being the continuous cumulative density function (cdf) of Y and $F(x)$ being the common continuous cdf of $X_{1}, X_{2}, ..., X_{k}$. Bhattacharyya and Johnson ^[23] introduced the reliability in a multi-component stress-strength (MSS) model, which is given by

Although complete sample cases have been used to derive statistical inference on multi-component stress-strength models in the reliability literature, this subject has received little attention under censored data, particularly progressive Type-II censored sample. The following sources are closely relevant to the structure of our research.

Kohansal ^[24] recently discussed estimating dependability in a multi-component stress-strength model, where data are observed using progressive Type-II censoring for the Kumaraswamy distribution. For the general class of inverse exponentiated distributions, Kizilaslan ^[25] studied the classical and Bayesian estimate of reliability in a multi-component stress-strength model based on complete data. Gunasekera ^[26] evaluated the reliability of a multi-component system using progressively Type-II censored sample with uniformly random removals. When the common parameter is known, the author obtains different interval inferences. Estimation is evaluated in the cases in which the common parameter is either known or unknown. Gunasekera ^[26] performed inference for known common parameter cases. In the Bayesian situation, we develop a uniformly-minimum-variance-unbiased estimator and use the Tierney-Kadane approximation technique. It should be noted that Gunasekera ^[26] did not consider the case of an unknown common parameter.

Progressive Type-II censoring is commonly used in life-testing trials to analyze data under time and cost restrictions. Type-I and Type-II are the two most common censoring techniques. When using Type-I censoring, a test is almost finished at a specific time, however when using Type-II censoring, the test is finished after a certain amount of failure times have been logged. These censoring systems prevent live units from being removed between studies. The following is a description of the progressive Type II censoring scheme. On the life test, the experimenter first arranges $N$ independent and identical units. When the first failure occurs, say at $t_{(1)},$ $r_{1}$ units are eliminated at random from the remaining $N-1$ surviving units. When the second failure occurs at $t_{(2)},$ $r_{2}$ units are eliminated at random from the remaining $N-r_{1}-2$ surviving units. When the $n^{th}$ failure occurs at time $t_{n}$, the experiment ends, and the $r_{n} = N-n-\sum_{i = 1}^{n-1}r_{i}$ surviving units are eliminated from the test. For various uses of this censoring in lifespan analysis, see Balakrishnan and Aggarwala ^[27] and Balakrishnan and Cramer ^[28]. For some useful implications on this censoring scheme, see Raqab and Madi ^[29], Wu et al. ^[30] and Rastogi and Tripathi ^[31].

In this study, we attempt to estimate $R_{s, k}$ when the underlying distribution is the power Lomax distribution (POLO). According to Rady et al. ^[32], the power Lomax POLO distribution is obtained by the power transformation $X = Y^{\frac{1}{\beta }}$, where the random variable $X$ has pdf in Eq ($1.2$). The pdf of the POLO distribution is defined by

The corresponding cumulative distribution function (CDF) and survival function of POLO distribution are given by

and

where $\alpha \;$ and $\beta \;$ are shape parameters and $\lambda \;$is a scale parameter.

Due to physical constraints in these sectors, such as limited power supply, maintenance resources, and/or system design life, we suggest the POLO distribution for modeling reliability and life-testing data sets. The survival function studies the possibility of breakdowns of organisms, technical units, and other systems failing beyond a certain point in time. The hazard rate is used to measure a unit's lifetime over the duration of its lifetime distribution. The hazard rate (HRF) is a significant criterion for determining lifetime distributions since it measures the probability of failing or dying based on the age attained.

When both stress and strength follow the POLO distribution, the major goal of this study is to estimate $R_{s, k}$ using both classical and Bayesian techniques. This paper studies the feasibility of estimating $R_{s, k}$ by using maximum likelihood under progressive Type-II censoring and constructing an asymptotic confidence interval when all the parameters are unknown, these are discussed in Section 1. In Section 2, the maximum product of spacing methodology is used to obtain the explicit estimator of $R_{s, k}$. Section 3 includes the construction of boot-p and boot-t confidence intervals. In Section 4, the Bayes estimates are determined under a squared error loss function (SELF) and a linear-exponential loss function (LINEX) using gamma informative priors. The Markov-Chain Monte-Carlo (MCMC) method is obtained for Bayesian computation. Also, the Bayesian credible and HPD credible intervals are constructed. A simulation study and real data set are analyzed in Sections 5 and 6 respectively. Finally, conclusion is presented in Section 7.

2.   Estimation of $R_{s, k}$ when $\alpha \;$ is unknown

Let $X_{1}, X_{2}, ..., X_{k}$ and $Y$ be an iid random samples taken from the general class of power Lomax POLO($\alpha _{1}, \beta, \lambda$) and POLO($\alpha _{2}, \beta, \lambda$) distributions, respectively, with a common shape parameter $\beta \;$and scale parameter $\lambda$. Under this setup, the reliability $R_{s, k}$ is as follows

where $u = (1+\frac{y^{\beta }}{\lambda })^{-\alpha _{2}}.$

since $B(., .)$ is the standard Beta function and $k$ and $i$ are integers.

2.1.   Maximum likelihood estimation of $R_{s, k}$

To achieve the desired MLE of $R_{s, k}$, the MLEs of $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$ are evaluated assuming the progressive Type-II censoring scheme. Suppose $N$ systems are employed in a life-testing experiment, with a progressive Type-II censored sample $\{ X_{i1}, X_{i2}, ..., X_{ik} \}$, $i = 1, 2, \ldots, n$ is generated from the general class of power Lomax POLO($\alpha _{1}, \beta, \lambda$), where the progressive censoring scheme is $\{ K, k, r_{1}, \ldots, r_{k} \}$. Consider a progressively censored sample $\{ Y_{1}, Y_{2}, ..., Y_{n} \}$ obtained from another broad class of power Lomax POLO($\alpha _{2}, \beta, \lambda$) using the censoring scheme $\{ N, n, S_{1}, \ldots, S_{n} \}$.

Then likelihood function of $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$ is obtained as

where the constants $c_{1}$and $c_{2}\;$are given by

As a result, the likelihood function is expressed as:

The log-likelihood function is as follows:

The likelihood equations are constructed with respect to the variable of interest by calculating the derivatives of Eq ($1.5$) in the following forms

and

The parameters $\alpha _{1}$ and $\alpha _{2}$ MLEs are derived from the solutions of Eqs ($1.6$) and ($1.7$), respectively:

and

By incorporating $\hat{\alpha}_{1}$ and $\hat{\alpha}_{2}$ into Eq ($1.2$), the MLE of $R_{s, k}$ becomes

2.2.   Asymptotic confidence intervals

We use the asymptotic distribution of MLE $\hat{R}_{s, k}$ to construct an asymptotic confidence interval for the multicomponent reliability $R_{s, k}$, also need to observe an asymptotic distribution of $\hat{\theta} = (\hat{ \alpha}_{1}, \hat{\alpha}_{2}, \hat{\beta}, \hat{\lambda})$. In this regard, let $E[I(\theta)]$ denote the expected Fisher-information matrix, where

The elements of this matrix are obtained as

and

The asymptotic variances (AV) of $\hat{\alpha}_{1}$ and $\hat{\alpha}_{2}$ are calculated from Fisher information as given below

and

The MLE of $R_{s, k}$ is asymptotically normal with mean $R_{s, k}$, and a corresponding asymptotic variance given by:

for more details one may refer to Rao ^[33] It should be noticed that we obtain $R_{s, k}$ and its derivatives for $(s, k) = (1, 3)$ and $(2, 4)$, independently to avoid the difficulty in deriving $R_{s, k}$.

Therefore, $100(1-\gamma)\%$ confidence interval of $R_{s, k}$ is constructed as given below

where $z_{\gamma /2}$ denotes the upper $\gamma /2$th quantile of the standard normal distribution and $\hat{V}(\hat{R}_{s, k})$ is the MLE of $V(R_{s, k})\;$ which is obtained by replacing $(\alpha _{1}, \alpha _{2}, \beta, \lambda)$ in $V(R_{s, k})$ by their corresponding MLEs.

3.   Maximum product of spacing estimation

Ng et al. ^[34] presented the MPS approach. MPS technique determines the parameter values that makes the observed data as uniform as possible, with respect to a given quantitative measure of uniformity and based on a progressively Type-II censored sample

Cheng and Amin ^[35] defined as the geometric mean of the spacing as

where

such that $\sum D_{i} = 1$, depending on MPS method that was introduced by Cheng and Amin ^[36] and progressive Type-II censored scheme that was discussed by Balakrishnan and Aggarwala ^[27] and Ng et al. ^[37]. For more application of MPS on complete samples see Abu El Azm et al. ^[38], Sabry et al. ^[39] and Singh et al. ^[40].

Then MPS of $\alpha _{1}, \alpha _{2}, \beta$, and $\lambda$ are obtained as

The natural logarithmic likelihood functions are

We partially differentiate Eq (2.3) with respect to the parameters $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$, then equate them to zero to obtain the normal equations for the unknown parameters. The estimators for $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$ can be found by solving the equations below.

and

The MPS $\hat{\alpha}_{1(MPS)}$, $\hat{\alpha}_{2(MPS)}$, $\hat{\beta}_{(MPS)}$ and $\hat{\lambda}_{(MPS)}$ can be obtained by solving simultaneously the likelihood equations

The Eqs (2.4)–(2.7), on the other hand, must be solved numerically using a nonlinear optimization approach. By incorporating $\hat{\alpha}_{1(MPS)}$ and $\hat{ \alpha}_{2(MPS)}$ into Eq (2.2), the MPS estimator of $R_{s, k}$ becomes

4.   Bootstrap confidence intervals (CIs)

According to point estimation, a parametric bootstrap interval informs us a lot about the population values of the quantity of our interest. Furthermore, CIs based on asymptotic results clearly make errors for small sample size. To determine the bootstrap CIs of $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$, two parametric bootstrap methods are explained. The percentile bootstrap (Boot-p) CIs, introduced by Efron ^[41], and the CIs known as the bootstrap-t (Boot-t), which was presented by Hall ^[42]. Boot-t was developed using a studentized 'pivot' and it requires a variance estimator for the MLE of $\alpha _{1}, \alpha _{2}, \beta,$ and $\lambda$.

4.1.   Parametric Boot-p

Step 1: Generate a bootstrap sample of size $nk$, $\left\{ x_{i1}^{\ast }, x_{i2}^{\ast }, ..., x_{ik}^{\ast }\right\} \;$from $\{x_{i1}, x_{i2}, ..., x_{ik}\}, i = 1, 2, ..., n, \;$and generate a bootstrap sample of size $n$, $\{ y_{1}^{\ast }, y_{2}^{\ast }, ..., y_{n}^{\ast } \}$ from $\{y_{1}, y_{2}, ..., y_{n}\}$. Compute the bootstrap estimate of $R_{s, k}$, say $\hat{R}_{s, k}^{\ast }$, using Eq (1.2).

Step 2: Repeat Step $1$, $NBoot$ times.

Step 3: Let $G_{1}(z) = P(\hat{R}_{s, k}^{\ast }\leq z)$ be the cumulative distribution function of $\hat{R}_{s, k}^{\ast }$. Define $\hat{R} _{s, k(boot-p)}^{\ast } = G_{1}^{-1}(z)$ for given $z$. The approximate bootstrap-p $100(1-\gamma)\%$ CI of $\hat{R}_{s, k}$, is given by

4.2.   Parametric Boot-t

Step 1: From the samples $\{x_{i1}, x_{i2}, ..., x_{ik}\}, i = 1, 2, ..., n$ and $\{y_{1}, y_{2}, ..., y_{n}\}, \;$then compute $\hat{R}_{s, k}$.

Step 2: The same as the parametric Boot-p in Step $1$.

Step 3: Compute the $T^{\ast }$ statistic defined as

where $V\left(\hat{R}_{s, k}^{\ast }\right)$can compute as in Eq (1.15).

Step 4: Repeat Steps 1–3, $NBoot\;$times.

Step 5: Let $G_{2}(z) = P(T^{\ast }\leq z)$ be the cumulative distribution function of $T^{\ast }$ for given $z$. Define $\hat{R}_{s, k(boot-t)}^{\ast } = \hat{R}_{s, k}+G_{2}^{-1}(z)\sqrt{\widehat{\sigma ^{2}\left(\hat{R} _{s, k}^{\ast }\right) }}.$ Then, the approximate bootstrap-t $100(1-\gamma)\%$ CI of $R_{s, k}$, is given by

5.   Bayes estimation

In this section, Bayesian estimates are obtained for the parameters that are assumed to be random, and the uncertainties in the parameters are described by a joint prior distribution, which has been developed before the collected failure data. The Bayesian approach is highly useful in reliability analysis because it may incorporate previous knowledge into the analysis. Bayesian estimates of the unknown parameters $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$, as well as some lifetime parameter $R_{s, k}$ under the SELF and LINEX loss function are developed. It is assumed here that the parameters $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$ are independent and follow the gamma prior distributions,

where all the hyperparameters $a_{i}\ $and $b_{i},$ $i = 1, 2, 3, 4$ are assumed to be known non negative numbers. To determine the elicit hyper-parameters of the independent joint prior (5.1), we can use ML estimates and variance-covariance matrix of MLE method. By equating mean and variance of gamma priors, the estimated of hyper-parameters can be written as

where, $L$ is the number of iteration and $\Omega$ is a vector of parameters.

Combining the likelihood function in Eq ($1.4$) with the priors in Eq ($4.1$), resulted with the posterior distribution of the parameters $\alpha _{1}, \alpha _{2}, \beta$ and $\lambda$ indicated by $\pi ^{\ast }\left(\alpha _{1}, \alpha _{2}, \beta, \lambda \mid \underline{\rm{x}}, \underline{\rm{y}}\right)$, which can be expressed as

A commonly used loss function is the SELF, which is a symmetrical loss function that assigns equal losses to overestimation and underestimation. If $\phi$ is the parameter to be estimated by an estimator $\hat{\phi}$, then the square error loss function is defined as:

Therefore, the Bayes estimate of any function of $\alpha _{1}, \alpha _{2},$ $\beta$ and $\lambda$, say $g\left(\alpha _{1}, \alpha _{2}, \beta, \lambda \right)$ under the SELF can be obtained as

where

Varian ^[43] considered the LINEX loss function $L\left(\triangle\right)$ for a parameter $\phi$ is given by

This loss function is suitable for situations where overestimation of is more costly than its underestimation. Zellner ^[44]. discussed Bayesian estimation and prediction using LINEX loss. Hence, under LINEX loss function in Eq ($4.3$), the Bayes estimate of a function $g\left(\alpha _{1}, \alpha _{2}, \beta, \lambda \right)$ is

where

The multiple integrals in Eqs ($4.3$) and ($4.6$) can not be obtained analytically. Thus, the MCMC technique can be used to generate samples from the joint posterior density function in Eq ($4.2$). In order to be able to implement the MCMC technique, we consider the Gibbs within the Metropolis-Hasting samplers procedure. The Metropolis-Hasting and Gibbs sampling are two useful MCMC methods that have been widely used in statistics.

The joint posterior density function of $\alpha _{1}, \alpha _{2},$ $\beta$ and $\lambda$ is obtained as follows:

Under the SELF and LINEX loss function, the Bayesian estimation of $R_{s, k}$ is the mean of the posterior function in Eq (4.7), which can be written as shown below

The integral given in Eq (4.8) is obviously impossible to be calculated analytically. As a result, the Bayesian estimator of $R_{s, k}$, specifically the Gibbs sampling methods, is obtained using this approach. The next subsection get across the specifics of these strategies.

5.1.   Gibbs sampling

The Gibbs sampling method is employed, which is a sub type of Monte-Carlo Markov Chain (MCMC) method, to create the Bayesian estimate of $R_{s, k}$ and the related credible interval. The idea behind this method is to use posterior conditional density functions to generate posterior samples of parameters of interest. The posterior density function of the parameters of interest is produced by Eq ($4.7$). The posterior conditional density functions of $\alpha _{1}, \alpha _{2},$ $\beta$ and $\lambda$ can be expressed as follows using this equation:

and

The conditional density function of $\alpha _{1}, \alpha _{2},$ $\beta$ and $\lambda$ cannot be obtained in the form of the well-known density functions, as shown by Eqs (4.9)–(4.12). In this case, we can utilize the Metropolis-Hasting (MH) technique, developed by Metropolis et al. ^[45], to create random-samples from the posterior density of $\alpha _{1}, \alpha _{2},$ $\beta$ and $\lambda$ using a normal proposal distribution.

The steps of Gibbs sampling are described as follows:

(1) Start with initial guess $\left(\alpha _{1}^{(0)}, \alpha _{2}^{(0)}, \beta ^{(0)}, \lambda ^{(0)}\right).$

(2) Set $l = 1.$

(3) Using the following M-H algorithm, generate $\alpha _{1}^{(l)}, \alpha _{2}^{(l)}, \beta ^{(l)}$ and $\lambda ^{(l)}\;$ from

$\pi _{1}^{\ast }\left(\alpha _{1}^{(l)}\mid \alpha _{2}^{(l-1)}, \beta ^{(l-1)}, \lambda ^{(l-1)}, \underline{\rm{x}}, \underline{\rm{y}}\right)$, $\pi _{2}^{\ast }\left(\alpha _{2}^{(l)}\mid \alpha _{1}^{(l)}, \beta ^{(l-1)}, \lambda ^{(l-1)}, \underline{\rm{x}}, \underline{\rm{y}}\right)$ $, \; \pi _{3}^{\ast }\left(\beta ^{(l)}\mid \alpha _{1}^{(l)}, \alpha _{2}^{(l)}, \lambda ^{(l-1)}, \underline{\rm{x}}, \underline{\rm{y}}\right) \;$ and $\pi _{4}^{\ast }\left(\lambda ^{(l)}\mid \alpha _{1}^{(l)}, \alpha _{2}^{(l)}, \beta ^{(l)}, \underline{\rm{x}}, \underline{\rm{y}}\right) \;$ with the normal proposal distributions

where $V\left(\alpha _{1}\right), \; V\left(\alpha _{2}\right), \; V\left(\beta \right)$ and $V\left(\lambda \right)$ can be obtained from the main diagonal in the inverse Fisher information matrix.

(4) Generate a proposal $\alpha _{1}^{\ast }\;$from $N\left(\alpha _{1}^{(l-1)}, V\left(\alpha _{1}\right) \right), \; \alpha _{2}^{\ast }\;$from $N\left(\alpha _{2}^{(l-1)}, V\left(\alpha _{2}\right) \right), \; \beta ^{\ast }\;$from $N\left(\beta ^{(l-1)}, V\left(\beta \right) \right) \;$and $\lambda ^{\ast }\;$from $N\left(\lambda ^{(l-1)}, V\left(\lambda \right) \right).$

(i) Evaluate the acceptance probabilities

(ii) Generate a $u_{1}$, $u_{2}, u_{3}\;$and $u_{4}$ from a uniform $(0, 1)$ distribution.

(iii) If $u_{1} < \eta _{\alpha _{1}}$, accept the proposal and set $\alpha _{1}^{(l)} = \alpha _{1}^{\ast }$, else set $\alpha _{1}^{(l)} = \alpha _{1}^{(l-1)}$.

(iv) If $u_{2} < \eta _{\alpha _{2}}$, accept the proposal and set $\alpha _{2}^{(l)} = \alpha _{2}^{\ast }$, else set $\alpha _{2}^{(l)} = \alpha _{2}^{(l-1)}$.

(iiv) If $u_{3}$ $< \eta _{\beta }$, accept the proposal and set $\beta ^{(l)} = \beta ^{\ast }$, else set $\beta ^{(l)} = \beta ^{(l-1)}$.

(v) If $u_{4} < \eta _{\lambda }$, accept the proposal and set $\lambda ^{(l)} = \lambda ^{\ast }$, else set $\lambda ^{(l)} = \lambda ^{(l-1)}.$

$(5)$ Compute $R_{s, k}^{(l)}$ at ($\alpha _{1}^{(l)}, \; \alpha _{2}^{(l)}, \; \beta ^{(l)}, \; \lambda ^{(l)}$).

$(6)$ Set $l = l+1.$

$(7)$ Repeat Steps (3)–(6), $N$ times and obtain $\alpha _{1}^{(l)}, \; \alpha _{2}^{(l)}, \; \beta ^{(l)}, \; \lambda ^{(l)}\;$and $R_{s, k}^{(l)}, l = 1, 2, ...N.$

$(8)$ To compute the CRs of $\alpha _{1}, \alpha _{2},$ $\beta$, $\lambda$ and $R_{s, k}$, $\psi _{k}^{\left(l\right) }, k = 1, 2, 3, 4, 5\;, \left(\psi _{1}, \psi _{2}, \psi _{3}, \psi _{4}, \psi _{5}\right) = \left(\alpha _{1}, \alpha _{2}, \beta, \lambda, R_{s, k}\right)$ as $\psi _{k}^{\left(1\right) } < \psi _{k}^{\left(2\right) }... < \psi _{k}^{\left(N\right) },$ then the $100(1-\gamma)\%$CRIs of $\; \psi _{k}$ is

The first $M$ simulated variants are discarded in order to ensure convergence and remove the affection of initial value selection. Then the selected samples are $\psi _{k}^{\left(i\right) }\;, j = M+1, ...N$, for sufficiently large $N$.

Based on the SELF, the approximate Bayes estimates of $\psi _{k}$ is given by

the approximate Bayes estimates for$\; \psi_{k}$, under LINEX loss function, from Eq ($4.6$) is

6.   Simulation

In this section random samples are generated from POLO distribution using the R-coding. The simulation experiment is carried out to determine the reliability coefficient and compare the suggested methods.

6.1.   Simulation study

The performance of the parameters and $R_{s, k}$ is compared using different sample sizes based on Monte Carlo simulation, where k = 5, and s = 2, 3, and 4. A total of 5, 000 random samples of size are $n_1 = 10, n_2 = 15, n_3 = 15, n_4 = 10, n_5 = 12$ created from the stress and strength populations and the sample size of censored sample are chosen as $\left(m_1 = 7, m_2 = 10, m_3 = 10, m_4 = 8, m_5 = 9\right)$, and $\left(m_1 = 9, m_2 = 13, m_3 = 13, m_4 = 9, m_5 = 11\right)$. This section examines some empirical data derived from Monte-Carlo simulations to see how the proposed methods perform with different sample sizes. For $(a_j, b_j); j =, ..., 4$, we may use the estimate and variance-covariance matrix of the MLE approach to elicit hyper-parameters of the independent joint prior. The estimated hyper-parameters are calculated by equating the mean and variance of gamma priors. For the random variables generating, the values of the parameters $\alpha_1, \alpha_2, \lambda$, and $\beta$ are chosen as follows:

Tables 1–6 show the simulation results of MLEs, MPS, Bayesian estimates, and interval estimations of $R_{s, k}$. All of the results are calculated using a total of 5000 simulated samples. The simulation methods are compared using the criteria of parameters estimation, the comparison is performed by calculating the Bias, the mean of square error (MSE), the length of asymptotic and bootstrap confidence intervals (L.CI) and coverage probability (CP) for each estimation method. In simulation results Tables 1–3, for each sample-size $m_i$, scheme (S), and estimator, the first four values represent the average bias and MSE of the parameters model, and the next three values represent the estimated risk for the corresponding stress-strength reliability when $k = 5$ and $s = 2, 3, 4$, respectively. In simulation results for CI Tables 4–6, for each sample size $m_i$, and scheme (S), in MLE, and MPS estimators, the first four values represent the average length of asymptotic CI (L.CI), CP, Boot-p (BP), and Boot-t (BT) of the parameters model, and the next three values represent the average length of delta CI of risk for the corresponding $R_{s, k}$ when $k = 5$ and $s = 2, 3, 4$, respectively. While, in Bayesian estimation, the average length of credible CI (L.CI) of the parameters model and risk for the corresponding stress-strength reliability when $k = 5$ and $s = 2, 3, 4$, respectively.

Numerical simulations, on the other hand, make it impossible to see in a basic sense how estimated dangers decrease with sample size. For probability, product spacing, and Bayesian estimates, we see this trend. In terms of estimated risks, the Bayesian estimates of $R_{s; k}$ perform significantly better than the MLE and MPS. We notice that the Bayesian estimate's predicted risks under the LINEX loss are often lower than those under the SELF. The Bayes estimates and their estimated risks are sometimes near to each other based on calculated findings. The average length of HPD intervals is found to be shorter than that of asymptotic confidence intervals. When the sample size is increased, the lengths of both intervals shrink. But, the bootstrap CI has the shortest length of CI.

7.   Application of real data

For demonstration purposes, the analysis of a pair of real data sets is shown. The idea is to figure out how we can create conditions in which there is a lot of droughts. We assert that there will be no excessive drought if the water capacity of a reservoir in an area in August for at least two years out of the next five years is more than the amount of water achieved in December of 2019. It's also feasible that in this case, rather than entire samples from both groups, one sees censored samples. To achieve this purpose, we used the monthly water capacity of the Shasta reservoir in California, as well as the months of August and December from 1975 to 2016. http://cdec.water.ca.gov/cgi-progs/queryMonthly/SHA contains the data. Some writers have previously used these data, including Nadar and Kizilaslan ^[46], and Kızılaslan and Nadar ^[47].

In the whole data case, assuming k = 5 and s = 2, $Y_1$ represents December 1975 capacity while $X_{11}, ..., X_{15}$ represents August capacities from 1976 to 1980. Also, $Y_2$ represents December 1981 capacity, while $X_{21}, ..., X_{25}$ represents August capabilities from 1982 to 1986. n = 7 data for $Y$ are acquired by continuing this approach till 2016. The corrected data are listed in Table 7.

To begin, we make sure that the POLO distribution can be utilized to examine the data set in Table 7. We get the MLEs of unknown parameters in Table 8. The Kolmogorov-Smirnov distance (KSD) values are also reported, along with the appropriate p-values. We can see from this table that the POLO distribution fits the data pretty well. For the data sets X and Y, the empirical cdf distribution, estimated pdf with histogram plot, and the PP-plot are given in Figures 1–3, respectively. In Figure 4, we plot the empirical cdf distribution, estimated pdf with histogram plot, and the PP-plot for $X = (x_1, x_2, x_3, x_4, x_5)$. These figures confirmed that the data have been fitted for POLO distribution.

Two distinct progressively censored samples are produced from the previous data sets for illustration purposes. In Table 9, the MLEs, MPS, and the Bayesian estimation of unknown parameters for the model have been obtained for this scheme. In Table 10, the estimation of reliability in a multi-component stress-strength model have been obtained. Figure 5 shows the trace and density plots for all parameters in the MCMC trace. Also it shows the trace and density plots for all parameters in an MCMC trace. The posterior density of MCMC results for each parameter is shown, which demonstrates a symmetric normal distribution that is identical to the proposed distribution. The convergence of the MCMC results is confirmed in Figure 5.

When the two schemes are compared using Bayesian and non-Bayesian, it is found that estimators in scheme 1 have lower standard errors than estimators in scheme 2. Also, it is found that estimators in scheme 1 have high reliability than estimators in scheme 2 and the whole scheme.

8.   Conclusions

The study discusses the multi-component stress-strength model. The reliability has been investigated where both the stress and strength variables follow the POLO distribution. To calculate the multi-component stress-strength reliability $R_{s; k}$, we apply both classical and Bayesian methods. We compute the Bayesian estimates of parameters and $R_{s; k}$ under symmetric and asymmetric loss functions by using the MH algorithm. We compute the classical estimates as MLE ad MPS of the parameter of the model and $R_{s; k}$ by using the Newton-Raphson (NR) and Markov Chain Monte Carlo algorithms. Based on the simulation analysis, we observe that the predicted risks of the proposed $R_{s; k}$ estimators show good behavior when the sample size increases. In general, as the sample size increases, the average length of the intervals decreases, therefor The average length of higher posterior density intervals is found to be shorter than that of asymptotic confidence intervals. Based on tabulated numerical results we find that the predicted risks of Bayesian estimation are often lower than the risks of the classical approaches. For illustrating the applicability of the multi-component stress strength model under POLO distribution we have examined a real life data taken from the monthly water capacity of the Shasta reservoir in California, and by using Kolmogorov-Smirnov distance (KSD)and its corresponding p-values we conclude that the above model fit the data very well. This study can further be extended by using different censoring schemes and different lifetime models.

Acknowledgments

This work was supported by the Deanship of Scientific Research, Vice Presidency for Graduate Studies and Scientific Research, King Faisal University, Saudi Arabia (Project No. GRANT421), King Faisal University (KFU), Al Ahsa, Saudi Arabia. The Authors, therefore, acknowledge technical and financial support of the Deanship of Scientific Research at KFU.

Conflict of interest

The authors declare no conflict of interest.