Estimation for inverse Weibull distribution under progressive type-Ⅱ censoring scheme

1.   Introduction

Kaller and Kamath ^[1] proposed two-parameter inverse Weibull distribution (IWD) to simulate the degradation of mechanical components of diesel engines. After that, IWD is considered an appropriate model to analyze lifetime data. For example, Abhijit and Anindya ^[2] found that IWD is superior to the normal distribution when using ultrasonic pulse velocity to measure concrete structures. Elio et al. ^[3] proposed a new model generated by appropriate mixing of IWD for modeling under extreme wind speed conditions. Langlands et al. ^[4] observed that breast cancer mortality data can be modeled and analyzed by IWD. Beyond that, IWD is widely used in reliability research. For example, Bi and Gui ^[5] considered the estimation of stress-strength reliability of IWD. They proposed an approximate maximum likelihood estimation for point and confidence interval estimations. Bayesian estimator and highest posterior density (HPD) confidence interval were derived using Gibbs sampling. Based on adaptive type-Ⅰ progressive hybrid censored scheme, Azm et al. ^[6] studied the estimation of unknown parameters of IWD when the data were competing risks data. The maximum likelihood estimation and Bayesian estimation were discussed. The asymptotic confidence intervals, the bootstrap confidence intervals and the HPD confidence intervals were derived. Then, two sets of real data were studied to illustrate maximum likelihood estimation and Bayesian estimation. Alslman and Helu ^[7] assumed that the two components were independent and identically distributed and considered the estimation of stress-strength reliability for IWD. Its estimators were derived by maximum likelihood estimation and maximum product of spacing method and compared using computer simulation. Shawky and Khan ^[8] assumed that stress and strength both followed IWD, and they focused on the multi-component stress-strength model. The estimation of reliability was obtained by maximum likelihood estimation. Monte Carlo simulation results indicate that the proposed estimating methods are effective.

The probability density function (PDF), the cumulative distribution function (CDF) and the reliability function of IWD are respectively given by

and

Here, $\eta$ is rate parameter and $\lambda$ is shape parameter. For convenience, the PDF (1.1) of IWD will be denoted by ${\text{IW}}(\eta, \lambda)$.

The assessment of product reliability often relies on the collection of life data, which can be obtained through a life test. This test involves observing whether a group of test samples fail during the test and recording their corresponding failure time. If the life test continues until all the test samples fail, the failure time can be recorded for all the samples, resulting in complete data. However, if the test stops before all the test samples fail, the data collected are called censored data. With the continuous advancement of science and technology, products are becoming more reliable and have longer lifespans. Collecting complete data in a life test can be expensive, making reliability analysis based on censored data a popular research topic among scholars. A progressive type-Ⅱ censored sample can be expressed as follows: consider an experiment where $n$ units are subjected to a life test at time zero, and the experimenter decides beforehand the number of failures to be observed, denoted by $r$. Upon observing the first failure time ${T_1}$, ${Q_1}$ out of the remaining $n - 1$ surviving units are randomly selected and removed. At the second observed failure time ${T_2}$, ${Q_2}$ out of the remaining $n - 2 - {Q_1}$ surviving units are randomly selected and removed. This process continues until the r-th failure is observed at time ${T_r}$, and the remaining ${Q_r} = n - r - {Q_1} - {Q_2} -... - {Q_{r - 1}}$ surviving units are all removed. The sample $T = ({T_1}, {T_2}, ..., {T_r})$ is referred to as a progressively type-Ⅱ censored sample of size $r$ from a sample of size $n$ with censoring scheme $Q = ({Q_1}, {Q_2}, ..., {Q_r})$.

IWD is one of the commonly used lifetime distributions in reliability estimation ^[5,6,7,8]. Most of the estimation methods are maximum likelihood estimation and Bayesian estimation, and there is a lack of research on inverse moment estimation. Therefore, this paper aims to provide three methods to compute the point estimations and construct generalized confidence intervals of unknown parameters.

The rest of this paper is arranged as follows: In Section 2, the maximum likelihood estimators (MLEs) are obtained by Newtown-Raphson method. The Lindley approximation is proposed to derive the Bayesian estimators (BEs) in Section 3. The inverse moment estimators (IMEs) and the construction of the generalized confidence intervals (GCIs) are discussed in Section 4. In Section 5, Monte Carlo simulation is conducted to evaluate the effect of these methods. Section 6 gives a set of real data as a demonstration. Finally, Section 7 gives the conclusions of this paper.

2.   Maximum likelihood estimation

In this section, we will discuss the MLEs of $\eta$, $\lambda$ and $R(t)$. Due to the nonlinear nature of the likelihood equations, the Newton-Raphson method is considered to solve likelihood equations numerically.

Let $T = ({T_1}, {T_2}, ..., {T_r})$ be the collected progressive type-Ⅱ sample under the censoring scheme $Q = ({Q_1}, {Q_2}, ..., {Q_r})$. Denote $t = ({t_1}, {t_2}, ..., {t_r})$ as the observation of $T$. We can easily get the likelihood function $l(\eta, \lambda; t)$ as follows:

where $\delta = n(n - {Q_1} - 1)(n - {Q_1} - {Q_2} - 2)...(n - {Q_1} - {Q_2} -... - {Q_{r - 1}} - r + 1)$.

Then the log-likelihood function $L(\eta, \lambda; t)$ is

Thus, the partial derivatives of $L(\eta, \lambda; t)$ with respect to $\eta$ and $\lambda$ are respectively given by

Denote the MLEs of $\eta$ and $\lambda$ as ${\hat \eta _{ML}}$ and ${\hat \lambda _{ML}}$ respectively, and they are the solutions of likelihood equations (2.5). According to the invariance of the maximum likelihood estimation, the MLE ${\hat R_{ML}}(t)$ is obtained by replacing the parameters with ${\hat \eta _{ML}}$ and ${\hat \lambda _{ML}}$. Since (2.3) and (2.4) are non-linear, the Newtown-Raphson method is considered to solve likelihood equations numerically. The elements of the Jacobi matrix are given from (2.6) to (2.9).

Define

and

where $\theta = (\eta, \lambda)$. The steps involved in Newton-Raphson iteration method for obtaining the MLEs of $\eta$ and $\lambda$ are given below.

Step 1. Pick an arbitrary starting estimate ${\theta _0}$, and desired precision $\varepsilon = {10^{ - 5}}$.

Step 2. Update ${\theta _0}$ as ${\theta _{new}} = {\theta _0} - {[H'({\theta _0})]^{ - 1}}H({\theta _0})$, where ${[H'({\theta _0})]^{ - 1}}$ is the inverse matrix of $H'({\theta _0})$.

Step 3. If $|{\theta _{new}} - {\theta _0}| \leqslant \varepsilon$, then ${\hat \theta _{ML}} = {\theta _{new}}$.

Step 4. If $|{\theta _{new}} - {\theta _0}| > \varepsilon$, set ${\theta _0} = {\theta _{new}}$ and return to Step 2.

Step 5. Repeat from Step 2 to Step 4 until the condition in Step 3 is achieved.

3.   Bayesian estimation

Statistical inference requires three kinds of information: overall information, sample information and prior information. The statistical inference based on the first two kinds of information is called classical statistics, and the statistical inference based on comprehensive consideration of the three kinds of information is called Bayes statistics. Prior information already exists before sampling, which mostly comes from experience and historical data. The distribution obtained by processing prior information is called a prior distribution. We all know that a random variable can be described by a certain distribution. Bayesian scholars believe that an unknown parameter can be regarded as a random variable. In other words, an unknown parameter can also be described by a certain distribution, that is, a prior distribution. Kundu and Howlader ^[9] considered the estimation of parameters of IWD using the Bayesian approach under the squared error loss function when the sample is a type-Ⅱ censoring sample. Akgul et al. ^[10] discussed Bayes estimation of the step-stress partially accelerated life test model with type-Ⅰ censored sample for the IWD. The point estimators were obtained using the Lindley approximation and Tierney–Kadane approximation. The credible intervals were constructed using the Gibbs sampling method. Helu and Samawi ^[11] considered the Bayesian inferences based on IWD based on progressive first-failure censored data. The point estimators were derived under three loss functions. Additionally, the estimators were calculated by Lindley approximation. Based on generalized adaptive progressively hybrid censored sample, Lee ^[12] considered the estimation of uncertainty measure for IWD. For the BE and HPD confidence interval, the Tierney-Kadane approximation and importance sampling technique were proposed.

In this section, the BEs of $\eta$, $\lambda$ and $R(t)$ are derived under the symmetric entropy (SE) loss function, scale squared error (SSE) loss function and LINEX loss function. Since there is a complex ratio of two integrals in BEs, Lindley approximation is proposed to solve this problem.

(i) The SE loss function is defined as (Xu et al. ^[13]):

(ii) The SSE loss function is defined as (Song et al. ^[14]):

(iii) The LINEX loss function is defined as (Varian ^[15]):

where $\hat \theta$ is the estimator of unknown parameter $\theta$ and $d$ is a nonnegative integer.

As an important part of Bayesian estimation, the selection of prior distribution will directly affect the final Bayesian estimation. It usually follows two rules making full use of prior information, such as empirical and historical data and being convenient for computational use. The most widely used prior distributions are mainly non-informative prior distribution, conjugate prior distribution and hierarchical prior distribution. The gamma prior belongs to the conjugate prior, which makes the computation process of estimation results easier.

Assume the prior distributions of $\eta$ and $\lambda$ are gamma prior. $\eta$ follows ${\text{Gamma}}({a_1}, {b_1})$ and $\lambda$ follows ${\text{Gamma}}({a_2}, {b_2})$.

Based on (3.4) and (3.5), the joint prior distribution is

Using (2.1) and (3.6), the posterior distribution is

where $K = {[\int_0^{ + \infty } {\int_0^{ + \infty } {{\eta ^{r + {b_1} - 1}}{\lambda ^{r + {b_2} - 1}}{e^{ - {a_1}\eta - {a_2}\lambda }}\prod\limits_{i = 1}^r {t_i^{ - \lambda - 1}{e^{ - \eta t_i^{ - \lambda }}}{{(1 - {e^{ - \eta t_i^{ - \lambda }}})}^{{Q_i}}}} } } d\eta d\lambda ]^{ - 1}}$.

Based on (3.7), the margin posterior distribution of $\eta$ is

The margin posterior distribution of $\lambda$ is

Lindley ^[16] proposed an approximation algorithm to calculate the ratio of two integrals in the form:

Here is a continuous function in $\eta$ and $\lambda$, $L(\eta, \lambda; t)$ as shown in (2.2) and $J(\eta, \lambda)$ is the logarithm of joint prior distribution (3.6). This ratio usually occurs in BE, which is why the Lindley approximation is often used to calculate the Bayesian estimator.

The expression (3.10) can be approximated by (3.11) under regularity conditions or with a large sample size. Here $A$, $B$ and $C$ are given by (3.12). ${L_{\eta \eta \eta }}$ denotes the third derivative of log-likelihood function (2.2) for $\eta$, and ${\hat L_{\eta \eta \eta }}$ represents the value of ${L_{\eta \eta \eta }}$ at $\eta = {\hat \eta _{ML}}$. ${\psi _{ij}}$ is the element of inverse matrix of $- {L_{ij}}$ ($i, j = \eta, \lambda$), and ${\hat \psi _{ij}}$ represents the value of ${\psi _{ij}}$ at $\eta = {\hat \eta _{ML}}$ and $\lambda = {\hat \lambda _{ML}}$. Other terms are defined as the same as the above rules. The detailed expressions are presented in (3.13) to (3.18).

3.1.   Bayesian estimation under SE loss function

Lemma 1. Suppose that $T$ is a random sample. The BE ${\hat \theta _{SE}}$ of unknown parameter $\theta$ under the SE loss function (3.1) for any prior distribution $\pi (\theta)$ is

where $E(\theta |T)$ and $E({\theta ^{ - 1}}|T)$ denote the posterior expectations of $\theta$ and ${\theta ^{ - 1}}$.

Proof. Based on the SE loss function (3.1), the Bayesian risk of ${\hat \theta _{SE}}$ is

To minimize $R({\hat \theta _{SE}})$, only need to minimize $E({S_1}(\theta, {\hat \theta _{SE}})\left| T \right.)$. Denote ${h_1}({\hat \theta _{SE}}) = E({S_1}(\theta, {\hat \theta _{SE}})\left| T \right.)$ for convenience.

Because

and the derivative is

The BE ${\hat \theta _{SE}}$ can be obtained by solving the equation ${h'_1}({\hat \theta _{SE}}) = 0$.

According to Lemma 1, the BEs ${\hat \eta _{SE}}$, ${\hat \lambda _{SE}}$ and ${\hat R_{SE}}(t)$ under the SE loss function can be written as:

and

From the marginal posterior distribution (3.8), the BE (3.20) may be written as

From the marginal posterior distribution (3.9), the BE (3.21) may be written as

From the posterior distribution (3.7), the BE (3.22) may be written as

It is obvious that (3.23)–(3.25) cannot be evaluated explicitly. Thus, the Lindley approximation is used to approximate them.

(i) When $W(\eta, \lambda) = \eta$, there are

Submitting $W(\eta, \lambda) = \eta$ and (3.26) in the expression (3.11), the posterior expectation $E(\eta |T)$ may be written as

(ii) When $W(\eta, \lambda) = {\eta ^{ - 1}}$, there are

Submitting $W(\eta, \lambda) = {\eta ^{ - 1}}$ and (3.28) in the expression (3.11), the posterior expectation $E({\eta ^{ - 1}}|T)$ may be written as

The BE ${\hat \eta _{SE}}$ of rate parameter $\eta$ under the SE loss function can be obtained by submitting (3.27) and (3.29) in the expression (3.20).

Using Lindley approximation, the BEs ${\hat \lambda _{SE}}$ and ${\hat R_{SE}}(t)$ under SE loss function are obtained as similar to the above steps.

3.2.   Bayesian estimation under SSE loss function

Lemma2. Suppose that $T$ is a set of simple random samples. The BE ${\hat \theta _{SSE}}$ of unknown parameter $\theta$ under the SSE loss function (3.2) for any prior distribution $\pi (\theta)$ is

Proof. The Bayesian risk of ${\hat \theta _{SSE}}$ based on SSE loss function (3.2) is

Denote ${h_2}({\hat \theta _{SSE}}) = E({S_2}(\theta, {\hat \theta _{SSE}})\left| T \right.)$, and

The derivative of ${h_2}({\hat \theta _{SSE}})$ is

Therefore, the BE ${\hat \theta _{SSE}}$ under the SSE loss function is derived by solving equation ${h'_2}({\hat \theta _{SSE}}) = 0$.

According to Lemma 2, the BEs ${\hat \eta _{SSE}}$, ${\hat \lambda _{SSE}}$ and ${\hat R_{SSE}}(t)$ under SSE loss function are presented in (3.31) and (3.33) respectively.

and

From the marginal posterior distribution (3.8), the BE (3.31) may be written as

From the marginal posterior distribution (3.9), the BE (3.32) can be written as

From the posterior distribution (3.7), the BE (3.33) can be written as

Then, Lindley approximation is used to approximate (3.34) to (3.36).

(i) When $W(\eta, \lambda) = {\eta ^{1 - d}}$, there are

Putting $W(\eta, \lambda) = {\eta ^{1 - d}}$ and (3.37) into (3.11), the posterior expectation $E({\eta ^{1 - d}}|T)$ can be written as

(ii) When $W(\eta, \lambda) = {\eta ^{ - d}}$, there are

Putting $W(\eta, \lambda) = {\eta ^{ - d}}$ and (3.39) into (3.11), the posterior expectation $E({\eta ^{ - d}}|T)$ can be written as

Hence the BE ${\hat \eta _{SSE}}$ under the SSE loss function are obtained by submitting (3.38) and (3.40) in (3.31).

Using Lindley approximation, the BEs ${\hat \lambda _{SSE}}$ and ${\hat R_{SSE}}(t)$ under SSE loss function can be obtained by the similar steps.

3.3.   Bayesian estimation under LINEX loss function

Lemma 3. Suppose that $T$ is a set of simple random samples. The BE ${\hat \theta _L}$ of unknown parameter $\theta$ under the LINEX loss function (3.3) for any prior distribution $\pi (\theta)$ is

Proof. The Bayesian risk of ${\hat \theta _L}$ based on LINEX loss function (3.3) is

Denote ${h_3}({\hat \theta _L}) = E({S_3}(\theta, {\hat \theta _L})\left| T \right.)$, and

The derivative of ${h_3}({\hat \theta _L})$ is

Therefore, the BE ${\hat \theta _L}$ under the LINEX loss function is derived by solving equation ${h'_3}({\hat \theta _L}) = 0$.

It follows from Lemma 3 that the BEs under LINEX loss function are

Subsequently, (3.42)–(3.44) can be written as

Next, the explicit expressions of these BEs are obtained by using the Lindley approximation. When $W(\eta, \lambda) = {e^{ - a\eta }}$, there are

According to Lindley's formula (3.11), the posterior expectation $E({e^{ - a\eta }}|T)$ can be written as

The BE ${\hat \eta _L}$ under LINEX loss function is derived by substituting (3.49) into (3.42). The BEs ${\hat \lambda _L}$ and ${\hat R_L}(t)$ can be acquired using a comparable method to the aforementioned steps, and therefore will not be reiterated here.

4.   Estimation based on generalized pivotal quantity

In Sections 2 and 3, the MLEs and BEs have been derived. However, we cannot obtain the explicit forms of MLEs and BEs easily by using both methods. In addition, it is necessary to select the appropriate initial values when using the Newtown-Raphson iteration method. Therefore, the generalized pivot quantity is constructed for deriving IMEs and GCIs in this Section. Compared with the maximum likelihood estimation and Bayes estimation, the inverse moment estimation is much simpler in calculation. It only needs some mathematical transformations and finally solves the equations. Wang ^[17] proposed a new method and named it as inverse moment estimation method in 1992. Additionally, the method was applied to parameter estimation of Weibull distribution. After that, inverse moment estimation has been widely used and studied. For example, Luo et al. ^[18] used the inverse third-moment method when forecasting a single time series using a large number of predictors in the presence of a possible nonlinear forecast function. Qin and Yuan ^[19] proposed an ensemble of IMEs to explore the central subspace. Based on progressive censored data, Gao et al. ^[20] proposed the pivotal inference for inverse exponentiated Rayleigh distribution. The point estimators were derived using the method that combined pivotal quantity with inverse moment estimation.

In this section, the IMEs and GCIs of $\eta$, $\lambda$ and $R(t)$ are obtained by constructing the generalized pivot quantity.

4.1.   Inverse moment estimation

Definition 1. Assume that $T$ is a random variable and $t$ is the observation of $T$, and ${\theta _1}$ is an interest parameter and ${\theta _2}$ is a nuisance parameter. A function $G(T; t, {\theta _1}, {\theta _2})$ is called a generalized pivotal quantity if it satisfies the following conditions:

(1) Given $t$, the distribution of $G(T; t, {\theta _1}, {\theta _2})$ is unrelated to both ${\theta _1}$ and ${\theta _2}$.

(2) Given $t$, the observation $G(t; t, {\theta _1}, {\theta _2})$ of generalized pivotal quantity $G(T; t, {\theta _1}, {\theta _2})$ is unrelated to ${\theta _2}$.

First, let

The distribution of ${X_i}$ is

Let ${\text{Exp}}(1)$ be the standard exponential distribution. It is obvious that ${X_i} \sim {\text{Exp}}(1)$, and ${X_r} < {X_{r - 1}} < ... < {X_1}$. Let

Then, ${S_1}, {S_2}, ..., {S_r}$ are independent of each other and ${S_i} \sim {\text{Exp}}(1)$. Denote $U = 2\sum\limits_{i = 1}^r {{S_i}}$ and $V = 2{S_1}$, so $U$ follows ${{\mathtt{χ}}^{\text{2}}}$ distribution with $2r - 2$ degrees of freedom and $V$ follows ${{\mathtt{χ}}^{\text{2}}}$ distribution with 2 degrees of freedom. Finally, let

Therefore, ${G_1}$ and ${G_2}$ are independent, ${G_1}$ follows ${\text{F}}$ distribution with $2$ and $2r - 2$ degrees of freedom and ${G_2}$ follows ${{\mathtt{χ}}^{\text{2}}}$ distribution with $2r$ degrees of freedom. According to Definition 1, ${G_1}$ is the generalized pivotal quantity of $\lambda$, but ${G_2}$ is neither a generalized pivotal quantity of $\eta$ nor $\lambda$. Since $(r - 1){(r - 2)^{ - 1}}$ is the mean of ${\text{F}}(2, 2r - 2)$ and $2r$ is the mean of ${{\mathtt{χ}}^{\text{2}}}(2r)$, Theorem 1 can be derived.

Theorem 1. Let $T = ({T_1}, {T_2}, ..., {T_r})$ be a progressive type-Ⅱ censored sample following ${\text{IW(}}\eta {\text{, }}\lambda {\text{)}}$, and let $Q = ({Q_1}, {Q_2}, ..., {Q_r})$ be the censoring scheme. Denote $t = ({t_1}, {t_2}, ..., {t_r})$ as the observation of $T$. The IME ${\hat \eta _G}$ of $\eta$ and the IME ${\hat \lambda _G}$ of $\lambda$ are determined by the following expressions. The IME ${\hat R_G}(t)$ is obtained by replacing the parameters with ${\hat \eta _G}$ and ${\hat \lambda _G}$.

4.2.   Generalized confidence interval

This section will discuss the GCIs by generalized pivotal quantity.

Lemma 3. Suppose that a set of constants ${k_i}$ ($i = 1, 2, ..., r$) satisfy $0 < {k_1} < {k_2} < ... < {k_r}$ and denote $G(\lambda) = \frac{{r(r - 1)k_r^{ - \lambda }}}{{\sum\limits_{i = 1}^r {k_i^{ - \lambda }} - rk_r^{ - \lambda }}}$.

(i) $G(\lambda)$ decreases monotonically when $\lambda > 0$;

(ii) The equation $G(\lambda) = k$ has only one solution, which $k > 0$ and $k$ is a constant.

Proof. (i) The derivative of $G(\lambda)$ is

According to $0 < {k_1} < {k_2} < ... < {k_r}$, it can be derived $0 < \ln {k_1} < \ln {k_2} < ... < \ln {k_r}$. That is $\ln {k_i} - \ln {k_r} < 0$ ($i = 1, 2, ..., r - 1$). Therefore, $G(\lambda)$ decreases monotonically when $\lambda > 0$.

(ii) Suppose that the equation $G(\lambda) = k$ has two unequal solutions, ${\lambda _1}$ and ${\lambda _2}$ respectively. Based on $G({\lambda _1}) = G({\lambda _2})$, there is

That is

Here, ${(\frac{{{k_i}}}{{{k_r}}})^{ - \lambda }}$ is monotone because of $\frac{{{k_i}}}{{{k_r}}} \geqslant 1$. Hence the above expression is inconsistent with the supposition, and the equation $G(\lambda) = k$ has only one solution.

Theorem 2. Let $T = ({T_1}, {T_2}, ..., {T_r})$ with observation $t = ({t_1}, {t_2}, ..., {t_r})$ be a progressive type-Ⅱ censored sample following ${\text{IW(}}\eta {\text{, }}\lambda {\text{)}}$, and let $Q = ({Q_1}, {Q_2}, ..., {Q_r})$ be the censoring scheme. ${{\text{F}}_\omega }(2, 2r - 2)$ and ${\mathtt{χ}}_\omega ^2(2r)$ denote the upper quantile of ${\text{F}}(2, 2r - 2)$ and ${{\mathtt{χ}}^2}(2r)$ respectively with $\omega = \frac{{1 - \sqrt {1 - \gamma } }}{2}$. The $100(1 - \gamma)\%$ GCIs of $\eta$ and $\lambda$ are given as follows:

where $\varphi ({t_1}, ..., {t_r}, {{\text{F}}_{1 - \omega }}(2, 2r - 2))$ is the solution of equation (4.8) and $\varphi ({t_1}, ..., {t_r}, {{\text{F}}_\omega }(2, 2r - 2))$ is the solution of Eq (4.9).

Proof. From Section 4.1, there are ${G}_{1}~\text{F}(2, 2r-2)$ and ${G_2} \sim {{\mathtt{χ}}^2}(2r)$. ${G_1}$ and ${G_2}$ are independent, so

The ${{\text{F}}_{1 - \omega }}(2, 2r - 2) < {G_1} < {{\text{F}}_\omega }(2, 2r - 2)$ may be written as

According to Lemma 3, ${{\text{F}}_{1 - \omega }}(2, 2r - 2)$ and ${{\text{F}}_\omega }(2, 2r - 2)$ are positive constants, the above inequation is equivalent to

${\mathtt{χ}}_\omega ^2(2r) < {G_2} < {\mathtt{χ}}_{1 - \omega }^2(2r)$ is equivalent to

5.   Monte Carlo simulation

In this Section, the proposed estimation methods are compared using MATLAB. For point estimates, the mean squared errors (MSEs) are calculated by Eq (5.1). For interval estimates, coverage probability (CP) is used to reflect the performance of GCIs. The simulation is carried out under true value $({\eta _{real}}, {\lambda _{real}}) = (2, 2)$ and different $n$, $r$ and $Q$, and the trials are $N$ at 1000 times. The hyper-parameter of the prior distribution is $({a_1}, {b_1}) = (2, 1.8)$ and $({a_2}, {b_2}) = (1.5, 2)$, and the parameters of SSE loss function and LINEX loss function are $d = 4$ and $a = 4$ respectively. The MSEs of $\eta$, and $\lambda$ are shown in Tables 1 and 2 respectively, and the MSEs of $R(t)$ are shown in Table 3 with $t = 3$. The CP values is shown in Table 4.

Balakrishnan and Sandhu ^[21] proposed an algorithm to produce progressive type-Ⅱ censored sample from any continuous distribution. The specific steps applied to the IWD are as follows:

Step 1. Generate samples ${w_1}, {w_2}, ..., {w_r}$ from $U(0, 1)$, where ${w_1}, {w_2}, ..., {w_r}$ are independent.

Step 2. Let ${v_i} = w_i^{{{(i + {Q_r} + {Q_{r - 1}} +... + {Q_{r - i + 1}})}^{ - 1}}}$ and ${u_i} = 1 - {v_r}{v_{r - 1}}...{v_{r - i + 1}}$.

Step 3. Let ${t_i} = {F^{ - 1}}({u_i}; {\eta _{real}}, {\lambda _{real}})$, where $F(\cdot)$ is the CDF (1.2) of IWD, $i = 1, 2, ..., r$.

Then, ${t_1} < {t_2} < ... < {t_r}$ are progressive type-Ⅱ censored data from ${\text{IW(}}\eta {\text{, }}\lambda {\text{)}}$ with a censoring scheme $Q$. Calculation results can be found in Tables 1–3.

From Table 1 to Table 3, we have the following conclusions:

(i) Obviously, considering the same $n$ and $r$, the censoring scheme $Q$ has a great influence on MSE.

(ii) Considering the same $n$, $r$ and $Q$, the BE of $\eta$ under SE loss function is the better than the MLE and IME.

(iii) Considering the same $n$, $r$ and $Q$, the BE $\lambda$ under SE loss function is better than others. However, the BE of $\lambda$ under SSE loss function is close to the BE under LINEX loss function.

(iv) Considering the same $n$, $r$ and $Q$, for the reliability $R(t)$, MSEs of MLEs and BEs under SE and SSE loss functions are relatively close, while, MSEs of BEs under LINEX loss function are larger than others.

From Table 4, We know that CP values for $\eta$ and $\lambda$ are all close to 95%.

6.   Real data analysis

There is a set of real data from Dumonceaux and Antle ^[22], which represents the maximum flood level (in millions of cubic feet per second) of the Susquehanna River at Harrisburg, Pennsylvania over 20 four-year periods (1890–1969) as:

0.265, 0.269, 0.297, 0.315, 0.324, 0.338, 0.379, 0.379, 0.392, 0.402, 0.412, 0.416, 0.418, 0.423, 0.449, 0.484, 0.494, 0.613, 0.654, 0.740.

According to the data, we can plot the empirical CDF and the CDF of the IWD, as shown in Figure1. In the IWD, we using the BEs under SE loss function as the value of parameter, i.e. $\eta = 0.0336, \lambda = 2.0431$. According to Figure 1, we can see that the IWD can well model the data. Therefore, we can conclude that this distribution is valid.

Now, the real data with censoring scheme $({Q_1}, {Q_2}, ..., {Q_{10}}) = (1, 1, ..., 1)$ are as follows:

0.265, 0.297, 0.324, 0.379, 0.392, 0.412, 0.418, 0.449, 0.494, 0.654.

Before proceeding with the estimation, it is necessary to establish the existence and uniqueness of the maximum likelihood estimate. However, proving this can be a complicated process due to the nonlinearity of the system of Eq (2.5). For this reason, we visualize it through Figure 2, where L1 represents $\frac{{\partial L(\eta, \lambda; t)}}{{\partial \eta }}$ in Eq (2.3) and L2 represents $\frac{{\partial L(\eta, \lambda; t)}}{{\partial \lambda }}$ in Eq (2.4).

From Figure 2, we know that the two curves intersect at only one point, indicating the presence of a unique MLE.

The estimates and generalized confidence intervals that obtained by using these proposed methods are shown in Table 5.

7.   Conclusions

In this paper, we have investigated the point estimation and interval estimation of parameters based on progressive type-Ⅱ censored sample for IWD. First, the Newton-Raphson iteration method is used to solve the likelihood equations of parameters and obtain their MLEs. Then, the BEs are derived based on SE and SSE loss functions, respectively. Finally, the IMEs are derived by generalized pivotal quantity. Additionally, the GCIs are also constructed by generalized pivotal quantity. Monte Carlo simulation is used to present the effect of the above estimators. The simulation results revealed that the estimators derived using Bayesian approach perform better than other methods in terms of MSE. Moreover, a real set of data is analyzed and the results coincide with simulation. Monte Carlo simulation results indicate that Bayesian estimation under the SE loss function works best among all the methods mentioned in this paper.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This research was funded by National Natural Science Foundation of China (No. 71661012).

Conflict of interest

The authors declare there is no conflict of interest.