Statistical analysis of dependent competing risks model from Gompertz distribution under progressively hybrid censoring.

Previous studies have mostly considered the competing risks to be independent even when the interpretation of the failure modes implies dependency. This paper studies the dependent competing risks model from Gompertz distribution under Type-I progressively hybrid censoring scheme. We derive the maximum likelihood estimations of the model parameters, and then the asymptotic likelihood theory and Bootstrap method are used to obtain the confidence intervals. The simulation results are provided to investigate the effects of different dependence structures on the estimations of parameters. Finally, one data set was used for illustrative purpose.

Background

The competing risks model involves multiple failure modes when only the smallest failure time and the associated failure mode are observed. This model is widely studied in the medical, actuarial, biostatistics and so on, under the assumption of independent competing risks. It is common that a failure is associated with one of the several competing failure modes. Previous studies have mostly considered the competing failure modes to be independent even when the interpretation of the failure modes implies dependency. Such as, in the study of colon cancer, the failure causes were cancer recurrence or death, obviously, such failure causes were dependent [see Lin et al. (1999)]. The competing risks model assuming independence among competing failure modes has been widely studied [see, e.g., Crowder (2001)]. Kundu et al. (2004) analyzed the progressively censored competing risks data, Sarhan (2007) analyzed the competing risks models with generalized exponential distributions, Cramer and Schmiedt (2011) studied the progressively censored competing risks data with Lomax distribution, other related works see, Bunea and Mazzuchi (2006); Balakrishnan and Han (2008); Pareek et al. (2009); Xu and Tang (2011), and so on.

The competing risks model under the assumption of dependent competing failure modes has been considered in the early work by Elandt-Johnson (1976). Afterwards, a number of corresponding works have been devoted to the dependent competing risks model. Zheng and Klein (1995) considered the dependence structure between failure modes is represented by an assumed Archimedean copula. Other works see Escarela and Carriere (2003); Kaishev et al. (2007).

In this paper, we present a dependent competing risks model from Gompertz distribution under Type-I progressively hybrid censoring scheme (PHCS). The Gompertz distribution is one of classical mathematical models and was first introduced by Gompertz (1825), which is a commonly used growth model in actuarial and reliability and life testing, and plays an important role in modeling human mortality and fitting actuarial tables and tumor growth. This distribution has been widely used, see, Ali (2010); Ghitany et al. (2014).

The Type-I PHCS was first proposed by Kundu and Joarder (2006) [see also Childs et al. (2008)]. This censoring scheme has been widely used in reliability analysis, see, Chien et al. (2011); Cramer and Balakrishnan (2013). It can be defined as follows: suppose n identical units are put to life test with progressive censoring scheme , the experiment is terminated at time , where and m are fixed in advance. At the time of the first failure of the remaining units are randomly removed, at the time of the second failure of the remaining units are randomly removed and so on. If the mth failure time occurs before time , all the remaining units are removed and the terminal time of the experiment is . On the other hand, if the mth failure time does not occur before time and only J failures occur before time , where . Then all the remaining units are removed and the terminal time of the experiment is . We denote the two cases as

Case I

Case II

The rest of the paper is organized as follows. “Model description” section provides the model description, “Maximum likelihood estimations (MLEs)” section presents the maximum likelihood estimations of the model parameters. The confidence intervals are provided in “Confidence intervals” section. “Simulation and data analysis” section presents the simulation and data analysis. Conclusion appears in “Conclusion” section.

Model description

It is assumed that the Gompertz distribution with shape parameter λ and scale parameter θ has the following probability density function (PDF), cumulative distribution function (CDF) and survival function respectively, where . We denote the Gompertz distribution by .

Suppose variables are independent and follows , , . Define , , then the distributions of are and , respectively.

Theorem 1

The joint survival function of is

Proof

Corollary 1

The joint PDF of can be written as

Proof For the cases and , can be easily obtained by . For the case , by the full probability formula, we have the fact thatwhere So from (4), we have

So we have . □

Figure 1 presents the surface plot of for different values of , from Fig. 1, we can see that the is unimodal. Define , and the distribution of is . indicates that are independent. Therefore, can be regarded as the dependence structure between .

Fig. 1

The surface plot of the joint PDF of with different values of when . a , b , c , d

Competing risks model

Consider two competing failure modes with latent lifetimes in the experiment under Type-I PHCS, the failure of an individual is caused by any single one of the two failure modes, obviously, the actual lifetime span is . Let r denotes the number of failures that occur before time τ, τ* denotes the terminal time. Then, at time all the remaining units are removed and the experiment is terminated, where , , rm = 0 in Case I and r = J, in Case II.

For the competing risks model under Type-I PHCS, are the observed failure data, where are order statistics, takes any integer in the set . For , . denotes the number of failures caused by the two competing failure modes, denotes the number of failures caused by competing failure mode , where .

Maximum likelihood estimations (MLEs)

The likelihood function for the two competing risks model under Type-I PHCS can be written aswhere So the likelihood function can be written as By setting the first partial derivative of about to zero, we get From (7), (8) and (9), the estimates of are given by

Substituting into and ignoring the constant, we obtain the profile log-likelihood function of λ as

Lemma 1

The profile log-likelihood function is concave.

Proof Denote , where . Therefore, we get where , .

by the Cauchy–Schwarz inequality, therefore , which implies that the second derivative of is negative, so is concave. □

From Lemma 1, we know that is unimodal and it has a unique maximum. Since is unimodal, most of the standard iterative procedure can be used to find the MLE. So we propose to use the following simple algorithm. Substituting into (10), the MLE of satisfies the following equation,where

Using the method of a simple iterative scheme proposed in the literature by Kundu (2007), we can solve the shape parameter from (13). Start with an initial guess of , say , then obtain and proceed in this way to obtain . Stop the iterative procedure when , some pre-assigned tolerance limit. Once we obtain , the MLEs of can be obtained from (11) as .

Confidence intervals

Observed fisher information

In this section, we will construct the asymptotic confidence intervals (ACIs) for the parameters using the asymptotic likelihood theory. The observed Fisher information matrix is denoted by

,

where the elements of which are negative second partial derivatives of . Denote as the approximate asymptotic variance–covariance matrix of the MLEs of and as the estimation of , we get

By the asymptotic distribution of MLEs, follows as approximately standard normal distribution. Therefore, the two-sided ACIs for are given by

where is the quantile of a standard normal distribution.

Bootstrap sample

Step1. Given and progressive censoring scheme , compute the MLEs based on the original Type-I progressively hybrid censored sample .

Step2. Based on , , generate a Type-I progressively hybrid censored sample .

a1. Generate a random sample from Uniform distribution , where are order statistics. Let , are order statistics followed Uniform distribution .

a2. We obtain the failures before time and the terminal time .

If , , ;

If , , , where J is obtained from the inequality

, for , we set .

Step3. Based on and , we obtain the MLEs .

Step4. Repeat steps 2–3 N times, we obtain N estimates . Arrange them in ascending order to obtain the bootstrap sample .

The two-sided percentile bootstrap confidence intervals (Boot-P CIs) for parameters

Simulation and data analysis

Simulation

In this section, we presented some simulation results to evaluate the performance of all the methods proposed in the previous sections for different sample size n, different effective sample size m and different dependence structure .

Consider two competing failure modes, the initial values for parameters are . Take the dependence structure , where indicates that the two competing failure modes are independent Generate the Type-I PHC samples from the Gompertz distribution for competing failure mode according to the algorithm proposed by Balakrishnan and Sandhu (1995). Take the terminal time , and n = 20, 30, 50, m = 4, 6, 8, 10, 15, the pre-fixed scheme are

To compute the MLEs of , we have used the iterative procedure described in “Maximum likelihood estimations (MLEs)” section and stopped the iterative procedure when the difference between two consecutive iterates is less than . Before going to compute the MLEs, we plot the profile log-likelihood function of λ in Fig. 2. Figure 2 shows that the profile log-likelihood function of λ is unimodal, the MLE of λ is close to 0.6, so we start the iteration with the initial guess that .

Fig. 2

Profile log-likelihood function of λ

Repeat 10,000 times for each given n, m, and censoring scheme, the average mean squared errors (MSEs) and the average absolute relative bias (RABias) and the coverage percentage of the ACIs and Boot-P CIs are shown in Tables 1, 2 and 3.

Table 1

n = 20,

m	θ 0	θ 0		θ 1		θ 2		λ
		MSEs	ACI	MSEs	ACI	MSEs	ACI	MSEs	ACI
		RABias	Boot-P	RABias	Boot-P	RABias	Boot-P	RABias	Boot-P
4	0	0.221	0.963	0.8171	0.900	0.7271	0.912	0.2516	0.904
			0.908	0.7144	0.869	0.837	0.836	0.7974	0.902
	0.3	0.1299	0.958	0.7194	0.914	0.6848	0.857	0.2402	0.980
		0.7855	0.934	0.6744	0.877	0.7961	0.852	0.7832	0.934
	0.8	0.4077	0.969	0.626	0.929	0.6395	0.940	0.2315	0.965
		0.7278	0.886	0.6189	0.907	0.7621	0.917	0.7672	0.978
	1.2	0.9094	0.974	0.6771	0.927	0.629	0.931	0.2245	0.963
		0.7431	0.896	0.5798	0.912	0.7307	0.923	0.7516	0.986
	1.6	1.7076	0.961	0.7727	0.918	0.639	0.915	0.2184	0.947
		0.7724	0.866	0.5649	0.905	0.7023	0.891	0.7409	0.991
8	0	0.1081	0.925	0.4369	0.899	0.5514	0.903	0.217	0.941
			0.898	0.5018	0.879	0.704	0.853	0.7241	0.897
	0.3	0.1057	0.957	0.3504	0.915	0.4752	0.947	0.215	0.949
		0.7145	0.925	0.4321	0.891	0.6533	0.869	0.7125	0.924
	0.8	0.2489	0.952	0.3168	0.957	0.4051	0.974	0.2038	0.987
		0.5349	0.967	0.3832	0.934	0.5686	0.967	0.6984	0.954
	1.2	0.5723	0.917	0.3634	0.905	0.365	0.943	0.1952	0.942
		0.563	0.968	0.3784	0.951	0.5278	0.977	0.6822	0.952
	1.6	1.1113	0.906	0.5161	0.899	0.3765	0.902	0.1962	0.937
		0.6053	0.879	0.4238	0.893	0.5122	0.904	0.6882	0.936

Table 2

n = 30,

m	θ 0	θ 0		θ 1		θ 2		λ
		MSEs	ACI	MSEs	ACI	MSEs	ACI	MSEs	ACI
		RABias	Boot-P	RABias	Boot-P	RABias	Boot-P	RABias	Boot-P
6	0	0.1572	0.907	0.6745	0.907	0.6719	0.898	0.2554	0.933
			0.903	0.6554	0.911	0.7931	0.895	0.7989	0.892
	0.3	0.0848	0.958	0.5776	0.898	0.6121	0.914	0.2484	0.929
		0.6838	0.921	0.5969	0.899	0.753	0.875	0.7953	0.928
	0.8	0.3152	0.971	0.4578	0.957	0.5416	0.971	0.239	0.968
		0.6385	0.868	0.517	0.897	0.6897	0.901	0.7812	0.943
	1.2	0.7678	0.980	0.4317	0.968	0.5079	0.934	0.234	0.988
		0.6835	0.937	0.4716	0.913	0.6526	0.915	0.7744	0.983
	1.6	1.5004	0.929	0.4752	0.918	0.4819	0.927	0.2291	0.951
		0.728	0.898	0.4655	0.826	0.6298	0.877	0.7648	0.889
10	0	0.1101	0.917	0.4546	0.913	0.5544	0.899	0.2253	0.914
			0.914	0.5234	0.895	0.7262	0.879	0.7441	0.927
	0.3	0.0733	0.929	0.3553	0.924	0.485	0.908	0.2227	0.931
		0.6668	0.920	0.4471	0.894	0.6565	0.897	0.7432	0.930
	0.8	0.2195	0.972	0.2674	0.961	0.3879	0.962	0.2173	0.947
		0.5113	0.939	0.3632	0.869	0.5776	0.929	0.7338	0.946
	1.2	0.5649	0.968	0.2733	0.977	0.3314	0.984	0.2108	0.967
		0.5739	0.965	0.3404	0.915	0.513	0.978	0.7224	0.938
	1.6	1.135	0.943	0.3163	0.929	0.3108	0.953	0.2088	0.905
		0.6222	0.865	0.3525	0.886	0.4856	0.912	0.72	0.894

Table 3

n = 50,

m	θ 0	θ 0		θ 1		θ 2		λ
		MSEs	ACI	MSEs	ACI	MSEs	ACI	MSEs	ACI
		RABias	Boot-P	RABias	Boot-P	RABias	Boot-P	RABias	Boot-P
10	0	0.1443	0.916	0.5431	0.931	0.6052	0.909	0.25	0.933
			0.897	0.588	0.865	0.7577	0.842	0.794	0.934
	0.3	0.0577	0.914	0.4382	0.948	0.535	0.914	0.2467	0.941
		0.5855	0.897	0.5133	0.878	0.7088	0.845	0.7933	0.928
	0.8	0.2409	0.929	0.311	0.961	0.4439	0.968	0.2375	0.948
		0.5498	0.935	0.4103	0.920	0.627	0.896	0.7867	0.956
	1.2	0.6486	0.968	0.2746	0.967	0.3905	0.955	0.2349	0.967
		0.6286	0.963	0.3673	0.936	0.5754	0.911	0.7782	0.967
	1.6	1.2801	0.941	0.2878	0.949	0.3568	0.947	0.2311	0.958
		0.672	0.855	0.3549	0.864	0.5369	0.866	0.7719	0.936
15	0	0.1052	0.928	0.396	0.937	0.5194	0.934	0.2193	0.928
			0.892	0.4954	0.899	0.6969	0.802	0.7281	0.927
	0.3	0.051	0.933	0.2819	0.946	0.4332	0.929	0.2161	0.941
		0.5615	0.921	0.4028	0.878	0.6385	0.863	0.7245	0.936
	0.8	0.1705	0.967	0.1874	0.968	0.3368	0.964	0.2095	0.973
		0.4524	0.959	0.298	0.936	0.5408	0.938	0.722	0.957
	1.2	0.4865	0.972	0.1725	0.971	0.282	0.978	0.2068	0.968
		0.5337	0.951	0.2746	0.941	0.4781	0.942	0.7157	0.949
	1.6	1.0174	0.928	0.2367	0.944	0.2386	0.929	0.2048	0.927
		0.5948	0.836	0.302	0.855	0.424	0.914	0.7125	0.934

From Tables 1, 2 and 3, the observations can be made. For fixed sampling scheme, sample size n and dependence structure , the MSEs and RABias decrease as the effective sample size m increase.

For fixed sampling scheme, sample size n and effective sample size m, as the dependence structure of competing failure modes become stronger, the MSEs and RABias get smaller, while the MSEs and RABias with are bigger, which shows that the performance of the MLEs depends on the strength of dependence. This also shows that the dependence structure is very important in the competing risks model.

For fixed sampling scheme, n, m and dependence structure , the ACIs are stable than the Boot-P CIs, they can maintain their coverage percentages at the pre-fixed normal level.

Data analysis

Using the procedures above, we generate the Type-I PHC samples when with initial value for parameters as , and the dependence structure , the censoring scheme as . The simulated data is listed in Table 4. The MLEs and 95 % ACIs and Boot-P CIs are shown in Table 5. The trace plot of the MLE for parameter using the iterative procedure is shown in Fig. 3, which shows that the estimate of converges to a value after about 1000 iterations.

Table 4

The simulated data

i	1	2	3	4	5	6	7	8	9	10
t_i	0.0035	0.0181	0.0435	0.0813	0.0860	0.1286	0.1483	0.1484	0.1929	0.4449
α_i	2	2	0	0	2	1	0	1	1	2

Table 5

MLEs and 95 % CIs of the parameters

Para.	True value	MLE	ACI	Boot-P CI
θ 0	0.8	0.8934	(0.3777, 2.1645)	(0.2811, 0.9764)
θ 1	1.2	0.6627	(0.1987, 1.5241)	(0.1728, 1.3569)
θ 2	1	0.8136	(0.1167, 1.7438)	(0.2718, 1.1114)
λ	0.6	0.6935	(0.1911, 2.9921)	(0.4962, 0.7265)

Fig. 3

Trace plot of MLE for λ using the iterative procedure

Conclusion

This paper proposed the dependent competing risks model from Gompertz distribution under Type-I PHCS. We obtained the MLEs and ACIs and Boot-P CIs for the parameters. Simulations showed that the ACIs are more stable than the Boot-P CIs and that the dependence structure is important in the competing risks model. For a given sample size, the performance of the MLEs declined with increasing dependence, which suggests that greater dependence will require a larger sample size to achieve a particular level of precision in estimation.