生存分析课程设计乘积极限法拟合样本的生存模型

1、乘积极限法拟合样本的生存模型摘要随着社会的进步和医学研究要求的不断提高，生存分析的应用范围不再仅仅 是字面上所理解的“生存分析”，更代表了一种处理缺失数据的基本分析思想。 它的研究内容主要包括两个方面:描述生存过程以及分析影响生存过程的因素。 由于生存分析在处理缺失数据上具有无法替代的作用，因此在临床试验中应用非 常广泛。随着统计软件的不断发展,生存分析的理论和应用将会越来越广泛和深 入，更多的符合生物医学实践的模型的建立将会越来越可行。生存分析是研究生 存现象和响应时间数据及其统计规律的一门学科。该学科在生物学、医学、保险 学、可靠性工程学、人口学、社会学、经济学等方面都有重要应用。目前这方

2、面 的教材，国内还不太多，而且大多局限于生存分析的某一特定应用领域。在现有 的几本教科书中，或者偏重于工程学，或者偏重于生物学和医学。本文主要通过 简单的例子对成绩极限法和Nelson-Aalen进行介绍和简单的比较。关键词：MATLAB；生存分析；乘积极限；Nelson-Aalen法PRODNCT-LIMITSURVIVAIMODELFITTINGSAMPLESABSTRACTAlong with social progress and the continuous improvement of medical research requirements, the scope of appl

3、ication of survival analysis is no longer just a literal understanding of survival analysis, but also represents a process of fundamental analysis censored thought. Its research mainly includes two aspects: process description and analysis of the factors affecting survival process. Because survival

4、analysis in dealing with censored data have irreplaceable, and therefore is widely used in clinical trials. With the continuous development of statistical software, the theory and application of survival analysis will be more extensive and in-depth, more in line with the practice established biomedi

5、cal model will be more feasible. Survival analysis is the study of the phenomenon of survival and response time data and the statistical laws of a subject. The discipline has important applications in biology, medicine, insurance, reliability engineering, demography, sociology, and economics. Curren

6、tly teaching in this area, not too much domestic and mostly confined to the survival of a specific application domain analysis. In the prior few textbooks, or emphasis on engineering, or emphasis on biology and medicine. In this paper, a simple example of a performance limit law and Nelson-Aalen int

7、roduction and a simple comparison.Key words: MATLAB ;Survival analysis; Product limit ;Nelson-Aalen method目录 TOC o 1-5 h z HYPERLINK l bookmark13 o Current Document 1问题提出21.1问题的重述21.2问题的提出2 HYPERLINK l bookmark19 o Current Document 2问题分析22.1相关知识简介22.1.1生存分析问题22.1.3符号说明22.1.3乘积极限法知识介绍32.1.4Nelson-Aal

8、en 法介绍42.1.5软件知识52.2程序编写前分析6 HYPERLINK l bookmark115 o Current Document 3问题假设6 HYPERLINK l bookmark123 o Current Document 4问题的求解64.1乘积极限求解74.2 Nelson-Aalen 法求解74.3两种求解方法的比较9 HYPERLINK l bookmark162 o Current Document 参考文献10 HYPERLINK l bookmark170 o Current Document 附录111问题提出1.1问题的重述现有一个由10只实验鼠组成的样本

9、，直到实验鼠全部死亡，以天为单位记 录发生死亡的时间如下：2，3，4，5，7，8，8，9，9，12。利用要求的方法求估计S(10)；用要求的方法估计A(l。)；又知道S=心。1.2问题的提出运用乘积极限法估计s(10)；用 Nelson-Aalen 法估计 A(10)；S*，由2)可得s(10)的又一估计值，与1)作比较，并说明原因。2问题分析2.1相关知识简介2.1.1生存分析问题生存分析是研究生存现象和响应时间数据及其统计规律的一门学科。该学科 在生物学、医学、保险学、可靠性工程学、人口学、社会学、经济学等方面都有 重要应用。随着社会的进步和医学研究要求的不断提高，生 存分析的应用范围不

10、再仅仅是字面上的“生存分析”，更代表了一种处理删失数据的基本分析思想。 它的研究内容主要包括两个方面:描述生存过程以及分析影响生存过程的因素。 2.1.3符号说明七 为在第j个子区间终点处同时发生死亡的人数；。为恰好在第j个死亡点之前的生存人数；qj为在第j个子区间内生存人数的死亡概率；S（t）为生存分布函数;a)为累计危险函数；。为第m个时间点；为样本含量2.1.3乘积极限法知识介绍乘积极限法，即kaplan-meier法。适用范围：小样本未分组资料、大样本 未分组资料，可利用概率乘法定理来计算生存率。s。)的估计：假设对一个有个对象的样本中的个体都从，=0开始观察，并在对所有的 对象死亡之

11、前终止研究，于是得到不同的时间处样本中终止的人数。在临床研究中，习惯上采用在每一死亡点进行分割的方法。如果r是恰好第j个 死亡点之前的风险集合(一般而言。可发生dj次死亡，但多数情况下令dj T)，a d.qj则由式j rj求得qj的估计量是在终止于第j个死亡点的那个区间开始时生 存，在那个区间内死亡的概率。二H因此，s(tm TOC o 1-5 h z v6)q (-)可以估计出由t = 0到。的生存概率。对于于匕, m+1”I”的估计量与相同，因此样本在。于t之间没有发生死亡事件，从而有 HYPERLINK l bookmark31 o Current Document TZTr d s(

12、)=H Tj ,t t t ,m = 1,2, mm+1j=】Ij )特别地，当t 时，kS (t)= 1。当tm - t tm+1时kS (t)的方差可以由下式近似表示Vari、：(t) S(t)2 IEj=1q jp raa ()如果用七代替qj，P代替P，S(t)代替S()VarS (t) A S O2芝-rj=1 j则得 d F，如果不存在终止者(即完整数据的情况)，则对于所有的j，+1 = rj - dj，显然S ()=rmj=1rm 11r1，此处r=n为初始样本，rm+1为恰好在第m+1个死亡点之前的生存人数(或称 为第m个死亡点之后的生存人数)，S心 就由生存到那一时刻的生存人

13、数与n 的比例来确定。2.1.4Nelson-Aalen 法介绍由累积危险函数定义 Xt) = jtX (y)dy = -InS(t)中定义S(t) =e-八(t),其中 0A(t)称为累积危险率函数。因此也可以通过首先估A(t)来估计S(t)，定义S (t) = e-t)于是-dj) rj入(t) = - ln S (t)又S(t)用(rmj=1将其带入入(t) = - ln S (t)，得大(t)=-In FIr -d I Ij J=-E ln(1j=1d-),t t t , m = 1,2,j因为-ln(1- d-j-) = d + 2()2+.，忽略二次项和更高次数的项，得jjj近似累

14、积危险率估计量为A (t)=于七 r jT j于是，得到生存函数的近似-弟S (t) = e j=1 rj该近似量称作Nelson-Aalen估计量2.1.5软件知识MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据 可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括 MATLAB和Simulink两大部分。MATLAB是matrix&laboratory两个词的组合，意为矩阵工厂（矩阵实验室）。 是由美国mathworks公司发布的主要面对科学计算、可视化以及交互式程序设计 的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态

15、系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研 究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决 方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C、Fortran）的 编辑模式，代表了当今国际科学计算软件的先进水平。MATLAB和Mathematical Maple并称为三大数学软件。它在数学类科技应用 软件中在数值计算方面首屈一指。MATLAB可以进行矩阵运算、绘制函数和数据、 实现算法、创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、 控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。MATLAB的基本

16、数据单位是矩阵，它的指令表达式与数学、工程中常用的形 式十分相似，故用MATLAB来解算问题要比用C，FORTRAN等语言完成相同的事情 简捷得多，并且MATLAB也吸收了像Maple等软件的优点，使MATLAB成为一个强大的数学软件。在新的版本中也加入了对C, FORTRAN, C+, JAVA的支持。2.2程序编写前分析由题意可知10只实验鼠发生死亡的时间如下：2, 3, 4, 5, 7, 8, 8, 9, 9, 12。则死亡时间点为：2, 3, 4, 5, 7, 8, 9, 12。则得到：d = 0, d = 1, d = 1, d = 1, d = 1,d = 0, d = 1, d

17、= 2, d = 2, d = 0, d = 0, d = 1.123456789101112r = 10, r = r - d , j = 1,2,3, ,11a dq = jjrj又可以得到A (t )= r j= j引入新参数令b = 1 - q , j = 1,2,3,.,12注：为了利于循环程序的编写在没有死亡的时间点j上令d = 03问题假设1、实验鼠的健康状况正常且处于同一健康水平；2、实验鼠所处的环境相同；3、实验鼠的死亡时间记录准确无误；4、实验鼠的死亡没有外力和人为干扰；5、实验鼠相互之间没有影响。4问题的求解4.1乘积极限求解 运行成绩极限法求解源程序（源程序见附录）得到

18、:Command WindowtoCoLunns1 through 10a.LDOO11L1.0. 12500. ：li29D. 16670. iOOO0. 0067ColuminE 11 through 121.0000Colunns 1 throagh 101.00000. 9000D. SS990. 8760Cl 96711. 00000.60000.3333l.OGOCiColuniis 11 throuah 12h 0000Columns 1 thTouah 10:L 00000. 90000. 8000O.OOO0.60000.60000. 50000.30000.100Q1110

19、0010Columns 11 through 120,1000由运行结果只 S（10）= 0.1000 其中 b = 1 -q , j = 1,2,3,.,124.2 Nelson-Aalen 法求解运行Nelson-Aalen法求解源程序（源程序见附录）得到:d11 10122001T1010997e 6911 1qthrough 100. LOGO 0. 1111through 121. 0000Q.12500. U29a0. 1657o. iaaa0u56670E -Cd Lunns L thjoujh 1000.10000. 21110. 33610. 47900.47900.6456

20、1.04561.71231.7123CoLlhuis LLthrough iz|1.71232.7123S 5 LthTGUgh IQi.aooo0.90490.90970. 71450.61940.619 +0.5243D. 36150.1004D. 1804CdLunns 11thTough L20.18040.D664结果如图所示：其中B的值即为A（t）的值S（10）= 0.18044.3两种求解方法的比较有以上结果知Nelson-Aalen法的估计值比乘积极限法的估计值偏大。原因如下:乘积极限法公式如下所示，较为精确未取近似值% rj)而相应的Nelson-Aalen法的计算式如下s ()=Wj=irm十1ri ，因-md)工+将)2 +.，忽略二次项和更高次数的项而使得 jrrjjddln(1)铝，rrj j又因为近似累积危险率估计量为估计偏小入(t) = Y d r j=1 j生存函数的近似为减函