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

乘积极限法拟合样本的生存模型

摘要

随着社会的进步和医学研究要求的不断提高,生存分析的应用范围不再仅仅是字面上所理解的“生存分析”,更代表了一种处理缺失数据的基本分析思想。它的研究内容主要包括两个方面:描述生存过程以及分析影响生存过程的因素。由于生存分析在处理缺失数据上具有无法替代的作用,因此在临床试验中应用非常广泛。随着统计软件的不断发展,生存分析的理论和应用将会越来越广泛和深入,更多的符合生物医学实践的模型的建立将会越来越可行。生存分析是研究生存现象和响应时间数据及其统计规律的一门学科。该学科在生物学、医学、保险学、可靠性工程学、人口学、社会学、经济学等方面都有重要应用。目前这方面的教材，国内还不太多，而且大多局限于生存分析的某一特定应用领域。在现有的几本教科书中，或者偏重于工程学，或者偏重于生物学和医学。本文主要通过简单的例子对成绩极限法和Nelson-Aalen进行介绍和简单的比较。

关键词：MATLAB；生存分析；乘积极限；Nelson-Aalen法

PRODNCT-LIMITSURVIVAIMODELFITTINGSAMPLES

ABSTRACT

Along with social progress and the continuous improvement of medical research requirements, the scope of application 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 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 biomedical 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. Currently 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 introduction and a simple comparison.

Key words:MATLAB ;Survival analysis; Product limit;Nelson-Aalen method

目录

1 问题提出 (2)

1.1 问题的重述 (2)

1.2问题的提出 (2)

2 问题分析 (2)

2.1相关知识简介 (2)

2.1.1生存分析问题 (2)

2.1.3符号说明 (2)

2.1.3乘积极限法知识介绍 (3)

2.1.4Nelson-Aalen法介绍 (4)

2.1.5软件知识 (5)

2.2程序编写前分析 (6)

3 问题假设 (6)

4 问题的求解 (6)

4.1 乘积极限求解 (7)

4.2 Nelson-Aalen法求解 (7)

4.3 两种求解方法的比较 (9)

参考文献 (10)

附录 (11)

1 问题提出

1.1 问题的重述

现有一个由10只实验鼠组成的样本，直到实验鼠全部死亡，以天为单位记录发生死亡的时间如下：2，3，4，5，7，8，8，9，9，12。

利用要求的方法求估计(10)S ；用要求的方法估计(10)Λ；又知道()

()t S t e

-Λ=。

1.2问题的提出

1.运用乘积极限法估计(10)S ；

2.用Nelson-Aalen 法估计(10)Λ；

3.()

()t S t e

-Λ=，由2）可得(10)S 的又一估计值，与1）作比较，并说明原因。

2 问题分析

2.1相关知识简介 2.1.1生存分析问题

生存分析是研究生存现象和响应时间数据及其统计规律的一门学科。该学科在生物学、医学、保险学、可靠性工程学、人口学、社会学、经济学等方面都有重要应用。随着社会的进步和医学研究要求的不断提高,生 存分析的应用范围不再仅仅是字面上的“生存分析”,更代表了一种处理删失数据的基本分析思想。它的研究内容主要包括两个方面:描述生存过程以及分析影响生存过程的因素。 2.1.3符号说明

j d 为在第j 个子区间终点处同时发生死亡的人数；

j

r 为恰好在第j 个死亡点之前的生存人数； j

q 为在第j 个子区间内生存人数的死亡概率； ()

S t 为生存分布函数；

()

t Λ 为累计危险函数；

m t 为第m 个时间点；

n 为样本含量

2.1.3乘积极限法知识介绍

乘积极限法，即kaplan-meier 法。适用范围：小样本未分组资料、大样本未分组资料，可利用概率乘法定理来计算生存率。

()

S t 的估计：

假设对一个有n 个对象的样本中的个体都从0t = 开始观察，并在对所有的

对象死亡之前终止研究，于是得到不同的时间t 处样本中终止的人数。

在临床研究中，习惯上采用在每一死亡点进行分割的方法。如果j r

是恰好第j 个死亡点之前的风险集合（一般而言。可发生

j

d 次死亡，但多数情况下令

1

j d =），

则由式

^

j j j

d q r =

求得

j

q 的估计量是在终止于第j 个死亡点的那个区间开始时生

存，在那个区间内死亡的概率。 因此，

()^

^

111m

m j j

m j j j j r d S t q r ==??

-??=-=

? ? ?????∏∏。

可以估计出由0t = 到m t 的生存概率。对于

[]1,m m t t t +∈，()S t 的估计量与()^

m S t 相

同，因此样本在m t 于t 之间没有发生死亡事件，从而有

()^

11,,1,2,,m

j j

m m j j

r d S t t t t m r +=??

-=≤<=

? ???∏ 。

特别地，当1t t < 时，()^

1

S t =。

当1m m t t t +≤<时

()

^

S t 的方差可以由下式近似表示

(){}()^

2

1

m

j

j j j

j

q Var S t r S t p

r =≈????

?∑，

如果用

^

j

q 代替

j

q ，

^

j

p 代替

j

p ，

()

^

S t 代替()S t ，则得

(){}

()(

)^

2

1

m j j j j

j j

d Var S t r S t r r d

=≈????-∑，

如果不存在终止者（即完整数据的情况），则对于所有的1,j j j

j r r d +=-，

显然

()^

11

11m

j m m j j r r S t r r ++=??== ? ???∏，

此处1r n = 为初始样本，1m r + 为恰好在第1m +个死亡点之前的生存人数（或称为第m 个死亡点之后的生存人数），()

m S t 就由生存到那一时刻的生存人数与n

的比例来确定。

2.1.4Nelson-Aalen 法介绍

由累积危险函数定义0()()ln ()t

t y dy S t λ∧==-?中定义()()t S t e -∧=,其中

()t Λ称为累积危险率函数。因此也可以通过首先估()t Λ来估计()S t ,定义

()

()t S t e

∧

∧

-∧=

于是

()ln ()t S t ∧

∧

∧=-

又^

1

()(

)m

j j

m j j

r d S t r =-=∏

将其带入()ln ()t S t ∧∧

∧=- ,得

()^

111ln ln(1),,1,2,,m m j j

j

m m j j j

j r d d t t t t m r r +==????-∧=-=--≤<=??

? ???????

∑∏ 因为21ln(1)(),2j j

j

j

j j

d d d r r r --

=

++???忽略二次项和更高次数的项，得

ln(1),j j j

j

d d r r --

≈

近似累积危险率估计量为

1

()m

j j j

d t r ∧

=∧=∑

于是，得到生存函数的近似

1

()m

j

j j d

r S t e

=-

∧

∑=

该近似量称作Nelson-Aalen 估计量

2.1.5软件知识

MATLAB 是美国MathWorks 公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB 和Simulink 两大部分。

MATLAB 是matrix&laboratory 两个词的组合，意为矩阵工厂（矩阵实验室）。是由美国mathworks 公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C 、Fortran ）的编辑模式，代表了当今国际科学计算软件的先进水平。

MATLAB 和Mathematica 、Maple 并称为三大数学软件。它在数学类科技应用软件中在数值计算方面首屈一指。MATLAB 可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。

MATLAB 的基本数据单位是矩阵，它的指令表达式与数学、工程中常用的形式十分相似，故用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。 则得到：

12345678910111211^

0,1,1,1,1,0,1,2,2,0,0, 1.10,,1,2,3,,11j j j j j j

d d d d d d d d d d d d r r r d j d q r +==============-==

又可以得到

()^1

m

j j j

d t r =Λ=∑

引入新参数令

1,1,2,3,...,12

j j b q j =-=

注：为了利于循环程序的编写在没有死亡的时间点j 上令0j d =

3 问题假设

1、实验鼠的健康状况正常且处于同一健康水平；

2、实验鼠所处的环境相同；

3、实验鼠的死亡时间记录准确无误；

4、实验鼠的死亡没有外力和人为干扰；

5、实验鼠相互之间没有影响。

4问题的求解

4.1 乘积极限求解

运行成绩极限法求解源程序（源程序见附录）得到：

由运行结果只()100.1000S = 其中1,1,2,3,...,12j j b q j =-=

4.2 Nelson-Aalen 法求解

运行Nelson-Aalen法求解源程序（源程序见附录）得到：

结果如图所示：

Λ的值

其中B的值即为()t

()100.1804

S=

4.3 两种求解方法的比较

有以上结果知Nelson-Aalen 法的估计值比乘积极限法的估计值偏大。 原因如下：

乘积极限法公式如下所示，较为精确未取近似值

()^

11

11m

j m m j j r r S t r r ++=??== ? ???∏，

而相应的Nelson-Aalen 法的计算式如下 因21ln(1)(),2j j

j

j

j j

d d d r r r --

=++???,忽略二次项和更高次数的项，而使得 ln(1),j j j

j

d d r r --

≈

估计偏小

又因为近似累积危险率估计量为

1

()m

j j j

d t r ∧

=∧=∑

生存函数的近似为减函数

1

()m

j

j j d

r S t e

=-

∧

∑=

所以用此方法求的S 偏大

参考文献

[1] 李晓林，孙佳美．生命表基础[M]．北京：中国财政经济出版社，2006：1-133.

[2] 刘国卫.MATLAB程序设计教程(第二版) [M].北京:中国水利水电出版社，2010：136-153.

[3] 张文彤，邝春伟.SPSS统计分析基础教程（第二版）[M].北京:高等教育出版社 ,2011:1-417.

[4] 姜启源，谢金星，叶俊编．数学模型（第三版）[M]．北京：高等教育出版社，2005：1-202.

[5] 刘琼荪，何中市．数学实验（第一版）[M]．北京：高等教育出版社，2004.01：1-247.

附录S t的源程序

求()

clc

clear

d(1)=0;d(2)=1;d(3)=1;d(4)=1;d(5)=1;d(6)=0;

d(7)=1;d(8)=2;d(9)=2;d(10)=0;d(11)=0;d(12)=1; r(1)=10;

for i=1:1:11

r(i+1)=r(i)-d(i);

end

q=d./r

b=1-q

s=1;

for i=1:12

S(i)=b(i)*s;

s=S(i);

end

Nelson-Aalen法求解程序

clc

clear

d(1)=0;d(2)=1;d(3)=1;d(4)=1;d(5)=1;d(6)=0;

d(7)=1;d(8)=2;d(9)=2;d(10)=0;d(11)=0;d(12)=1; r(1)=10;

for i=1:1:11

r(i+1)=r(i)-d(i);

end

d

r

q=d./r

b=0;

for i=1:12

B(i)=q(i)+b;

b=B(i);

S(i)=exp(-B(i));

end

B

S