Statistical Analysis for Competing Risks' Model with Two Dependent Failure Modes from Marshall-Olkin Bivariate Gompertz Distribution.

The bivariate or multivariate distribution can be used to account for the dependence structure between different failure modes. This paper considers two dependent competing failure modes from Gompertz distribution, and the dependence structure of these two failure modes is handled by the Marshall-Olkin bivariate distribution. We obtain the maximum likelihood estimates (MLEs) based on classical likelihood theory and the associated bootstrap confidence intervals (CIs). The posterior density function based on the conjugate prior and noninformative (Jeffreys and Reference) priors are studied; we obtain the Bayesian estimates in explicit forms and construct the associated highest posterior density (HPD) CIs. The performance of the proposed methods is assessed by numerical illustration.

1. Introduction

It is extremely common that the failure of a product or a system contains several competing failure modes in reliability engineering; any failure mode will lead to the failure result. Competing risks' data contain the failure time and the corresponding failure mode, which can be modeled by the competing risks' model and has been commonly performed in many research fields, such as engineering and medical statistics. Previous studies have mostly assumed the competing failure modes to be independent; Wang et al. [1], Ren and Gui [2], and Qin and Gui [3] focused on the independent competing risks' model under progressively hybrid censoring from Weibull and Burr-XII distributions. Objective Bayesian analysis for the competing risks' model with Wiener degradation phenomena and catastrophic failures was studied by Guan et al. [4]. In practice, the independency relationship between different failure modes is a very special case; a more common situation is dependency. That is, the failure mechanisms are interactive and interdependent; the occurrence of one failure mode will affect the occurrence of other failure modes. For example, a ship fixed carbon dioxide fire extinguishing system can fail due to pressure gauge, distribution value, cylinder group, and so on; these failure modes are dependent because they all are related to the storage environment. Therefore, it is more reasonable to assume dependency among different competing failure modes. The competing risks' model considers the product or system with multiple dependent competing failure modes, any one of which will cause the occurrence of failure. The dependent competing risks' model has been extensively studied. Zhang et al. [5] and Zhang et al. [6] studied the dependent competing risks' model under accelerated life testing (ALT) by copula function to measure the dependence between different competing failure modes; the results indicate the copula construction method has good accuracy and universality. Wang and Yan [7] and Wu et al. [8] also studied this model under ALT and progressively hybrid-censoring scheme using Clayton copula and Gumbel copula, respectively. For other related works, see the works of Lo and Wike [9] and Fang et al. [10].

In addition to using copula function to handle the relationship between different competing failure modes, the bivariate or multivariate distribution also can be used to account for the correlation between different failure modes. The Marshall–Olkin distribution [11], which has many good properties, is the best-known bivariate distribution and has been discussed extensively; it has a parameter to describe the dependence structure. Li et al. [12], Kundu and Gupta [13], and Bai et al. [14] provided reviews on Marshall–Olkin–Weibull distribution; Kundu and Gupta [13] obtained the explicit forms of the unknown parameters when the shape parameter is known; when the shape parameter is unknown, they used the importance sampling to compute the Bayesian estimates of the unknown parameters. Bai et al. [14] discussed the statistical analysis for the accelerated dependent competing risks' model under Type-II hybrid censoring schemes. Guan et al. [15] studied objective Bayesian analysis for the Marshall–Olkin exponential distribution based on reference priors; they also found that some of the reference priors are also matching priors and the posterior distributions based on these priors are proper.

The Gompertz distribution is a widely used growth model which has been studied extensively; Ismail [16] studied the Bayesian analysis of Gompertz distribution parameters and acceleration factor in the case of partially accelerated life testing under Type-I censoring scheme. Ghitany et al. [17] considered a progressively censored sample from Gompertz distribution; they discussed the existence and uniqueness of the MLEs of the unknown parameters. The Gompertz distribution plays an important role in fitting clinical trials' data in medical science and can be used to the theory of extreme-order statistics. In this paper, we will study the dependent competing risks' model from the Marshall–Olkin bivariate Gompertz (MOGP) distribution, which is a bivariate distribution with Gompertz marginal distributions. We focus our attention on the statistical analysis of the model parameters, including classical likelihood inference, Bayesian analysis, and objective Bayesian analysis. Because the Bayesian analysis based on conjugate prior is sensitive to the hyperparameters, inappropriate choice of it will cause bad priors. Based on this reason, we propose the objective Bayesian analysis based on noninformative priors for comparison. The objective Bayesian inference has been studied by Guan et al. [14], Bernardo [18], and Berger and Bernardo [19] based on Reference and Jeffreys priors.

In the rest of this paper, we will present the model description and some properties. Section 3 presents the MLEs and associated bootstrap CIs. In Section 4, Bayesian estimates and associated HPD CIs based on conjugate prior, Jeffreys prior [20], and reference priors [18] are obtained, and these priors lead to proper posteriors which are proved. Section 5 presents some results obtained from simulation study and illustrative analysis. Section 6 gives some final concluding remarks.

2. Model Description

Suppose that f(t; λ, θ) is a Gompertz distribution; the density function and reliability function of it arewhere λ is shape parameter and θ is scale parameter.

Suppose U0,  U1,  and  U2 are three independent Gompertz variables with different scale parameters, that is, U0 ~ GP(λ, θ0), U1 ~ GP(λ, θ1), and U2 ~ GP(λ, θ2). Let T1=min(U0, U1) and T2=min(U0, U2); we obtain T1 ~ GP(λ, θ0+θ1) and T2 ~ GP(λ, θ0+θ2). Then, the pair of variables (T1, T2) follows the MOGP distribution denoted by (T1, T2) ~ MOGP(λ, θ0, θ1, θ2). When θ0=0, the two variables T1 and T2 are independent and T1 and T2 will be dependent when θ0 > 0; hence, θ0 can be regarded as a correlation coefficient between T1 and T2.

The joint PDF of (T1,  T2) can be written as

The surface plots of f(t1, t2; λ, θ0, θ1, θ2) are presented in Figure 1. From Figure 1, we can see that f(t1, t2; λ, θ0, θ1, θ2) is a unimodal function.

Figure 1

Surface plot of f(t1, t2; λ, θ0, θ1, θ2) with different values of λ,  θ0,  θ1,  θ2. (a) (λ,  θ0,  θ1,  θ2)=(3,  0.5,  2,  1). (b) (λ,  θ0,  θ1,  θ2)=(3,  1.5,  0.5,  2). (c) (λ,  θ0,  θ1,  θ2)=(1,  0.5,  0.5,  0.5). (d) (λ,  θ0,  θ1,  θ2)=(1,  0.2,  0.8,  0.6).

Put n identical products into test, and each product has two dependent failure modes with lifetimes T1,  T2, (T1, T2) ~ MOGP(λ, θ0, θ1, θ2). Then, the system lifetime is X=min(T1, T2) ~ MOGP(λ, θ0+θ1+θ2). Let δ0=I(T1=T2), δ1=I(T1 < T2), and δ2=I(T1 > T2), for l=1, ⋯, n, where I(·) is an indicator function. Then, we can compute n0=∑δ0, n1=∑δ1,  n2=∑δ2, and n=n0+n1+n2.

Theorem 1 .

For l=1, ⋯, n, δ0=I(T1=T2), δ1=I(T1 < T2), and δ2=I(T1 > T2), We have

Proof

For l=1, ⋯, n, we have δ0+δ1+δ2=1,

Therefore, (δ0, δ1, δ2) ~ Multinomial(1; θ0/(θ0+θ1+θ2), θ1/(θ0+θ1+θ2), θ2/(θ0+θ1+θ2)).

The likelihood function iswhere

Then, we obtain

3. Classical Inference

3.1. Maximum Likelihood Estimates (MLEs)

The MLEs of θ0,  θ1,  θ2,  and λ can be obtained by maximizing the logarithm of L(x; λ, θ0, θ1, θ2). Set the first partial derivation of log  L(x; λ, θ0, θ1, θ2) about θ0,  θ1,  θ2,  λ to 0, i.e.,

From (8), we get the MLEs of θ0,  θ1,  and θ2 as

Substituting into log  L(x; λ, θ0, θ1, θ2), we obtainwhich is the profile logarithm likelihood function of λ.

We can show that ∂2h(λ)/∂λ2 < 0, which implies that h(λ) is concave. Some iterative schemes can be used to find the MLE for λ, such as Newton–Raphson algorithm.

3.2. Bootstrap Confidence Intervals (CIs)

Since it is hard to construct the exact CIs for the unknown parameters, we consider the Bootstrap method to construct CIs for parameters θ0,  θ1,  θ2,  and λ. The Bootstrap method is a resampling method to estimate some statistical characteristics for the unknown parameters by taking samples from the original samples repeatedly; the obtained samples are called Bootstrap samples. This method has a great practical value since it does not need to assume the overall distribution or construct the pivot quantity. We generate the Bootstrap sample by the following three steps:

Step 1: for the fixed value of n and observed data (x1, x2, ⋯, x), we get the estimates based on the maximum likelihood method.

Step 2: for the values of n, , we generate the sample (x1, x2, ⋯, x). Then, get the MLEs .

Step 3: repeat Step 2 M times to obtain M sets of the values . Arrange them as follows to get the Bootstrap sample:

Based on the Bootstrap sample and by percentile Bootstrap (Boot-P) method, we construct the Boot-P CIs for θ0,  θ1,  θ2,  λ at 1 − γ confidence level as

4. Bayesian Inference and HPD CIs

4.1. Conjugate Prior

In this section, we suppose the shape parameter λ is known. Denote θ=θ0+θ1+θ2, which has a Gamma prior with hyperparameters a and b as

Due to θ0/θ+θ1/θ+θ2/θ=1, so given θ, (θ1/θ,  θ2/θ) follows a Dirichlet prior with hyper parameters c0,  c1, and c2, that is,

Therefore, the joint prior of θ0, θ1,  and θ2 becomeswhere c=c0+c1+c2.

4.2. Jeffreys Prior

According to Jeffreys [20], Jeffreys prior is proportional to the square root of the determinant of the Fisher information matrix. From (7), we obtain the Fisher information matrix of (θ0, θ1,  θ2) as

From Theorem 1, we have n=n · θ/(θ0+θ1+θ2),  i=0,  1,  2, so I(θ0, θ1,  θ2) can be written as

Thus, the Jeffreys prior is given by

Theorem 2 .

Based on the Jeffreys prior π2(θ0, θ1,  θ2), the joint posterior distribution of (θ0, θ1,  θ2) is proper.

From (6) and (7), we obtain the joint posterior distribution of (θ0, θ1,  θ2) based on π2(θ0, θ1,  θ2) as

Integrating π2(θ0, θ1,  θ2|x) with respect to θ0, θ1,  and θ2, we obtainwhere A=∑(e − 1) and B(·,  ·) is a beta function.

Thus, the joint posterior distribution of (θ0, θ1,  θ2) based on π2(θ0, θ1,  θ2) is proper.

4.3. Reference Priors

Bernardo [18] and Berger and Bernardo [19] proposed the reference prior which plays a vital role in the objective Bayesian inference. We set μ0 ≡ θ=θ0+θ1+θ2, μ1=θ0/θ, and μ2=θ1/θ; the transformation from (θ0,  θ1,  θ2) to (μ0,  μ1,  μ2) is one-to-one with the inverse transformation θ0=μ0μ1, θ1=μ0μ2, and θ2=μ0(1 − μ1 − μ2). The Jacobian matrix of the transformation has the form

The likelihood function (3) becomes

The Fisher information matrix of (μ0,  μ1,  μ2) can be written as

Theorem 3 .

Under the ordering groups {μ0,  (μ1,  μ2)} and {(μ1,  μ2), μ0}, the reference priors are the same, which is given by ; the corresponding reference prior for (θ0, θ1,  θ2) is

Under the ordering groups {μ0,  μ1,  μ2}, {μ0,  μ2, μ1}, {μ1,  μ0,  μ2}, and {μ1,  μ2,  μ0}, the reference priors are the same, which is given by ; the corresponding reference prior for (θ0, θ1,  θ2) is

Under the ordering groups {μ2,  μ0,  μ1} and {μ2,  μ1, μ0}, the reference priors are the same, which is given by ; the corresponding reference prior for (θ0, θ1,  θ2) is

The Fisher information matrix of (μ0,  μ1,  μ2) is

where ∑11=n/μ02 and .

The reference prior for the ordering groups {μ0,  (μ1,  μ2)} and {(μ1,  μ2), μ0} is the same as in [21], which is given by

The inverse of I1 is

According the notations in [18], we obtain h1=1/μ02, h2=1/μ1(1 − μ1), and h3=(1 − μ1)/(μ2(1 − μ1 − μ2)).

Choose the compact sets Ω={(μ0, μ1, μ2)|a0 < μ0 < b0,  a1 < μ1,  a2 < μ2,  μ1+μ2 < d}, such that a0,  a1,  a2⟶0, b0⟶∞, and d⟶1, as k⟶∞. Then, we have

where :

Then, we get the reference prior as

where (μ0,  μ1,  μ2) is an inner point of Ω.

Similarly, under the ordering group {μ0,  μ2,  μ1}, the reference prior is ω(μ0,  μ1,  μ2).

The Fisher information matrix of {μ1,  μ0,  μ2} is

The inverse of I2 is

Similarly, we obtain h1=1/μ1(1 − μ1), h2=1/μ02, and h3=(1 − μ1)/(μ2(1 − μ1 − μ2)).

Choose the compact sets Ω={(μ1, μ0, μ2)|a0 < μ1,  a1 < μ0 < b1,  a2 < μ2,  μ1+μ2 < d}, such that a0,  a1,  a2⟶0, b1⟶∞, and d⟶1, as k⟶∞. Then, we have

where ,

Let (μ1,  μ0,  μ2) be an inner point of Ω; we get the reference prior as

Similarly, under the ordering group {μ1,  μ2,  μ0}, the reference prior is ω(μ0,  μ1,  μ2).

The Fisher information matrix of {μ2,  μ1,  μ0} is

The inverse of I3 is

Then, we obtain h1=1/μ2(1 − μ2), h2=(1 − μ2)/(μ1(1 − μ1 − μ2)), and h3=1/μ02.

Choose the compact sets Ω={(μ2, μ1, μ0)|a0 < μ2,  a1 < μ1,  μ2+μ1 < d, a2 < μ0 < b2}, such that a0,  a1,  a2⟶0, b2⟶∞, and d⟶1, as k⟶∞. Then, we havewhere ,

Let (μ2,  μ1,  μ0) be an inner point of Ω, we obtain the reference prior as

Similarly, under the ordering group {μ2,  μ0,  μ1}, the reference prior is ω(μ0,  μ1,  μ2). According to the one-to-one transformation from (μ0,  μ1,  μ2) to (θ0,  θ1,  θ2), we can obtain the reference priors π2(μ0,  μ1,  μ2), π3(μ0,  μ1,  μ2), π4(μ0,  μ1,  μ2) from ω, ω, and ω, respectively.

Theorem 4 .

Based on the reference priors π3(θ0, θ1,  θ2) and π4(θ0, θ1,  θ2), the posterior distributions of (θ0, θ1,  θ2) are proper.

The joint posterior distributions of (θ0, θ1,  θ2) based on reference prior π3(θ0, θ1,  θ2) and π4(θ0, θ1,  θ2) are, respectively, as

Integrating π3(θ0, θ1,  θ2|x) and π4(θ0, θ1,  θ2|x) with respect to θ0, θ1,  and θ2, respectively, we obtain

Thus, the posterior distributions of (θ0, θ1,  θ2) based on π3(θ0, θ1,  θ2) and π4(θ0, θ1,  θ2) are proper.

4.4. Bayesian Estimates

The joint posterior distributions of (θ0, θ1,  θ2) based on π1,  π2,  π3,  and π4 are, respectively, aswherewhere w1=Γ(∑2c)b/Γ(a) and w2=∏2b/Γ(c).

Thus, we obtain

Similarly,where

We obtainwhere

We obtainwhere

Then, we have

From (9)–(12), we get the Bayesian estimates of parameters θ0,  θ1,  θ2,  and θ against squared error loss function based on π1,  π2,  π3,  and π4, respectively, which are listed in Table 1.

Table 1

Bayesian estimates of parameters based on different priors.

Prior	θ 0	θ 1	θ 2	θ
π 1	λ(n0+c0)(n+a)/(n+c0+c1+c2)(A+bλ)	λ(n1+c1)(n+a)/(n+c0+c1+c2)(A+bλ)	λ(n2+c2)(n+a)/(n+c0+c1+c2)(bλ+A)	(n+a)λ/A+bλ
π 2	nλ(2n0+1)/A(2n+3)	nλ(2n1+1)/A(2n+3)	nλ(2n2+1)/A(2n+3)	nλ/A
π 3	nλ(2n0+1)/2A(n+1)	nλ(2n1+2n2+1)(2n1+1)/4A(n+1)(n1+n2+1)	nλ(2n1+2n2+1)(2n2+1)/4A(n+1)(n1+n2+1)	nλ/A
π 4	nλ(2n0+2n2+1)(2n0+1)/4A(n+1)(n0+n2+1)	nλ(2n1+1)/2A(n+1)	nλ(2n0+2n2+1)(2n2+1)/2A(n+1)(n0+n2+1)	nλ/A

4.5. HPD Credible Intervals

The HPD credible intervals of parameters θ0,  θ1,  θ2,  and θ can be constructed by the Monte Carlo method studied by Chen and Shao [22].

Step 1: given the value of n and the observed data (x1, x2, ⋯, x), compute the Bayesian estimates of based on π1,  π2,  π3,  and π4, respectively.

Step 2: repeat Step 1 M times; we obtain M sets of the values based on π1,  π2,  π3,  and π4, respectively. Arrange them in the ascending order, we obtain

Step 3: compute the CIs at 1 − γ confidence level as

: the HPD CIs for θ, v=0,1,2, and θ are the shortest intervals among and w=1,  2, ⋯, M − (1 − γ)M, respectively.

5. Numerical Simulation and Illustrative Example

5.1. Simulation

Suppose the common shape parameter λ is known. The initial values for parameters (λ,  θ0, θ1,  θ2) are (3,  1,  2,  1). The initial values for the hyperparameters a,  b,  c0,  c1,  and c2 are all 0.001. Take the sample size n = 10, 20, 30, and 50. Generate the random samples (x1, x2, ⋯, x) from MOGP(λ, θ0, θ1, θ2) by the following steps:

Step 1: for a fixed value n, generate n samples u01, u02, ⋯, u0 from GP(λ, θ0), u11, u12, ⋯, u1 from GP(λ, θ1), and u21, u22, ⋯, u2 from GP(λ, θ2). Then, we obtain T1=min(u0, u1) and T2=min(u0, u2),  l=1,2, ⋯, n.

Step 2: compute (x, δ0, δ1, δ2),  l=1,2, ⋯, n, where x=min(T1, T2), δ0=I(T1=T2), δ1=I(T1 < T2), and δ2=I(T1 > T2).

Repeat the procedures 10,000 times; we get the values of the mean squared errors (MSEs) of the MLEs, the average lengths (ALs), and coverage probabilities (CPs) of the 95% Boot-P CIs, and the MSEs of the Bayesian estimates, the ALs, and CPs of the 95% HPD CIs, which are shown in Table 2–5. From the results in Table 2–5, we can make the following conclusions.

Table 2

MSEs, ALs, and CPs of θ0,  θ1,  θ2,  and θ (n = 10).

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MSE	0.4858	0.8030	0.4865	0.9374
		Boot-AL	2.2414	2.7146	2.2266	1.7920
		Boot-CP	0.9339	0.9294	0.9405	0.9321
Bayes	π 1	MSE	0.4850	0.8012	0.4857	0.9340
		HPD-AL	2.0388	2.5119	2.0425	1.9018
		HPD-CP	0.9663	0.9440	0.9645	0.9369
	π 2	MSE	0.4055	0.5903	0.4061	0.9374
		HPD-AL	1.7980	2.2183	1.8016	1.9034
		HPD-CP	0.9552	0.9399	0.9539	0.9335
	π 3	MSE	0.4678	0.5732	0.3909	0.9374
		HPD-AL	1.8797	2.2193	1.7850	1.9034
		HPD-CP	0.9481	0.9405	0.9569	0.9460
	π 4	MSE	0.3748	0.7042	0.3754	0.9374
		HPD-AL	1.7724	2.3192	1.7760	1.9034
		HPD-CP	0.9527	0.9468	0.9515	0.9405

Table 3

MSEs, ALs, and CPs of θ0,  θ1,  θ2,  and θ (n = 20).

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MSE	0.2505	0.4519	0.2523	0.6907
		Boot-AL	1.5795	1.9048	1.5807	1.2957
		Boot-CP	0.9488	0.9483	0.9412	0.9407
Bayes	π 1	MSE	0.2503	0.4512	0.2520	0.6893
		HPD-AL	1.4434	1.7573	1.4382	1.3635
		HPD-CP	0.9832	0.9692	0.9831	0.9415
	π 2	MSE	0.2335	0.3766	0.2350	0.6907
		HPD-AL	1.3512	1.6476	1.3462	1.3640
		HPD-CP	0.9746	0.9447	0.9762	0.9409
	π 3	MSE	0.2551	0.3662	0.2293	0.6907
		HPD-AL	1.3834	1.6486	1.3398	1.3640
		HPD-CP	0.9668	0.9506	0.9777	0.9598
	π 4	MSE	0.2201	0.4260	0.2216	0.6907
		HPD-AL	1.3399	1.6868	1.3348	1.3640
		HPD-CP	0.9614	0.9525	0.9613	0.9498

Table 4

MSEs, ALs, and CPs of θ0,  θ1,  θ2,  and θ (n = 30).

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MSE	0.1752	0.3345	0.1771	0.6049
		Boot-AL	1.2849	1.5451	1.2896	1.0510
		Boot-CP	0.9651	0.9516	0.9654	0.9415
Bayes	π 1	MSE	0.1751	0.3341	0.1770	0.6040
		HPD-AL	1.1710	1.4354	1.1727	1.1164
		HPD-CP	0.9919	0.9629	0.9901	0.9427
	π 2	MSE	0.1694	0.2922	0.1712	0.6049
		HPD-AL	1.1197	1.3745	1.1213	1.1167
		HPD-CP	0.9839	0.9723	0.9835	0.9418
	π 3	MSE	0.1814	0.2849	0.1679	0.6049
		HPD-AL	1.1377	1.3750	1.1177	1.1167
		HPD-CP	0.9783	0.9770	0.9851	0.9615
	π 4	MSE	0.1612	0.3228	0.1629	0.6049
		HPD-AL	1.1132	1.3966	1.1148	1.1167
		HPD-CP	0.9896	0.9581	0.9885	0.9638

Table 5

MSEs, ALs, and CPs of θ0,  θ1,  θ2,  and θ (n = 50).

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MSE	0.1158	0.2460	0.1161	0.5380
		Boot-AL	0.9947	1.1981	1.0018	0.8227
		Boot-CP	0.9829	0.9578	0.9831	0.9554
Bayes	π 1	MSE	0.1157	0.2458	0.1161	0.5375
		HPD-AL	0.9118	1.1075	0.9071	0.8677
		HPD-CP	0.9955	0.9822	0.9954	0.9724
	π 2	MSE	0.1150	0.2243	0.1154	0.5380
		HPD-AL	0.8874	1.0786	0.8828	0.8679
		HPD-CP	0.9884	0.9900	0.9870	0.9702
	π 3	MSE	0.1209	0.2196	0.1137	0.5380
		HPD-AL	0.8961	1.0791	0.8811	0.8679
		HPD-CP	0.9813	0.9933	0.9901	0.9721
	π 4	MSE	0.1105	0.2417	0.1109	0.5380
		HPD-AL	0.8842	1.0892	0.8796	0.8679
		HPD-CP	0.9938	0.9789	0.9931	0.9717

The MSEs of MLEs and Bayesian estimates decrease as the sample size increases. For given sample size n, the Bayesian estimates based on π1, π2,  and π4 are smaller than the MSEs of MLEs. The MSEs of Bayesian estimates of θ0 and θ2 based on π4 are smaller than that based on π1, π2,  and π3. The MSEs of Bayesian estimates of θ1 based on π3 are smaller than that based on π1, π2,  and π4. The MSEs of Bayesian estimates of θ based on π1 are smaller than that based on π2, π3,  and π4.

The CPs of Boot-P and HPD CIs are all close to 0.95. The ALs of Boot-P and HPD CIs decrease; the associated CPs increase when the sample size increases. The CPs of HPD CIs based on Bayesian estimates are larger than the CPs of Boot-P CIs based on MLEs.

5.2. Illustrative Analysis

5.2.1. Simulated Data

For illustrative purposes, with initial value for parameters (λ, θ0, θ1, θ2) as (3,1,2,1), we use the procedures mentioned above to generate U0,  U1,  and U2 from GP(3,1), GP(3,2), and GP(3,1), respectively. We then get T1=min(U0, U1) and T2=min(U0, U2); the latent lifetime of the system is min(T1, T2). The simulated data are listed in Table 6. The MLEs, Bayesian estimates, and associated 95% CIs for parameters θ0,  θ1,  θ2,  and θ are shown in Table 7. From Table 7, all the MLEs and Bayesian estimates of (θ0,  θ1,  θ2, θ) are close to the true value.

Table 6

The simulated data when n = 25.

(0.0019

1)

(0.0025

1)

(0.0062

1)

(0.0361

0)

(0.0651

2)

(0.0675

1)

(0.1108

2)

(0.1447

1)

(0.1509

1)

(0.1694

2)

(0.1737

0)

(0.1859

1)

(0.1900

2)

(0.2218

0)

(0.2307

2)

(0.2537

1)

(0.2558

2)

(0.2734

2)

(0.2750

0)

(0.3349

1)

(0.3528

1)

(0.3790

0)

(0.3824

2)

(0.3875

1)

(0.5336

1)

Table 7

Point estimates and 95% CIs of θ0,  θ1,  θ2,  and θ.

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MLE	0.8029	1.9270	1.2847	4.0146
		Boot-CI	(0.0408, 1.8099)	(0.1891, 2.8803)	(0.0642, 1.8178)	(0.7671, 4.4112)
Bayes	π 1	Bayes	0.8029	1.9267	1.2845	4.0141
		HPD CI	(0.0898, 1.6991)	(0.3473, 2.9083)	(0.0396, 1.7374)	(0.8141, 4.1876)
	π 2	Bayes	0.8332	1.8937	1.2877	4.0146
		HPD CI	(0.0694, 1.7457)	(0.3070, 2.8296)	(0.0364, 1.7697)	(0.7893, 4.4378)
	π 3	Bayes	0.8492	1.8841	1.2812	4.0146
		HPD CI	(0.0798, 1.7561)	(0.1251, 2.5045)	(0.0315, 1.4392)	(0.6368, 4.1407)
	π 4	Bayes	0.8189	1.9301	1.2656	4.0146
		HPD-CP	(0.0638, 1.4011)	(0.2448, 2.8554)	(0.0456, 1.7418)	(0.8474, 4.3484)

5.2.2. Real Data

Use the procedures mentioned above to a real dataset. Kundu and Gupta [13] analyzed the football data of UEFA Champions' League data which are presented in Table 1. From the data, T1 and T2 can be regarded as two dependent failure modes, and n0=7, n1=17,  and n2=13. This data have been fitted by Marshall–Olkin bivariate Gompertz distribution (see Wang et al. [23]).

The MLEs, Bayesian estimates, and associated 95% CIs for parameters θ0,  θ1,  θ2,  and θ are shown in Table 8. From Tables 7 and 8, Bayesian estimates under different priors are close to MLEs, and the lengths of 95% Boot-p CIs associated to MLEs are longer than the lengths of 95% HPD CIs associated to Bayesian estimates.

Table 8

Point estimates and 95% CIs of θ0,  θ1,  θ2,  and  θ.

Method		Para.	θ 0	θ 1	θ 2	θ
MLE		MLE	1.0882e-2	0.4664e-2	1.3214e-2	2.8760e-2
		Boot-CI	(0.7354e-3, 1.2770e-2)	(0.3244e-2, 2.4914e-2)	(0.4077e-3, 1.5990e-2)	(0.7324e-2, 3.5454e-2)
Bayes	π 1	Bayes	1.1882e-2	0.6879e-2	1.3758e-2	3.2520e-2
		HPD CI	(0.2305e-2, 1.1210e-2)	(0.3998e-2, 1.9825e-2)	(0.1861e-2, 1.4366e-2)	(1.0088e-2, 3.8182e-2)
	π 2	Bayes	1.0832e-2	0.4856e-2	1.3073e-2	2.8760e-2
		HPD CI	(0.8317e-3, 1.1676e-2)	(0.2315e-2, 1.8702e-2)	(0.8807e-3, 1.3700e-2)	(0.5246e-2, 3.2002e-2)
	π 3	Bayes	1.0974e-2	0.4817e-2	1.2969e-2	2.8760e-2
		HPD CI	(0.8499e-3, 1.1535e-2)	(0.3160e-2, 1.7152e-2)	(0.8814e-3, 1.4144e-2)	(0.6994e-2, 3.2140e-2)
	π 4	Bayes	1.0803e-2	0.4919e-2	1.3038e-2	2.8760e-2
		HPD-CP	(0.8949e-3, 1.1372e-2)	(0.3365e-2, 1.8866e-2)	(0.6968e-3, 1.4164e-2)	(0.7224e-2, 3.0467e-2)

6. Conclusion

This paper discussed the point estimates and CIs for the parameters of the dependent competing risks' model from MOGP distribution. We studied the appropriateness of the posteriors based on conjugate prior and Jeffreys and Reference priors, obtained the Bayesian estimates in closed forms, and constructed the associated HPD CIs. From the simulations results, the use of the Bayesian method can be recommended if the priors are available. The results of the illustrative analysis show that the proposed methods work well; from the lengths of CIs, we can conclude the Bayesian estimates are better than MLEs in general.