Python生存分析实战：从Kaplan-Meier到Cox模型，用lifelines库处理删失数据

1. 生存分析从概念到Python实战的桥梁如果你在数据分析、医学研究或者金融风控领域工作大概率听说过“生存分析”这个词。它听起来有点专业甚至带点医学统计的严肃感但它的核心思想其实非常直观研究某个事件发生前所经历的时间。这个“事件”可以是任何你关心的终点——病人的复发、设备的故障、客户的流失、贷款人的违约。传统的数据分析方法比如线性回归或者分类模型在处理这类数据时往往会遇到一个棘手的问题删失数据。简单说就是你的观察期结束了但有些个体还没发生那个“事件”。比如一项为期5年的临床试验结束了一部分患者依然健在你不知道他们最终会在何时去世。粗暴地忽略这些数据会损失信息而错误地标记为“未发生事件”又会引入偏差。生存分析就是专门为解决这个问题而生的统计方法。近年来随着Python在数据科学领域的统治地位生存分析的门槛大大降低。你不再需要依赖昂贵的专业统计软件如SAS、SPSS或者艰深的R语言编程。lifelines库的出现让用Python进行专业级的生存分析变得像调用scikit-learn拟合一个模型一样简单。这个库几乎囊括了从最基础的Kaplan-Meier曲线、Cox比例风险模型到更复杂的参数模型、竞争风险模型等一系列工具。它提供了清晰统一的API、优秀的可视化支持并且能很好地与pandas、numpy、matplotlib等生态无缝集成。对于数据科学家、量化分析师和临床研究员来说掌握lifelines意味着你拥有了一把解开时间-事件数据奥秘的瑞士军刀。这篇内容就是为你准备的lifelines实战入门指南。无论你是刚接触生存分析概念的新手还是熟悉理论但苦于没有合适Python工具的中级从业者都能从这里找到一条清晰的路径。我们将避开繁复的数学推导聚焦于如何用代码解决实际问题。我会带你从数据准备开始一步步完成描述性分析、模型建立、结果解读和可视化并分享那些官方文档里不会写的、我在实际项目中踩过的坑和总结的技巧。我们的目标很明确让你在读完这篇文章后能独立地使用lifelines库处理你自己的生存数据。2. 核心工具箱lifelines库全景解析与安装配置在深入代码之前我们有必要对lifelines这个工具箱有一个全景式的认识。它不是一个单一功能的脚本而是一个模块化、层次分明的生态系统。理解它的结构能帮助你在面对具体问题时快速找到正确的“工具”。2.1 lifelines的核心模块与功能地图lifelines库主要围绕几个核心的生存分析任务进行组织非参数估计这是生存分析的起点用于在不假设任何特定分布的情况下估计生存函数。核心类是KaplanMeierFitter。它最经典的应用就是绘制Kaplan-Meier生存曲线直观比较不同组如治疗组 vs 对照组的生存状况。lifelines还提供了NelsonAalenFitter用于估计累积风险函数。半参数回归模型这是生存分析中最常用、最强大的工具即Cox比例风险模型。在lifelines中对应的类是CoxPHFitter。它允许你探究多个协变量如年龄、性别、治疗方案对生存风险的影响同时不指定基线风险的分布形式灵活性极高。参数回归模型当你对事件发生时间的分布有一定先验知识时例如你知道失效时间可能服从指数分布或威布尔分布参数模型可能更有效率。lifelines提供了WeibullAFTFitter,LogNormalAFTFitter,LogLogisticAFTFitter等加速失效时间模型。它们直接对生存时间的分布进行建模结果有时更容易解释。竞争风险模型在现实世界中个体可能面临多种类型的“失败”。例如病人可能因癌症死亡我们关心的事件也可能因心脏病或其他原因死亡竞争事件。AalenJohansenFitter和CoxPHFitter在特定配置下可以用于处理这类问题。统计检验与模型评估lifelines.statistics模块提供了如log-rank检验用于比较两条或多条生存曲线是否有统计学差异。对于回归模型库内置了评估比例风险假设、模型拟合优度如Concordance Index以及进行残差分析的工具。实用工具与数据处理lifelines提供了datasets模块加载经典生存数据用于练习utils模块包含如datetimes_to_durations这样的函数帮助你将现实中的日期时间数据转换为生存分析所需的“持续时间”和“事件指示”格式。2.2 环境搭建与避坑指南安装lifelines非常简单通常一行命令即可pip install lifelines或者使用condaconda install -c conda-forge lifelines注意lifelines依赖numpy,pandas,scipy,matplotlib等科学计算栈。如果你使用的是纯净的Python环境pip install lifelines会自动安装这些依赖。但如果你在某个特定的虚拟环境或容器中工作有时版本冲突会导致问题。一个常见的坑是pandas版本过高可能与lifelines的某些内部函数不兼容。如果遇到奇怪的报错可以尝试指定稍旧一点的pandas版本例如pip install pandas1.5.3。安装完成后在Jupyter Notebook或Python脚本中导入常用组件import pandas as pd import numpy as np import matplotlib.pyplot as plt # 设置中文显示如果需要 plt.rcParams[font.sans-serif] [SimHei, Arial Unicode MS] # 用于Windows plt.rcParams[axes.unicode_minus] False # 导入lifelines核心模块 from lifelines import KaplanMeierFitter, CoxPHFitter, WeibullAFTFitter from lifelines.statistics import logrank_test from lifelines.datasets import load_rossi # 示例数据集 from lifelines.utils import datetimes_to_durations, concordance_index这里有一个实操心得我强烈建议在Jupyter Notebook中进行生存分析的探索性工作。生存分析的结果尤其是生存曲线和模型摘要具有极强的可视化属性。在Notebook中能够即时看到图表并交互式地调整参数对于理解数据和模型行为至关重要。lifelines的.plot()方法返回的就是matplotlib的轴对象你可以很方便地在此基础上进行自定义美化。3. 数据准备生存分析格式的标准化转换生存分析模型对输入数据的格式有严格且统一的要求。无论你原始的数据来自CSV、数据库还是API在送入lifelines之前都必须转换为“标准格式”。这一步是后续所有分析的基础也是最容易出错的地方。3.1 理解核心数据列时间、事件与协变量一个标准的生存分析数据集至少需要两列持续时间通常命名为T或duration。它表示从观察起点到发生“事件”或观察结束所经过的时间。时间的单位可以是天、月、年等但整个数据集必须统一。关键点这个时间必须是数值型。事件指示通常命名为E或event_observed。它是一个布尔值True/False或二进制数值1/0。1或True表示在T时刻我们关心的事件确实发生了。0或False表示在T时刻该个体被删失了——即观察期结束时事件仍未发生或者因其他原因如失访我们无法再追踪其结局。示例一个客户流失分析的数据集。客户A注册后第100天取消服务。则T100,E1。客户B注册后我们观察了200天直到分析日期他仍是活跃用户。则T200,E0。除了这两列其他所有你希望用来预测生存时间的变量都称为协变量比如用户的年龄、性别、消费金额、产品版本等。3.2 从现实时间数据转换一个典型场景现实中我们的数据往往记录的是具体的开始日期和结束日期或删失日期。lifelines.utils中的datetimes_to_durations函数是处理这类数据的利器。假设你的原始数据框df有以下列start_date: 观察开始日期如用户注册日、患者入院日end_date: 事件发生日期如流失日、死亡日如果未发生事件则为分析截止日期或最后联系日期。event_occurred: 是否发生事件是/否。from lifelines.utils import datetimes_to_durations # 将日期列转换为 pandas 的 datetime 类型 df[start_date] pd.to_datetime(df[start_date]) df[end_date] pd.to_datetime(df[end_date]) # 使用函数进行转换 # 参数说明 # start_times: 开始时间列 # end_times: 结束时间列 # event_observed: 事件是否发生的布尔列 T, E datetimes_to_durations( start_timesdf[start_date], end_timesdf[end_date], event_observeddf[event_occurred].map({是: True, 否: False}) # 转换为布尔值 ) # 将计算出的生存时间和事件指示符添加回原数据框 df[duration] T df[observed] E.astype(int) # 转换为0/1便于查看 print(df[[start_date, end_date, event_occurred, duration, observed]].head())重要提示datetimes_to_durations默认计算的是两个日期之间的天数差。如果你的业务逻辑是以月或年为单位需要手动转换。例如如果你想以“月”为单位可以df[‘duration_months’] T / 30.44近似值但更严谨的做法是使用pandas的DateOffset或直接计算月份差。务必确保时间单位的业务意义明确。3.3 数据清洗与探索性分析在建模前必须对数据进行清洗和探索。处理异常值和缺失值持续时间 0检查是否有duration 0的记录。这通常意味着数据录入错误例如结束日期早于开始日期。这类数据需要调查并修正或删除。协变量缺失lifelines的模型无法处理缺失值。你需要决定是删除缺失行、用中位数/众数填充还是使用更复杂的插补方法。对于关键协变量缺失过多可能要考虑放弃该变量。初步探索计算总体事件率df[‘observed’].mean()。这让你对数据的删失比例有个大致了解。绘制持续时间的直方图分别对事件组observed1和删失组observed0绘制看分布是否有明显差异。分类变量分布查看各个协变量在不同事件状态下的分布这能提供最初的线索。# 示例基本的数据探查 print(f总样本数: {len(df)}) print(f事件发生数: {df[observed].sum()}) print(f事件发生率: {df[observed].mean():.2%}) print(f中位生存时间天: {df[duration].median()}) # 查看分类变量与事件的关系 if group in df.columns: event_rate_by_group df.groupby(group)[observed].mean() print(\n各分组事件发生率:) print(event_rate_by_group)这个阶段花的时间越多后续建模就越顺畅。数据质量直接决定了模型结果的可信度。4. 生存函数的非参数估计Kaplan-Meier曲线实战当我们还没有引入任何协变量只想了解整体的生存模式或者想比较几个亚组例如不同治疗方案的生存差异时Kaplan-Meier估计器是我们的首选工具。它提供了一种直观、无偏的方式来估计生存函数。4.1 拟合与绘制整体生存曲线使用KaplanMeierFitter非常简单。# 初始化拟合器 kmf KaplanMeierFitter() # 拟合数据 # 参数duration时间 event_observed事件 kmf.fit(durationsdf[duration], event_observeddf[observed]) # 绘制生存曲线 fig, ax plt.subplots(figsize(10, 6)) kmf.plot_survival_function(axax) ax.set_title(整体生存曲线 (Kaplan-Meier)) ax.set_xlabel(时间 (天)) ax.set_ylabel(生存概率) ax.grid(True, linestyle--, alpha0.7) plt.show()这张图会显示一条从1.0100%存活开始随时间下降的曲线。曲线上的每一个“台阶”代表一个事件的发生。如何解读在时间点t曲线对应的纵坐标值S(t)表示一个个体从开始到时间t仍然“存活”未发生事件的概率估计。你还可以直接从拟合器中获取关键统计量# 计算中位生存时间 median_survival_time kmf.median_survival_time_ print(f中位生存时间: {median_survival_time:.2f} 天) # 查询特定时间点的生存概率 time_points [30, 90, 180] # 第3090180天 survival_probabilities kmf.predict(time_points) for t, prob in zip(time_points, survival_probabilities): print(f在时间 {t} 天生存概率约为 {prob:.2%})注意事项kmf.median_survival_time_返回的是生存概率下降到50%时对应的时间。如果曲线从未下降到0.5以下这个值会是inf无穷大。此时你可以报告其他分位数比如kmf.quantile(0.75)表示生存概率为75%的时间点。4.2 分组比较与Log-Rank检验生存分析中最常见的问题之一是“A组和B组的生存情况有显著差异吗” 我们可以通过为不同组分别拟合KM曲线并对其进行可视化比较和统计检验来回答。# 假设数据中有一个‘treatment’列值为‘Drug_A’或‘Drug_B’ ax plt.subplots(figsize(10, 6))[1] # 为每个治疗组分别拟合并绘制 groups df[treatment].unique() for name, grouped_df in df.groupby(treatment): kmf KaplanMeierFitter() kmf.fit(durationsgrouped_df[duration], event_observedgrouped_df[observed], labelname) kmf.plot_survival_function(axax) ax.set_title(不同治疗组的生存曲线比较) ax.set_xlabel(时间 (天)) ax.set_ylabel(生存概率) ax.legend() ax.grid(True, linestyle--, alpha0.7) plt.show()从图上可以直观看出哪条曲线“更高”生存情况更好。但视觉差异不一定代表统计显著。我们需要进行Log-Rank检验这是比较生存曲线最常用的非参数方法。from lifelines.statistics import logrank_test # 将数据按组分开 group_A df[df[treatment] Drug_A] group_B df[df[treatment] Drug_B] # 执行Log-Rank检验 results logrank_test( durations_Agroup_A[duration], durations_Bgroup_B[duration], event_observed_Agroup_A[observed], event_observed_Bgroup_B[observed] ) print(fLog-Rank检验结果:) print(f 检验统计量: {results.test_statistic:.4f}) print(f P值: {results.p_value:.6f}) if results.p_value 0.05: print( - P值 0.05拒绝原假设认为两组生存分布存在显著差异。) else: print( - P值 0.05没有足够证据认为两组生存分布存在显著差异。)实操心得Log-Rank检验的零假设是“两组生存分布相同”。一个很小的p值通常0.05意味着我们有足够证据认为两组生存情况不同。但要注意Log-Rank检验对全时间段的整体差异敏感如果两组曲线早期差异大但后期交叉检验效能可能会降低。此时需要结合图形和专业判断。4.3 为曲线添加置信区间与风险表专业的生存曲线图通常会包含置信区间以展示估计的不确定性以及底部的风险表显示在每个时间点仍处于风险中的个体数。fig, (ax1, ax2) plt.subplots(2, 1, figsize(12, 10), gridspec_kw{height_ratios: [3, 1]}, sharexTrue) # 在上方子图绘制带置信区间的生存曲线 kmf.fit(df[duration], df[observed], labelAll Patients) kmf.plot_survival_function(axax1, ci_showTrue) # ci_showTrue 显示置信区间 ax1.set_title(生存曲线 with 95% Confidence Interval) ax1.set_ylabel(生存概率) ax1.legend() ax1.grid(True, linestyle--, alpha0.7) # 在下方子图绘制风险表 # lifelines没有内置风险表绘图但我们可以手动计算并绘制 # 这里使用一个简化版本展示在几个关键时间点的风险人数 at_risk_times [0, 30, 90, 180, 365] risk_counts [] for t in at_risk_times: # 计算在时间t仍处于风险中的个体数即duration t count (df[duration] t).sum() risk_counts.append(count) ax2.bar(at_risk_times, risk_counts, width20, alpha0.7, colorskyblue) ax2.set_xlabel(时间 (天)) ax2.set_ylabel(风险人数) ax2.set_title(风险表 (Number at Risk)) ax2.set_xticks(at_risk_times) ax2.grid(True, axisy, linestyle--, alpha0.7) plt.tight_layout() plt.show()添加这些元素能使你的分析图表更具专业性和说服力尤其是在学术或商业报告场景中。5. 半参数建模核心Cox比例风险模型详解Kaplan-Meier曲线很棒但它只能描述不能解释。当我们想了解哪些因素会影响生存时间并且量化这些影响时就需要回归模型。Cox比例风险模型是生存分析回归的黄金标准。5.1 模型原理与假设的通俗理解你可以把Cox模型想象成一个“风险预测器”。它对每个个体在时间t的风险率hazard rate建模。风险率可以粗略理解为在活到时间t的条件下在接下来一个极短的单位时间内发生事件的概率。Cox模型的核心公式是h(t|X) h0(t) * exp(b1*X1 b2*X2 ... bp*Xp)h(t|X)给定协变量X的个体在时间t的风险率。h0(t)基线风险函数。它描述了所有协变量都为0时的风险随时间的变化。Cox模型的巧妙之处在于它不关心h0(t)的具体形式因此是“半参数”只关心协变量的影响。exp(b1*X1 ...)这部分量化了协变量对风险的影响。b是模型的系数需要我们从数据中估计。关键概念风险比对于二分类变量如性别男1女0系数b的意义是在其他变量不变的情况下男性相对于女性其风险比为exp(b)。如果exp(b) 1说明该变量是风险因素会增加事件发生的风险。如果exp(b) 1说明该变量是保护因素会降低风险。如果exp(b) 1说明该变量无影响。Cox模型有一个核心假设比例风险假设。它假定任意两个个体的风险比是恒定的不随时间改变。例如如果男性的风险是女性的2倍那么这个“2倍”的关系在整个研究期间都应该保持不变。验证这个假设至关重要我们稍后会讲到。5.2 使用lifelines拟合你的第一个Cox模型假设我们有一个数据集包含患者的生存时间、事件指示以及年龄、癌症分期、治疗方案等协变量。# 准备数据选择用于建模的列 # CoxPHFitter要求数据中必须包含持续时间列和事件指示列以及其他协变量。 model_df df[[duration, observed, age, stage, treatment]].copy() # 处理分类变量Cox模型需要数值输入。 # 对于有序分类如stage: I, II, III可以映射为数值1,2,3或进行独热编码。 # 对于无序分类如treatment: A, B必须进行独热编码。 model_df pd.get_dummies(model_df, columns[stage, treatment], drop_firstTrue, dtypeint) # drop_firstTrue 是为了避免多重共线性。例如对于stage如果创建了stage_II和stage_III # 那么stage_I就被作为了参考基线。 print(建模数据前5行:) print(model_df.head()) print(f\n数据形状: {model_df.shape})现在我们可以拟合Cox模型了。# 初始化CoxPHFitter cph CoxPHFitter() # 拟合模型 # 参数dataframe duration_col event_col cph.fit(model_df, duration_colduration, event_colobserved, show_progressTrue) # 打印模型摘要 cph.print_summary()print_summary()会输出一个非常详细的表格包含coef: 每个协变量的系数b。exp(coef): 风险比HR exp(b)。这是最重要的解读指标。se(coef): 系数的标准误。z: z统计量系数除以标准误。p: p值。用于检验该系数是否显著不为0。lower 0.95 / upper 0.95: 风险比95%的置信区间。如何解读 以age变量为例假设其exp(coef) 1.05p 0.05。exp(coef)1.05意味着年龄每增加一岁死亡或事件风险是原来的1.05倍即风险增加了5%。95% CI: 如果置信区间不包含1例如 [1.02, 1.08]进一步支持该效应是显著的。p 0.05表明年龄的影响具有统计学意义。5.3 模型诊断比例风险假设检验如果比例风险假设不成立Cox模型的估计可能是有偏的。lifelines提供了便捷的检验方法。# 执行比例风险假设检验 from lifelines.statistics import proportional_hazard_test # proportional_hazard_test 需要模型拟合后的结果和原始数据 results proportional_hazard_test(cph, model_df, time_transformrank) # time_transform通常用rank print(results.summary)检验结果会为每个变量提供一个p值。原假设是“满足比例风险假设”。因此如果某个变量的p值很小如0.05则拒绝原假设认为该变量违反了比例风险假设。如果p值较大则没有证据拒绝原假设可以认为该变量满足假设。如果假设被违反了怎么办分层对于违反假设的分类变量可以将其作为分层变量。这意味着模型为这个变量的每个水平估计一个不同的基线风险函数但协变量的系数风险比在各层间保持一致。在lifelines中使用strata参数cph.fit(..., strata[stage])。加入时间交互项对于连续变量可以将其与时间或时间的函数进行交互允许其效应随时间变化。例如cph.fit(model_df.assign(age_time df[‘age’] * df[‘duration’]), …)但这会使得模型解释变得复杂。考虑参数模型如果主要变量严重违反假设可以考虑使用参数AFT模型如Weibull模型它们没有比例风险假设的限制。5.4 模型预测与可视化拟合好的Cox模型可以用来进行预测。预测部分风险比给定一个新个体的特征我们可以预测他的风险评分即线性部分b1*X1…并与其他个体比较相对风险。# 假设一个新患者年龄50岁 stage_III1 (stage III期) treatment_B1 (使用B疗法) new_patient pd.DataFrame.from_dict({ age: [50], stage_II: [0], stage_III: [1], treatment_B: [1] }) # 注意DataFrame的列名和顺序必须与训练数据完全一致 # 对于在get_dummies中被drop的基线类别如stage_I, treatment_A我们赋值为0。 # 预测风险评分即线性预测值 partial_hazard cph.predict_partial_hazard(new_patient) print(f该患者的风险评分相对于基线: {partial_hazard.values[0]:.4f}) # 这个值本身没有绝对意义但可以用于比较。值越大风险越高。预测生存函数我们还可以预测这个新个体在不同时间点的生存概率。# 预测该患者在不同时间点的生存概率 survival_function cph.predict_survival_function(new_patient) # survival_function 是一个DataFrame索引是时间点列是预测的生存概率 # 绘制预测的生存曲线 ax survival_function.plot(figsize(10,6)) ax.set_title(新患者的预测生存曲线) ax.set_xlabel(时间 (天)) ax.set_ylabel(生存概率) ax.grid(True, linestyle--, alpha0.7) plt.show() # 获取中位生存时间预测 # 注意如果预测的生存曲线从未降到0.5以下中位生存时间会是无穷大(inf) predicted_median cph.predict_median(new_patient) print(f预测的中位生存时间: {predicted_median.values[0]:.2f} 天)可视化协变量效应我们可以固定其他变量观察某个连续变量如年龄变化对生存曲线的影响。# 创建一个虚拟数据集固定其他变量让年龄从40岁变化到70岁 # 假设参考患者stage_II0, stage_III0 (即stage I), treatment_B0 (即treatment A) ages_to_plot np.linspace(40, 70, 4) # 40, 50, 60, 70岁 plotting_df pd.DataFrame({ age: ages_to_plot, stage_II: 0, stage_III: 0, treatment_B: 0 }) ax plt.subplots(figsize(10,6))[1] for i, row in plotting_df.iterrows(): sf cph.predict_survival_function(row.to_frame().T) sf.T.plot(axax, labelfAge{int(row[age])}) ax.set_title(不同年龄患者的预测生存曲线 (其他因素固定)) ax.set_xlabel(时间 (天)) ax.set_ylabel(生存概率) ax.legend(title年龄) ax.grid(True, linestyle--, alpha0.7) plt.show()这种可视化能非常直观地展示关键风险因素的影响力。6. 进阶技巧与实战问题排查掌握了基础操作后在实际项目中你会遇到更复杂的情况。下面分享一些进阶技巧和常见问题的解决方法。6.1 处理连续变量的非线性关系样条函数在Cox模型中我们默认连续变量如年龄对风险的影响是线性的即log(HR)与年龄呈线性关系。但现实中这种关系可能是非线性的。例如年龄对风险的影响可能在中年和老年阶段不同。我们可以使用限制性立方样条来捕捉这种非线性。lifelines支持通过patsy公式来指定样条。首先需要安装patsy和scipy。# 使用patsy公式指定包含年龄样条的模型 from lifelines import CoxPHFitter import patsy # 为年龄创建3个节点的限制性立方样条项 # ‘cr’代表自然立方样条 model_df_spline model_df.copy() # 使用patsy.dmatrix创建样条基 age_spline_basis patsy.dmatrix(cr(age, df3), model_df_spline, return_typedataframe) # 将生成的样条列加入数据框并移除原始的年龄列 model_df_spline pd.concat([model_df_spline.drop(columns[age]), age_spline_basis], axis1) # 拟合包含样条的Cox模型 cph_spline CoxPHFitter() cph_spline.fit(model_df_spline, duration_colduration, event_colobserved) cph_spline.print_summary()现在模型不再输出一个age的系数而是输出多个样条基的系数。要解读年龄的影响最好通过绘图来观察预测的风险比随年龄变化的曲线。6.2 模型性能评估一致性指数对于分类或回归模型我们有准确率、R平方等指标。对于生存模型最常用的评估指标是一致性指数也称为C-index或C-statistic。它的含义类似于ROC曲线下的面积AUC它衡量的是模型预测的“风险排序”与实际观察到的生存时间排序之间的一致性。C-index为0.5表示模型没有预测能力随机猜测为1表示完美预测。# 计算训练数据上的C-index c_index_train cph.concordance_index_ print(f训练集C-index: {c_index_train:.4f}) # 如果你有独立的测试集可以这样计算 # 假设 test_df 是处理好的测试集数据 # c_index_test concordance_index(test_df[duration], -cph.predict_partial_hazard(test_df), test_df[observed]) # 注意concordance_index函数期望风险评分越高生存时间越短所以对预测的风险评分取负。一般来说在生物医学领域C-index在0.7以上可以认为模型有不错的区分能力。6.3 常见报错与排查技巧ConvergenceError: Convergence halted due to matrix inversion problems...原因这是最常见的问题之一。通常意味着数据存在完全分离或多重共线性。例如某个分类变量的某个类别下所有个体要么都发生了事件要么都没发生。排查检查分类变量的分布df[‘your_variable’].value_counts()。检查事件在不同类别中的分布pd.crosstab(df[‘your_variable’], df[‘observed’])。如果某个类别样本量极少或事件完全一致考虑合并类别或删除该变量。检查连续变量是否有异常大的值尝试标准化。解决使用penalizer参数添加一个微小的L2正则化这可以稳定矩阵求逆过程。cph CoxPHFitter(penalizer0.01) # 尝试一个小的惩罚系数如0.01或0.1 cph.fit(...)ValueError: The following columns are constant...原因数据中存在方差为0的常数列。解决在拟合前移除这些列model_df model_df.loc[:, model_df.std() 0]生存曲线不下降或中位生存时间为inf原因事件发生率太低或者大量删失发生得很早导致估计的生存概率在很长一段时间内都无法下降到0.5以下。解读这本身是一个有意义的发现说明群体生存情况很好。在报告中可以报告其他分位数如75%生存时间或平均生存时间需注意平均生存时间在大量删失下可能估计不准。预测新数据时报错提示列名不匹配原因pandas.get_dummies在训练集和测试集上可能生成不同的列如果某个类别只在其中一个集合中出现。解决在数据处理阶段先定义一个所有可能类别的完整集合然后使用pd.get_dummies(..., columns..., drop_firstTrue).reindex(columnspredefined_columns, fill_value0)来确保训练和测试数据框的列结构完全一致。6.4 一个完整的工作流示例最后让我们用一个模拟的、更接近真实场景的数据集串联起从数据加载到模型评估的完整流程。# 1. 加载lifelines内置数据集示例 from lifelines.datasets import load_rossi rossi load_rossi() print(rossi.head()) print(f事件率: {rossi[arrest].mean():.2%}) # 2. 数据准备此数据集已基本处理好 df rossi.copy() # 假设我们想研究‘fin’是否有经济资助和‘age’年龄对‘arrest’再逮捕的影响 # ‘week’是生存时间 # 3. 探索性分析按‘fin’分组绘制KM曲线 kmf_fin KaplanMeierFitter() ax plt.subplots(figsize(10,6))[1] for name, grouped_df in df.groupby(fin): kmf_fin.fit(grouped_df[week], grouped_df[arrest], labelffin{name}) kmf_fin.plot_survival_function(axax) ax.set_title(KM Curve by Financial Aid (fin)) plt.show() # 4. 拟合Cox模型 cph CoxPHFitter() cph.fit(df, duration_colweek, event_colarrest) cph.print_summary() # 5. 模型诊断 proportional_hazard_test(cph, df).print_summary() # 6. 预测示例一个没有经济资助的25岁个体 new_individual pd.DataFrame.from_dict({ fin: [0], age: [25], race: [1], wexp: [0], mar: [0], paro: [0], prio: [2] }) # 注意必须包含数据集中所有的协变量列 print(f\n预测风险评分: {cph.predict_partial_hazard(new_individual).values[0]:.3f}) print(f预测的1年内52周不被逮捕的概率: {cph.predict_survival_function(new_individual).loc[52].values[0]:.2%}) # 7. 评估 print(f\n模型C-index: {cph.concordance_index_:.4f})通过这个完整的流程你应该已经能够使用lifelines库独立地处理和分析自己的生存数据了。记住生存分析不仅仅是运行代码更重要的是理解数据背后的业务逻辑、合理解读统计结果并将发现用清晰的语言和可视化呈现出来。从Kaplan-Meier曲线到Cox模型你手中的工具已经可以解决大部分常见的生存分析问题。