广义指数分布下区间删失数据贝叶斯回归分析

李 群, 董小刚, 王纯杰, 赵 波

(长春工业大学 基础科学学院， 吉林 长春 130012)

广义指数分布下区间删失数据贝叶斯回归分析

李 群, 董小刚, 王纯杰*, 赵 波

(长春工业大学 基础科学学院， 吉林 长春 130012)

研究了在两参数广义指数分布下的区间删失寿命时间的贝叶斯回归分析模型。生存时间在服从广义指数分布的情况下，假定形状参数的先验分布来自伽马分布，建立了尺度参数与生存时间贝叶斯回归模型，从而得到生存时间的变化。选取MCMC算法对参数进行估计，并运用R软件进行了模拟。

广义指数分布; 区间删失; 贝叶斯回归; MCMC算法

0 引 言

在可靠性寿命试验中，两参数广义指数分布可简称广义指数分布或GE分布。作为指数分布的推广，由于广义指数分布对于删失时间数据有很好的分析效果，而且还可以作Gamma分布和Weibull分布的替代分布，因而在寿命试验和可靠性工程中有着重要的应用[1-2]。寿命数据分析已经成为航空、工程、医学和生物科学等多个领域中统计学家和实际工作者十分关心的一个问题，因此，对广义指数分布的研究有着十分重要的实际意义。同时，在生存分析中也经常研究感兴趣的时间与哪些因素有密切的关系，也会研究不同的药物类型中，哪种药物对于患者更有效果等。文中将通过建立贝叶斯回归模型来进行研究感兴趣的时间与相关因素的关系及影响[3]。

李荣[4]于2006年给出了一篇删失实验寿命的贝叶斯威布尔生存回归模型，建立了威布尔分布下关于参数λ的回归模型，并给回归系数赋予先验分布。在删失寿命实验的条件下，给出了贝叶斯威布尔回归模型的似然函数，基于Gibbs抽样得出参数的后验分布，利用WinBUGS软件包求解威布尔回归模型的贝叶斯估计的过程。朱惠明[5]等于2007年给出了删失试验寿命的贝叶斯生存极值回归模型，同样引入参数λ的协变量,并建立了贝叶斯回归模型，用MCMC方法和Gibbs抽样获得参数后验分布，同样利用WinBUGS软件包求解极值回归模型的贝叶斯估计的过程。Upadhyay[6]发表了基于Gibbs抽样下对数正态回归的后验分析，分别建立对数正态分布的均值、方差两个参数关于协变量影响的贝叶斯回归模型。Puja Makkar[7]给出了头颈癌在对数正态模型下的贝叶斯生存分析，在不知道先验信息的情况下，采用Gibbs抽样的方法得到参数的后验分布，并分析不同的治疗方案对患头颈癌患者寿命的影响。

广义指数分布是由Gupta R D和Kundu D于1999年提出的。此外Gupta R D[8-9]等给出了广义指数分布的一些统计推断的性质。Kundu D[10]等于2008年给出了广义指数的贝叶斯估计的相关理论。此外，郭环[11]研究了两参数广义指数分布的一些参数估计方法和优良性质，给出了在一定条件下两个参数的贝叶斯估计。但是上述文献均没有涉及广义指数分布的贝叶斯生存回归模型。

MCMC方法是一种简单易行、广泛应用的计算随机模拟方法。该方法的核心思想是构建一个概率转移矩阵，建立一个以分布π(x)为平稳分布的Markov链，得到π(x)的样本，通过随机抽样得到的样本就可以进行各种统计推断和估计[12]。MCMC方法中最常用的一种方法是Metropolis-Hastings，该方法最早由Metropolis于1953年给出的，后来Metropolis的算法由Hastings改进，合称为M-H算法[13-14]。M-H算法是MCMC的基础方法，由M-H算法演化出了许多新的抽样方法，包括目前在MCMC中最常用的Gibbs抽样也可以看做M-H算法的一个特例[15]。

1 符号描述和模型介绍

文中主要研究的是区间删失下的广义指数分布模型的建立及贝叶斯回归分析的应用。下面假设第i个个体满足以下关系：

假定每个个体都可以观测两次，其中U、V代表两个随机变量，并且以概率1满足U

文中采用的是广义指数分布对区间删失数据进行建模[8-9]。广义指数分布的密度函数为：

(1)

其分布函数为：

(2)

生存函数为:

(3)

风险函数为：

(4)

式中： α----形状参数;

λ----尺度参数。

形状参数为α，尺度参数为λ的广义指数分布记为GE(α,λ)。

其对应的全数据的似然函数为：

(5)

文中研究的是区间删失情况下的贝叶斯回归模型，则区间删失情况下的似然函数为：

(6)

故区间删失数据的对数似然函数可以表示为：

(7)

接下来建立贝叶斯层次模型如下：

(8)

α～gamma(1,0.001)

其中,λi指每个个体生存时间所服从的广义指数分布的尺度参数，βj,j=0,1,…,m的先验分布为正态分布，α的先验分布为gamma分布。

这样就可以建立起区间删失数据的广义指数分布贝叶斯回归模型。接下来可根据贝叶斯层次模型写出后验的联合密度函数,即后验似然函数[3,16]为：

(9)

故得到后验对数似然函数为：

(10)

接着，运用MCMC算法对参数进行估计。

2 数值模拟

用数值模拟过程来评价文中建立的模型性能，给出模拟步骤如下：

1)产生来自均匀分布U[-2,2]的N个独立同分布的x1,x2,…,xN。

2)设定β0=1，β1=1，α=1.5，并令λi=exp(β0+β1xi)。

3)产生N个服从广义指数分布的失效时间T,形状参数α=1.5，尺度参数λi=exp(β0+β1xi)。

4)产生N个服从参数为θ1=6的指数分布的第一次观测时间U,产生N个服从参数为θ2=0.2指数分布的第二次观测时间V,并满足U

5)比较U、V和失效时间T的大小关系，若TV，则令δ3=1，否则δ3=0。令δ2=1-δ1-δ3。

6)给出β和α的先验分布。并写出先验似然函数(LL)和后验似然函数(LP)。

7)应用MCMC算法估计参数β和α。

按照上述算法步骤，循环500次计算出待估参数β和α的均值和方差。样本量设定为N分别为200、300、500，模拟结果见表1。

表1 样本量为200,300,500的模拟结果

由表1 可以看出，在样本量不同,且左删失比例约为0.2，右删失比例约为0.4的情况下，模拟参数的估计值较真值偏差较小，能够给出对应参数较好的估计结果，并且精度会随着样本量的增加而增加，样本标准差也会随着样本量的增加而减小。由此可见，该模型用来进行后验估计是可行的。在算法的选择上也可采用其他的算法进行估计。

3 实证分析

对一个实际数据例子进行研究分析，选取的数据是1976年到1980年之间在波士顿进行乳腺癌早期治疗的回顾性研究数据。该数据由Finkelstein和Wolfe在1985年展现出来，数据是由94位病人组成，其中分为给予放射性治疗组(RT)和放射性疗法加辅助性化学治疗组(RCT)。放射治疗组共计46位病人，放疗加辅助化疗组共有48位病人[17]。

在研究过程中，病人每4～6个月随访一次，然而，病人的实际访问时间不同，每个病人的两次随访时间也是不同的。在就诊过程中医生会根据乳腺收缩程度来评估病人情况。这项研究的目的是为了比较这两组治疗方式对患者的治疗效果，看放疗辅助化疗方法是否可以提高患者的无复发率和总的生存率。但是有一些实验和临床证据表明，化疗加剧了正常组织对放射治疗的急性反应。这个数据包含了关于乳腺收缩的信息，但是没有精确的观测时间。这里有38例患者在研究期内没有明显的乳腺收缩，所以这部分观测设定为右删失数据，即这个区间观测没有右侧端点。对于其他患者，观测时间的时间间隔代表着在这段时间内发生过乳腺收缩。观测时间的左端点是从第一次诊所就诊时间开始，到最后一次就诊时发现乳腺收缩截止。例如，观测到的(6,10]表示在第6个月随访时患者未出现乳腺收缩，但是在下一次随访，即第10个月时，患者出现了乳腺收缩。乳腺收缩情况出现在第6个月至第10个月两次随访之间，但精确的时间未知。所以我们用区间的删失时间数据来描述乳腺收缩。将这组数据进行详细地分析估计，观测数据见表2。

在进行数据分析的过程中，若第i个病人属于放射治疗组，定义协变量xi=0；若第i个病人属于放疗辅助化疗组，定义协变量xi=1，并且假定乳腺癌发作时间服从广义指数分布。估计结果见表3。

通过表3的实验结果可以求出

λ=exp(-15.873x)

可以判断出：当病人属于放射治疗组时，λ=1；当病人属于放疗辅助化疗组时，0<λ<1。从而根据生存函数可以判断出,放疗辅助化疗方法可以提高患者的无复发率和总的生存率。

表2 乳腺癌观测数据

表3 乳腺癌数据估计结果

4 结 语

在贝叶斯框架下建立了服从广义指数分布的生存时间的尺度参数同影响生存时间的相关因素之间的回归模型，并给出后验似然函数，采用MCMC方法对后验似然函数进行求解最大值，同时解出了待估参数。并对该模型进行了模拟，模拟效果较好。并将该方法应用到乳腺癌数据例子中，结果表明，放射疗法辅助化疗方法对于提高患者的总的生存率有着一定的效果。

[1]GuptaRD,KunduD.Exponentiatedexponentialfamily;analternativetogammaandweibulldistributions[J].BiometricalJournal,2001,43(1):117-130.

[2] Gupta R D, Kundu D. Generalized exponential distribution: different method of estimations[J]. Journal of Statistical Computation and Simulation,2001,69(4):315-337.

[3] Ibrahim J G, Chen M H, Sinha D. Bayesian survival analysis[M]. New York: John Wiley & Sons Ltd.,2005.

[4] 李荣,朱慧明.删失试验寿命的贝叶斯威布尔生存回归模型[J].统计与决策,2006(24):20-22.

[5] 朱慧明,李荣,方博文.删失试验寿命的贝叶斯生存极值回归模型[J].系统工程与电子技术,2007,29(11):1988-1990.

[6] Upadhyay S K, Peshwani M. Posterior analysis of lognormal regression models using the Gibbs sampler[J]. Statist. Papers,2008,49:59-85.

[7] Puja Makkar, Puneet K, Srivastava R S Singh, et al. Bayesian survival analysis of head and neck cancer data using lognormal model[J]. Communication in Statistics-Theory and Methods,2014,43(2):392-407.

[8] Gupta R D, Kundu D. Generalized exponential distributions[J]. Australian and New Zealand Journal of Statistics,1999,41(2):173-188.

[9] Gupta R D, Kundu D. Generalized exponential distributions: statistical inferences[J]. Technical Report, The University of New Brunswick, Saint John.,1999,41(3):111-115.

[10] Kundu D, Gupta R D. Generalized exponential distribution: bayesian estimations[J]. Computational Statistics & Data Analysis,2008,52(4):1873-1883.

[11] 郭环.两参数广义指数分布的参数估计与数值模拟[D].武汉：华中科技大学,2013.

[12] 黄小艳.MCMC方法分析[J].中国市场,2015(14):185-186.

[13] Metropolis N, Rosenbluth A W, Rosenbluth M N, et al. Equation of state calculations by fast computing machines[J]. The Journal of Chemical Physics,1953,21(6):1087-1092.

[14] Hastings W K. Monte carlo sampling methods using markov chains and their applications[J]. Biometrika,1970,57(1):97-109.

[15] Geman S, Geman D. Stochastic relaxation, Gibbs distributions, and the Bayesian restoration of images[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence,1984(6):721-741.

[16] Bernardo J, Smith A. Bayesian theory[M]. West Sussex: John Wiley & Sons.,2000.

[17] Finkelstein D M, Ra W. A semiparametric model for regression analysis of interval-censored failure time data[J]. Biometrics,1985,41(4):933-945.

Bayesian survival regression analysis of interval censored data with generalized exponential Model

LI Qun, DONG Xiaogang, WANG Chunjie*, ZHAO Bo

(School of Basic Sciences, Changchun University of Technology, Changchun 130012， China)

Bayesian regression analysis model of interval censored lifetime under two-parameter Generalized Exponential is studied. Provided that the lifetime comes from generalized exponential distribution, and the prior distribution of shape parameter derives from the gamma distribution, the Bayesian regression model influenced by scale parameter and survival time is established to obtain the variation of lifetime. MCMC algorithm is used to estimate the parameters, and R software is used for simulation.

generalized exponential distribution; interval censored; bayesian regression; MCMC algorithm.

2016-07-19

国家自然科学基金青年基金项目(11301037)； 国家自然科学基金资助项目(11571051)； 吉林省教育厅“十三五”规划项目(2016317)

李 群(1991-)，女，汉族，山东菏泽人，长春工业大学硕士研究生，主要从事生存分析方向研究,E-mail:liqun91@live.com. *通讯作者：王纯杰(1978-)，女，汉族，辽宁辽阳人，长春工业大学副教授，博士，主要从事数理统计、生存分析方向研究,E-mail:wangchunjie@ccut.edu.cn.

10.15923/j.cnki.cn22-1382/t.2016.6.16

O 212.4

A

1674-1374(2016)06-0597-06