实验用家蝇寿命分布的两个参数

实验用家蝇寿命分布的两个参数

家蝇寿命模型建立家蚕是中医学和生物实验中最常用的实验动物之一。许多新药物的效果首先通过对家蚕的实验初步证实。在这些实验中，有必要处理家蚕的生命数据。在文献中，由于不知道家蚕生命分布的分布类型，所以有两种常见的方法。其中之一是作为正态分布来处理的。而采用非参数方法进行处理。正态分布家蚕的生命数据，并武断地进行处理。例如，在两份长度测试中使用t-验证（见），这是不可取的方法。t-验证可以在大样本中确保第一次错误的概率，但验证的效率不一定得到保证。在小样本中，甚至第一次错误的概率，而是以参数方法来确保每个人都失去的信息。任何一种学习方法都是从所有空白的地方开始的，这尤其是对像家蚕这样的实验不太可能少另一方面，从大量研究事实来看，动物、植物和材料的死亡必须基于一定的机制。因此，生命分布的类型相对固定，因此可以使用这种方法。我们直接基于数据来寻找家蝇寿命的分布类型.数据来源于某种药物对家蝇寿命影响的实验数据,实验分为对照组,处理组Ⅰ(低剂量用药组)和处理组Ⅱ(高剂量用药组).在每个组里将雌雄分别处理,因此共有6个组.本文的是这样安排的,先用Kaplan-Meier方法给出了生存函数的非参数估计及相应的死亡力度函数的估计,然后据此导出家蝇的参数生存模型并对模型中的参数给出了生物学的解释,之后讨论删失样本下模型参数的估计方法,最后通过图形直观地描述了模型对数据的拟合状况,文章的结尾是一些注释.生存函数的描述设X是表示寿命的随机变量,则相应的生存函数定义为Sx(t)=P(X>t).生存函数的估计采用Kaplan-Meier乘积限估计的方法(见例如).寿命试验的数据记录为:(X1,δ1),…,(Xn,δn),其中δk=0或1视xk是截尾寿命还是真实寿命而定.为得到估计,将数据(X1,δ1),…,(Xn,δn)整理为(t1,n1,d1),…,(tk,nk,dk),其中tj表示死亡时间(假设t1≤t2≤…≤tk),dj是在ti时刻死亡的个体数目,nj表示tj时刻面临风险的个体数目(包括在tj时刻死亡的个体数目).下面的六个图形是我们的六组实验数据的生存函数的Kaplan-Meier估计(横轴为时间,单位是天),图中的“”表示在该处有删失的数据.从这些图中看出,所得六组实验的生存函数已经具有基本的形状,曲线的跃度较小,因此,我们可以认为,所取样本量用于寻找生存函数的分布类型是合适的.除生存函数外,死亡力度是生存分析中一个具有基本重要性的概念,其定义为死亡力度在时间点(ti+ti+1)/2处的估计为:这一估计的直观解释是:ni个个体在时间区间(ti,ti+1)内平均单位时间内死亡的频率.根据频率估计概率的思想,这样的估计是直观合理的.下面是对数死亡力度的图(横轴是时间,单位是天),图中“+”为对数死亡力度对时间的散点,直线是这些散点对时间拟合的回归直线.为0时的gomprtz模型从图2中可以很明显地看出,死亡力度函数的对数可以用线性函数来拟合,即:根据知道,生存函数为其中λ=eb,这是正半轴上的Gompertz分布.当a=0时,理解为α趋于0时的极限,即指数分布S(t)=e-λx.此外,还应有α≥0,否则会有S(+∞)>0,这是不合理的.相应的密度函数:.同样,当α=0时,理解为α趋于0时的极限,即指数分布f(t)=λe-λx.很显然,在Gompertz模型中,因此,在此模型中,参数λ是初始的死亡力度,家蝇的性别,羽化时的体重等因素都可能对λ的值发生影响;而参数α则是死亡力度在单位时间内增加的百分比,也就是说,在Gompertz生存模型中,在单位时间内,个体的死亡力度以固定的百分比增加,死亡加速因子是常数.我们在这里应该注意到,习惯上,Gompertz模型常写成位置-刻度模型的形式.但是在讨论半轴上的Gompertz模型时,这样的写法是不可取的,这会丢掉指数分布这一重要特例.换句话说,正半轴上的Gompertz模型不再具有位置-刻度的形式.刻失数据下gomprtz的二乘估计根据正半轴上Gompertz分布的两个参数的生物学解释,寻找这两个参数的矩估计是不现实的,我们这里叙述删失数据下Gompertz分布的极大似然估计和简单最小二乘估计.(1)极大似然最大解关于Gompertz模型的极大似然估计,在文献中有比较详细的叙述,正半轴上Gompertz模型的极大似然估计与所述相似.为了后文说明及引用方便,我们将其简要叙述如下.设样本为X1,…,Xr,Xr+1,…,Xn,其中,前r个为观察到的真实寿命,后n-r个为删失寿命.这时模型的参数的对数似然函数为:相应的似然方程为:如果则似然方程的解为,=0.否则,从第一个似然方程得代入第二个似然方程,得:该方程的解只能利用数值方法得到(比如说Newton-Raphson方法).从该方程得到α的极大然似估计为,相应就得到λ的极大似然估计为关于极大似然估计,可以证明:(1)方程(3)的解存在且唯一;(2)>0的充分必要条件是(2)全实轴gomprtz模型的最小二乘估计正如上面所述,似然方程的求解只有通过计算机的帮助才能实现,其次,似然方程的NewtonRaphson方法求解依赖于良好的迭代初值的选取,因此,发展一种简单的估计参数的方法是必要的,这可以得到良好的用来求解似然方程的初值.在全实轴的Gompertz模型中,对生存函数取两次对数就得到线性函数,据此导出参数的最小二乘估计.当讨论正半轴的Gompertz模型时,情形稍有不同,此时的最小二乘估计充分利用了对数死亡力度函数是时间的线性函数这一特性.根据最小二乘估计的公式有:其中,这样,得到参数的简单最小二乘估计为:下面的表给出了这两种方法在六组数据下的数值结果.从表中可以看出,两种方法估计的结果实际上相差不远.参数模型估计的拟合情况下面给出的图形表示了Kaplan-Meier估计与参数模型估计的拟合情况.在图中,阶梯函数是Kaplan-Meier估计的结果而光滑曲线是参数模型估计的结果,从图中可以看出,参数模型与数据的拟合非常良好.mprtz生存模型与苏氏最优模型的关系从数据出发,给出了实验用家蝇生存函数的分布类型,正是正半轴上的Gompertz分布.Gompertz分布的重要性是勿容置疑的,比如就有这方面的论述.但一般说来,在生存分析中引入Gompertz分布主要是用来估计Weibull分布的参数,而不是直接用来作为生存模型使用.我们在这里不但给出了正半轴的Gompertz生存模型的实例而且对其中的参数给出了生物学上的解释.从这一点来看,正半轴的Gompertz分布与全实轴上的Gompertz分布还是有着些许的不同.此外,正半轴的Gompertz分布与Weibull分布又有着某些相似性,这就是,正半轴的Gompertz生存模型的死亡力度函数的对数是年龄的线性函数,而Weibull生存模型的死亡力度函数的对数是对数年龄的线性函数.从这一点出发,假如说我们考虑这样的生存模型类,即死亡力度的对数可以表示为α+bg(x)的情形,其中g(x)是已知函数,α,b是未知参数,则指数生存模型,正半轴上的Gompertz生存模型以及Weibull生存模型皆属此类.在药效实验中,需要对家蝇施以一定剂量的药物,以观察该药对家蝇寿命的影响,这实际上是个假设检验问题.根据本文前面对参数λ和α给出的生物学解释,λ是生命初始时刻的死亡力度,这是不会被任何处理影响的量,因此所谓药物效果仅仅是作用在参数α上,这样来看,药物的效果是影响了死亡加速因子α,因而相应的假设检验问题就只需在参数α上施行就行了.由于篇