生存分析，顾名思义是指研究个体生存时间及其影响它们的因素。具有生存结果的研究类型包括临床试验、前瞻性和回顾性观察性研究以及动物实验。这里“个体的生存”可以推广抽象成某些关注的事件。 所以生存分析就成了研究某一事件与它的发生时间的联系的方法。
这个方法广泛的用在医学、生物学等学科上，近年来也越来越多人用在互联网数据挖掘中，例如去预测信息在社交网络的传播程度，或者去预测用户流失的概率。 相比于逻辑回归或常见的数据分析，生存分析有个很大的优点是可以处理删失数据（Censored Data， 只有个体的部分数据）。
生存研究可以包括估计生存分布，比较各种治疗或干预措施的生存分布，或阐明影响生存时间的因素。比如说医生想研究病人在服药后的健康状况，跟踪期为2015一整年。有些人接近2015年底才开始服药，那么他就只有很短的数据。如果在一年内知道病人死亡了，那么这个病人的数据是完全的；而健康的病人的数据到2015年末就是缺失的了，因为假如病人在一年后死亡了，我们并不知道这件事。或者说病人在一年内突然失联了，那么他的数据也是缺失，如图1 所示。如果直接删掉这个样本，对整体的数据信息是有损失的。考虑一个极端情况，如果数据集中的全部样本都是如上所述的情况，那怎么办呢？盲目的删除这些会提高模型的方差，以至于降低估计精度。这时逻辑回归就不太适用。
图1 在研究过程中删失数据的举例
生存分析在以下几个方面类似于线性和逻辑回归分析：
- （通常）有一个结果变量和一个或多个预测因子；
- 检验关于预测因子与结果变量关系的统计假设特别有趣；
- 调整混杂协变量至关重要；
- 通过分析残差和其他方法进行模型选择和检查假设是关键要求。
因此，我们应该熟悉经典假设检验的基本概念和回归分析的原则。熟悉分类数据分析方法，包括列联表、分层列联表、泊松和逻辑回归，也很重要。然而，生存分析与这些经典统计方法的不同之处在于，删失（Censoring）在几乎所有情况下都起着核心作用，而该主题的理论基础要复杂得多。
**起始事件(Initial Event)：**反应生存时间起始特征的事件，如疾病确诊、某种疾病治疗开始等。 **失效事件(Failure Event)：**在生存分析随访研究过程中，一部分研究对象可观察到死亡，可以得到准确的生存时间，它提供的信息是完全的，这种事件称为失效事件，也称之为死亡事件、终点事件。
终点事件和起始事件是相对而言的，它们都由研究目的决定，须在设计时明确规定，并在研究期间严格遵守，不能随意改变。
生存分析的响应通常被称为终点时间（Failure Time）、生存时间（Survival Time）或事件时间（Event Time）。
**生存时间（Survival Time）**广义上指某个起点事件开始到某个终点事件发生所经历的时间，度量单位可以是年、月、日、小时等，常用符号
删失（Censoring）是指开始时间或结束时间没有被精准观测，从而导致数据不完备的情况。当开始或结束的事件没有被精确地观察到时，就会出现。最常见的例子是右删失（Right Censoring），当最终端点只超过特定值的结果。
常见形式中，如果$T^$是一个表示失败时间的随机变量，$U$是一个表示删失事件发生时间的随机变量，我们观察到的是$T,,=,,min(T^,U)$和一个删失指标$\delta,,=,,I[T^*,<,U]$也就是说，根据$T$是删失时间还是观察到的失败时间，$\delta$分别为0或1。不常见的情况是，可能出现左删失（Left Censoring），已知事件发生在特定时间之前，或者区间删失（Interval Censoring），失败时间只发生在特定的时间间隔内。
生存数据包含以上数据，一般表示为${t_i,\delta_i,x_i}^n_{i=1}$，其中,
生存分析的目的就是估计生存分布，比较两个或两个以上的生存分布，或（更普遍地）评估许多因素对生存的影响。
以观察（随访）时间为横轴，以生存率为纵轴，将各个时间点所对应的生存率连接在一起的曲线图。 生存曲线是一条下降的曲线，分析时应注意曲线的高度和下降的坡度。平缓的生存曲线表示高生存率或较长生存期，陡峭的生存曲线表示低生存率或较短生存期。
生存分析方法依赖于生存分布，其中两种关键的指定方法是生存函数和风险函数。
**生存函数(Survival Function)**表示个体在t时刻还存活的概率，即事件在t时刻还未发生的概率。形式上， $$ S(t),=,pr(T,>,t),,0<t<\infty $$ 其中，$T$是生存时间。这个函数在时间为0时取值为1，随着时间的推移而减少（或保持恒定），当然永远不会降到0以下。正如这里的定义，它是连续的。 生存函数通常用风险函数来定义，即瞬时失效率。这个概率是，假设一个受试者已经存活到时间$t$，他或她在接下来的一小段时间段内失败，再除以该时间间隔的长度。在形式上，**风险函数(Hazard Function)**可以表示为 $$ h(t),=,lim_{\delta\to0}\frac {pr(t<T<t+\delta|T>t)}\delta $$ 这个函数也被称为强度函数(intensity function)或死亡率(force of mortality)。
风险和生存功能是指定生存分布的两种方式。为了说明两者之间的联系，首先考虑一个风险最初非常高的情况。这种危险可能适用于描述生命早期死亡率较高的动物的寿命。图2.1说明了该风险函数(a)和相应的生存函数(b)。接下来要考虑一下相反的，风险最初很低，随后生命后期增加。这样的风险将描述初始死亡风险较低的生物体。这一点如图2.1(c)和(d)中相应的生存率所示。
图2.1具有高初始风险(a&b)和低初始风险(c&d)的风险函数和生存函数
人口统计数据提供了风险和生存函数的另一个说明，如下一个示例所示，我们看到了高早期和高晚期风险函数的元素。
示例2.1 从1940年到2004年，每个日历年中按年龄划分的男性和女性的每日危险率包含在三维阵列“survexp.us”中，这是R包“survival”的一部分。这些风险率是通过使用Therneau 和 Offord 中描述的方法得出的美国生命表。图2.2显示了2004年美国男性和女性的寿命分布。这里用对数尺度绘制的风险图，显示了人类寿命的几个特征。首先，生命的最初几天和几周特别危险，风险在生命的第一个月后迅速消退。危险在青少年时期增加，然后达到最低水平，直到在中年时期开始稳步增加。众所周知，男性的死亡率高于女性。并显示了相应的生存函数。这个例子还表明，危险函数可能显示了在生存曲线中可能不明显的风险变化的细节。
图2.2 2004年美国男性和女性的风险函数和生存函数
> tm <- c(0,
1/365, # birth
7/365, # first day of life
28/365, # seventh day of life
1:110) # subsequent years
> hazMale <- survexp.us[,"male","2004"] # 2004 males
> hazFemale <- survexp.us[,"female","2004"] # 2004 females 除了生存和危险函数外，还有其他几种方法来定义生存分布。常用的生存分析之外通常使用的累积分布函数(CDF)， $$ F(t),=,pr(T\geq t),0<t<\infty $$ 这是生存函数的补充，就像生存函数一样，它也是连续的。在生存分析中，该函数被称为累积风险函数（不要与下面定义的累积风险相混淆）。概率密度函数(PDF)， $$ f(t),=,-\frac {d}{dt}S(t),=,\frac {d}{dt}F(t) $$ 是CDF的变化率，或负生存函数的变化率。风险函数与PDF和生存函数的关系为， $$ h(t),=,\frac {f(t)}{S(t)} $$ 也就是说，时间$t$时的风险是一个事件发生在时间$t$附近的概率除以受试者在时间$t$时生存的概率。累积风险函数定义为风险函数下直到时间t的区域，即， $$ H(t),=,\int^{t}_{0}h(u)du $$ 生存函数也可以根据风险函数来定义， $$ S(t),=,exp(-\int^{t}_0h(u)du)=exp\big(-H(t)) $$ 正是这种关系允许我们计算与危险函数对应的生存函数，如图2.1和图2.2。
有几种生存分布可用于建模生存数据。
指数分布：
风险函数$H(t),=,\int ^t_0h(u)du,=,\lambda t|^t_0,=,\lambda t$
因此，时间$t$ 处的累积风险只是图中阴影矩形的面积$ht$
图2.3 指数风险
生存函数$S(t),=,e^{-H(t)},=,e^{-\lambda t}$
性质：
风险率函数$h(t),=,\lambda$
概率密度函数$f(t),=,h(t)S(t),=,\lambda e^{-\lambda t}$
指数随机变量的平均值$E(T),=,\int^\infty_0S(t)dt,=,1/\lambda$ 方差$Var(T),=,\frac 1{\lambda^2}$
Weibull分布:
累计风险函数$H（t）,=,(\lambda t)^\alpha$
生存函数
性质：
风险率函数$h(t),=,\alpha\lambda(\lambda t)^{\alpha-1},=,\alpha\lambda^\alpha t^{\alpha-1}$
随机变量的平均值$E(T),=,\frac {\Gamma(1+1/\alpha)}\lambda$ 中位数
绘制韦伯风险函数
weibSurv <- function(t, shape, scale) pweibull(t, shape=shape, scale=scale, lower.tail=F)
curve(weibSurv(x, shape=1.5, scale=1/0.03), from=0, to=80, ylim=c(0,1), ylab=’Survival probability’, xlab=’Time’)weibHaz <- function(x, shape, scale) dweibull(x, shape=shape, scale=scale)/pweibull(x, shape=shape, scale=scale,
lower.tail=F)
curve(weibHaz(x, shape=1.5, scale=1/0.03), from=0, to=80, ylab=’Hazard’, xlab=’Time’, col="red")
图2.4 韦伯风险函数
>tt.weib <- rweibull(1000, shape=1.5, scale=1/0.03)
>mean(tt.weib)
[1] 31.35497
>median(tt.weib)
[1] 26.84281
>gamma(1 + 1/1.5)/0.03 # mean
[1] 30.09151
>(log(2)^(1/1.5))/0.03 # median
[1] 26.10733Gamma分布：
概率密度函数
绘制伽马风险函数
gammaHaz <- {function(x, shape, scale) dgamma(x, shape=shape, scale=scale)/pgamma(x, shape=shape, scale=scale, lower.tail=F)}
curve(gammaHaz(x, shape=1.5, scale=1/0.03), from=0, to=80, ylab=’Hazard’, xlab=’Time’, col="red")
图2.5 伽马风险函数
均匀分布：
生存函数$S(t),=,\frac{w-t}w$
性质：
风险率函数$h(t)=\frac 1{w-t}$
随机变量的期望$E(T)=\frac w2$ 方差$Var(T)=\frac{w^2}{12}$
Gompertz分布
风险率函数$h(t),=,Be^t,t\geq0,c>1,B>0$
生存函数$S(t),=,exp\big(\frac B{lnC}(1-c^t))$
Makeham分布
Makeham分布是对Gompertz分布的修正，Gompertz分布用幂函数对与年龄相关的危险率进行建模，但没有考虑到所有年龄段共有的一些死亡风险。
危险率函数$h(t),=,A+Be^t,t\geq0,c>1,B>0,A>-B$
生存函数$S(t),=,exp\big(\frac B{lnC}(1-c^t)-At)$
如果我们知道一个生存随机变量的风险函数，我们可以用方程推导出生存函数。对于某些参数族，这做得很简单。但如果风险函数更复杂，我们需要使用数值方法来计算积分。例如，再次考虑图2.1中人类危害函数。为了得到相应的生存图，我们首先计算一个差异向量“tm.diff”，然后使用“cumsum”函数找到累积风险函数，最后利用生存函数与累积风险的关系得到“survMale”和“survFemale”。
绘制这些生存函数与“tm”的结果
在下面的代码中，“tm.diff”是每个矩形的宽度，“survMale”和“survFemale”分别代表男性和女性的生存曲线。
>tm.diff <- diff(tm)
>survMale <- exp(-cumsum(hazMale*tm.diff)*365.24)
>survFemale <- exp(-cumsum(hazFemale*tm.diff)*365.24)
图2.6 2004年男性累积风险的计算说明
存储在“survexp.us”中的危险函数实际上是在1天、1周、1、月、1年以及之后每年评估的阶跃函数。这里显示的只是前两年，阴影区域代表了1.5年的累积风险
现在我们有了男性和女性的生存分布，我们可以计算出2004年美国男性和女性的平均死亡年龄，即图中各自生存曲线下的面积。
在下面的代码中，“tm.diff”是每个矩形的宽度，“survMale”和“survFemale”分别是男性和女性的矩形的高度。矩形的面积之和给出了积分的值：
>sum(survMale*tm.diff) # mean age of male death in 2004
[1] 73.8084
>sum(survFemale*tm.diff) # mean age of female death in 2004
[1] 78.90526Kaplan-Meier Estimator
我们已经看到，如果使用参数模型建模生存数据，则有各种各样的风险函数模型可以供选择。但是对于特定的应用程序应该使用哪个参数模型呢？当模拟人类或动物的生存时，很难知道要选择什么参数家族，而且通常没有一个可用的家族有足够的灵活性来模拟分布的实际形状。因此，在医疗和卫生应用中，非参数方法具有解释生物生存变化的灵活性，具有相当大的优势。
要了解生存函数的非参数估计量，其中最广泛使用的是乘积极限法，也被称为Kaplan-Meier估计，是存活到下一个失败时间的条件概率的失败时间的乘积。在形式上表示为， $$ \hat{S}(t),=,\prod_{t_i\leq t}(1-\hat{q_i}),=,\prod_{t_i\leq t}(1-\frac{d_i}{n_i}) $$ 其中$n_i$ 是在$t_i$ 时面临风险的受试者数量，$d_i$ 是当时失败的个体数量。
我们使用所谓的*“delta方法”*来得到方差。
乘积极限估计量的置信限 $$ var(\hat S(t)),\approx,[\hat S(t)]^2\sum_{t_i\leq t}\frac {d_i}{n_i(n_i-d_i)} $$ KM估计的部分渐进性质(Delta-Method) $$ V(\hat\lambda_i),=,\frac1{n_j}\hat\lambda_j(1-\hat\lambda_i),=,\frac{d_j(n_j-d_j)}{n^3_j} $$
这种方法不需要对生存时间的分布做出假定，但是却可以通过一个模型来分析生存时间的分布规律，以及危险因素对生存时间的影响，有以下几种方法可以用协变量来建模生存数据：
Parametric Families 参数族
一般的方法是选择我们已经讨论过的参数分布之一，并让该分布的参数依赖于协变量。一般来说，若有k个组，可以给每个组在一个家族中有自己的分布。这是一种可行的方法，但它并不完全节俭，也不容易被解释。
Accelerated Life 加速寿命
Proportional Odds 比例优势
Proportional Hazards 比例风险
对生存数据建模的另一种方法是假定协变量，这种方法的主要缺陷是其在所有持续时间内按比例地增加或减少风险。而cox比例风险模型可以同时评估几个因素对生存的影响，将允许我们将回归模型拟合到有删失的生存数据中，就像人们在线性和逻辑回归中可以做的那样。然而风险函数的简单关系转化为生存函数更复杂的关系，会使得生存建模更加不易。
Cox模型最初是由D.R. Cox所开发的，它可以被写为变量的风险对数的多元线性回归，而基线风险是随时间变化的“截距”项。为基线风险选择不同的参数形式会导致比例风险家族中出现不同的模型。
生存分析，顾名思义，其最广泛的应用是在医学领域，无论是临床领域对于某种慢性疾病的治疗做贡献，还是在公共卫生领域探寻感染流行性病的致死率，又或者是对于药物经济学的某种模型评价，甚至在动物医学领域无处不见这种分析方法的应用。目前更是被广泛拓展到工业、经济领域等其他领域中，来研究某个工业系统失效的时间或者某种产品客户的流失率。下面是生存分析几个主要的作用。
估计 ：
根据样本生存资料估计总体生存率及其它有关指标（中位生存期等）。具体在医学领域可以根据患者治疗后的生存时间资料，估计不同时间的生存率、生存曲线以及中位生存期等 。 比较 ：
对不同处理组生存率和生存时间进行比较，如比较不同疗法治疗某种疾病的生存率，以了解哪种治疗方案较优。
影响因素分析 ：
通过分析生存时间和协变量的相关性，探索影响生存时间长短的因素，或平衡某些因素影响后，研究某个或某些因素对生存率的影响。如了解影响患者痊愈的主要因素，包括患者的年龄、性别、病程、治疗方案等 。 预测 ：
具有不同因素水平的个体生存预测，如根据患者信息以及医生的治疗方案等预测该患者在某个时间的生存率。
[1] Applied Survival Analysis Using R. Dirk F.Moore, 2006.
[2] Aalen, O.: Nonparametric inference for a family of counting processes. Ann. Stat. 701–726 (1978)
[3] Barker, C.: The mean, median, and confidence intervals of the Kaplan-Meier survival estimate- computations and applications. Am. Stat. 63(1), 78–80 (2009)
[4]Parametric Survival Models. German Rodrıguez, 2001.