Python生存分析实战:从Kaplan-Meier到Cox模型,用lifelines处理时间-事件数据
2026/8/23 4:19:02 网站建设 项目流程

1. 项目概述:为什么是生存分析与lifelines?

在数据科学和统计分析领域,我们常常会遇到一类特殊的时间-事件数据。比如,你想知道一款新药上市后,患者从开始服药到病情复发平均需要多长时间;或者,你想评估一个营销活动后,客户从注册到流失的周期。这类问题的核心,不仅仅是事件“是否”发生,更重要的是事件“何时”发生。传统的分类模型(如逻辑回归)只能告诉你“会不会流失”,而生存分析(Survival Analysis)则能更进一步,回答“大概多久后会流失”以及“不同因素如何影响这个时间”。

生存分析起源于医学和工程可靠性领域,用于研究从某个起点(如诊断、产品出厂)到某个特定终点事件(如死亡、故障)发生的时间。它的强大之处在于能优雅地处理“删失数据”——那些在我们观察期结束时,终点事件尚未发生的数据。例如,在研究结束时,有些患者依然存活,有些客户尚未流失,这些数据并非无效,它们提供了“至少存活了这么久”的宝贵信息。忽视它们会导致严重的估计偏差。

对于Python用户而言,lifelines库无疑是进入生存分析世界最友好、功能最全面的门票。它不像R语言的survival包那样有着悠久的历史包袱,而是为Python的数据科学生态量身定制,与pandasnumpyscikit-learn等库无缝集成。lifelines封装了从非参数估计(如Kaplan-Meier曲线)到半参数模型(Cox比例风险模型),再到参数模型(如Weibull、Log-Normal)等一系列方法,并且提供了清晰、一致的API和优秀的可视化支持。

如果你手头有包含时间戳和事件状态的数据,并且想知道时间背后的故事,那么生存分析和lifelines就是你不可或缺的工具。无论你是医疗数据分析师、金融风控专家,还是用户增长策略师,掌握它都能让你对“时间”这个维度有更深刻、更量化的理解。

2. 核心概念与数据准备:理解你的“生存数据”

在动手写代码之前,我们必须把生存分析的核心概念和数据结构理清楚。这就像盖房子前先看明白图纸,能避免后续很多“返工”的坑。

2.1 生存分析的三要素

任何生存分析数据集,无论来自哪个领域,通常都包含以下三个核心要素:

  1. 生存时间(Duration/T):从起始时间点到终点事件发生,或者到观察截止所经过的时间。单位可以是天、月、年等。关键点:这个时间必须是连续或近似连续的。例如,客户的生命周期(天)、设备的无故障运行时间(小时)、疾病复发时间(月)。

  2. 事件状态(Event/E):一个指示终点事件是否发生的标志。通常用1表示事件发生(如死亡、流失、故障),用0表示删失(如研究结束时仍存活、客户仍活跃、设备仍在运行)。这是生存分析处理“不完全信息”的核心

  3. 协变量(X):可能影响生存时间的特征或变量。例如,患者的年龄、治疗方案;客户的性别、消费等级;设备的运行环境、生产批次等。我们的模型目标就是量化这些协变量如何影响生存概率。

2.2 用pandas准备你的数据

lifelinespandasDataFrame是绝配。你的数据应该被组织成一个DataFrame,其中至少有两列分别对应生存时间和事件状态。

假设我们有一个模拟的客户流失数据集:

import pandas as pd import numpy as np from lifelines.datasets import load_rossi # lifelines自带的一个经典数据集 # 加载一个示例数据集(关于罪犯再逮捕) rossi = load_rossi() print(rossi.head()) print(f"\n数据列名: {rossi.columns.tolist()}") print(f"‘week’列是生存时间, ‘arrest’列是事件状态(1=再被捕)")

对于你自己的数据,准备流程通常如下:

# 假设你有一个原始的客户订单/活跃日志表 `df_raw` # 1. 定义观察起点和终点 # 起点:客户首次购买日期 # 终点事件:客户超过90天无任何交互(视为流失) # 观察截止日期:2023-12-31 df_raw['first_purchase_date'] = pd.to_datetime(df_raw['first_purchase_date']) df_raw['last_interaction_date'] = pd.to_datetime(df_raw['last_interaction_date']) analysis_cutoff_date = pd.Timestamp('2023-12-31') # 2. 计算生存时间(单位:天) def calculate_duration(row): # 如果最后交互日期后90天超过了截止日期,则为删失(事件未发生) churn_deadline = row['last_interaction_date'] + pd.Timedelta(days=90) if churn_deadline > analysis_cutoff_date: # 删失:生存时间 = 截止日期 - 首次购买日期 duration = (analysis_cutoff_date - row['first_purchase_date']).days event = 0 else: # 事件发生:生存时间 = (最后交互日期+90天) - 首次购买日期 duration = (churn_deadline - row['first_purchase_date']).days event = 1 return pd.Series([duration, event]) df_survival[['T', 'E']] = df_raw.apply(calculate_duration, axis=1) # 3. 添加协变量 df_survival['age'] = df_raw['age'] df_survival['subscription_tier'] = df_raw['tier'] # 例如,1=基础,2=高级 df_survival['avg_monthly_spend'] = df_raw['monthly_spend'] # 检查数据 print(df_survival[['T', 'E', 'age', 'subscription_tier']].head()) print(f"\n事件发生率: {df_survival['E'].mean():.2%}") print(f"平均生存时间(含删失): {df_survival['T'].mean():.1f} 天")

注意:生存时间T必须大于0。如果你的数据中出现了0或负数,需要检查起点和终点的定义逻辑。常见的坑是起点和终点定义在同一天,导致生存时间为0。

2.3 数据质量的快速检查

在建模前,花几分钟做以下检查能省去大量调试时间:

  • 查看分布:用df_survival[‘T’].describe()df_survival[‘T’].hist()看看生存时间的分布是否合理,有无异常值(比如超过理论最大值的记录)。
  • 检查删失比例:事件率(E的均值)不宜过低(如<5%)或过高(如>95%),否则模型可能不稳定或缺乏信息量。通常20%-80%是比较理想的范围。
  • 协变量与时间的相关性:简单画个散点图,看看关键协变量(如avg_monthly_spend)与生存时间T是否有直观关系。这有助于后续理解模型结果。
  • 检查唯一性:确保没有重复的个体ID(如果你的数据是个人层面的)。

3. 生存函数的非参数估计:Kaplan-Meier与组间比较

当我们暂时不关心“为什么”,只想先看看“整体存活情况如何”或者“A组和B组谁的生存率更高”时,非参数方法是最佳起点。它不对生存时间的分布做任何假设,直接根据数据计算生存概率。

3.1 绘制整体的Kaplan-Meier生存曲线

Kaplan-Meier估计器是生存分析的基石。它给出了在任意时间点t,个体“存活”(即事件尚未发生)的概率估计值S(t)

from lifelines import KaplanMeierFitter import matplotlib.pyplot as plt # 初始化拟合器 kmf = KaplanMeierFitter() # 传入生存时间和事件状态进行拟合 kmf.fit(durations=df_survival['T'], event_observed=df_survival['E']) # 绘制生存曲线 fig, ax = plt.subplots(figsize=(10, 6)) kmf.plot_survival_function(ax=ax) ax.set_title('整体客户生存曲线(Kaplan-Meier估计)') ax.set_xlabel('生存时间(天)') ax.set_ylabel('生存概率 S(t)') ax.grid(True, linestyle='--', alpha=0.7) # 可以在曲线上标注中位生存时间 median_survival_time = kmf.median_survival_time_ if median_survival_time is not np.inf: ax.axvline(x=median_survival_time, color='red', linestyle=':', alpha=0.8, label=f'中位生存时间: {median_survival_time:.1f}天') ax.legend() plt.tight_layout() plt.show() # 打印关键摘要信息 print(f"中位生存时间: {median_survival_time}") print(f"在时间点 t=180 天的生存概率: {kmf.survival_function_at_times(180).iloc[0]:.3f}") kmf.print_summary()

print_summary()会输出一个详细的表格,包括估计的生存概率在不同时间点的置信区间,非常实用。

3.2 分组比较:Log-Rank检验

业务中更常见的问题是:“付费用户和免费用户的留存率有显著差异吗?” 这时我们需要比较不同组的生存曲线。lifelines提供了便捷的接口进行分组KM估计和Log-Rank检验(一种检验多条生存曲线是否相同的非参数方法)。

from lifelines.statistics import logrank_test # 假设我们根据‘subscription_tier’分组 tiers = df_survival['subscription_tier'].unique() tiers.sort() fig, ax = plt.subplots(figsize=(10, 6)) for tier in tiers: mask = df_survival['subscription_tier'] == tier kmf_tier = KaplanMeierFitter(label=f'订阅等级 {tier}') kmf_tier.fit(durations=df_survival.loc[mask, 'T'], event_observed=df_survival.loc[mask, 'E']) kmf_tier.plot_survival_function(ax=ax) ax.set_title('不同订阅等级的客户生存曲线对比') ax.set_xlabel('生存时间(天)') ax.set_ylabel('生存概率 S(t)') ax.grid(True, linestyle='--', alpha=0.7) plt.tight_layout() plt.show() # 执行Log-Rank检验(例如,比较等级1和等级2) mask_tier1 = df_survival['subscription_tier'] == 1 mask_tier2 = df_survival['subscription_tier'] == 2 results = logrank_test( df_survival.loc[mask_tier1, 'T'], df_survival.loc[mask_tier2, 'T'], event_observed_A=df_survival.loc[mask_tier1, 'E'], event_observed_B=df_survival.loc[mask_tier2, 'E'], labels=['Tier 1', 'Tier 2'] ) results.print_summary()

Log-Rank检验的结果会给出一个p值。如果p值小于0.05(或你设定的显著性水平),我们就可以拒绝“两组生存曲线相同”的原假设,认为两组的生存时间分布存在统计学上的显著差异。

实操心得:KM曲线和Log-Rank检验是向非技术背景同事(如产品经理、业务方)展示结果的神器。一张图加上“两组差异显著(p<0.01)”的结论,比任何复杂的模型系数都更有说服力。但要注意,KM曲线只能展示单一分组变量的影响,无法同时控制其他因素。

4. 半参数模型:Cox比例风险模型入门

Kaplan-Meier曲线告诉我们“是什么”,而Cox比例风险模型(Cox Proportional Hazards Model)则试图解释“为什么”。它是生存分析中最常用、最核心的回归模型,属于半参数模型——它对风险函数的基础形状不做假设(非参数部分),但假设协变量对风险的影响是乘性的且不随时间改变(参数部分)。

4.1 模型拟合与解读

“风险”(Hazard)可以理解为在某一时间点,个体瞬间发生事件的概率。Cox模型的形式是:h(t|X) = h0(t) * exp(β1*X1 + β2*X2 + ...)其中h0(t)是基准风险函数(未知),exp(βi)就是风险比(Hazard Ratio, HR)。

from lifelines import CoxPHFitter # 1. 初始化并拟合模型 cph = CoxPHFitter() # 确保数据中'T'和'E'列名正确,其他列将被视为协变量 cph.fit(df_survival, duration_col='T', event_col='E') # 2. 查看模型摘要 cph.print_summary()

print_summary()的输出是理解模型的关键,我们逐列解读:

  • coef: 系数 β。正值表示该变量增加风险(缩短生存时间),负值表示降低风险(延长生存时间)。
  • exp(coef): 风险比 HR = exp(β)。这是更直观的指标。
    • HR > 1: 该变量是风险因素。例如,HR=1.5表示拥有该特征(或每增加一单位)的个体,其事件发生的风险是参照组的1.5倍。
    • HR < 1: 该变量是保护因素。例如,HR=0.7表示风险降低到参照组的70%。
  • se(coef): 系数的标准误,衡量估计的精度。
  • z: 检验统计量 z = coef / se(coef),用于计算p值。
  • p: p值。通常认为p<0.05时,该协变量对风险有显著影响。
  • -log2(p): p值的另一种表示,值越大越显著。

4.2 模型假设检验:比例风险假设

Cox模型最重要的前提假设是“比例风险”(Proportional Hazards, PH),即任意两个个体的风险比不随时间改变。如果假设不成立,模型的估计可能是有偏的。lifelines提供了两种主要方法来检验PH假设:

# 方法1:基于Schoenfeld残差的统计检验 cph.check_assumptions(df_survival, p_value_threshold=0.05, show_plots=True)

运行check_assumptions会输出一个详细的报告。对于每个协变量,它会进行统计检验并给出p值。如果某个变量的p值小于阈值(如0.05),则提示该变量可能违反了PH假设。同时,它会绘制Schoenfeld残差随时间变化的图,理想情况下,残差应该随机分布在0附近,没有明显的趋势。

如果PH假设被违反怎么办?

  1. 分层(Stratification):对于分类变量,可以在拟合时指定其为分层变量。这意味着模型为这个变量的每个层级估计一个不同的基准风险函数h0(t),但协变量的系数β在不同层间保持一致。这适用于重要的分类预测因子(如肿瘤分期)但不关心其具体HR值的情况。
    cph_strat = CoxPHFitter() cph_strat.fit(df_survival, duration_col='T', event_col='E', strata=['subscription_tier'])
  2. 引入时间交互项:如果连续变量(如年龄)违反PH假设,可以将其与时间(或时间的函数)交互,允许其效应随时间变化。这增加了模型复杂度。
    # 在数据中创建一个与时间交互的项(例如,年龄*log(t)) # 这通常需要更深入的理论指导
  3. 考虑参数模型或加速失效时间模型:如果主要变量严重违反PH假设,可以考虑放弃Cox模型,转而使用参数模型(如Weibull AFT模型)。

注意事项:在实践中,轻微的PH假设违反有时是可以接受的,尤其是当样本量很大时。check_assumptions是一个重要的诊断工具,但不应成为教条。需要结合统计结果和业务意义综合判断。

4.3 可视化模型结果

让模型结果更直观:

# 1. 绘制系数森林图(Forest Plot) fig, ax = plt.subplots(figsize=(8, len(cph.params)/2)) cph.plot(ax=ax) ax.set_title('Cox模型风险比森林图') plt.tight_layout() plt.show() # 森林图直观展示了每个变量的风险比(HR)及其置信区间。如果区间包含1,则效应不显著。 # 2. 绘制部分依赖图:展示某个变量对生存概率的影响 fig, axes = plt.subplots(1, 2, figsize=(14, 5)) # 假设我们关注‘avg_monthly_spend’和‘age’ cph.plot_partial_effects_on_outcome(covariates='avg_monthly_spend', values=[0, 50, 100, 150], ax=axes[0]) axes[0].set_title('月均消费对生存概率的影响') axes[0].set_ylabel('生存概率 S(t)') cph.plot_partial_effects_on_outcome(covariates='age', values=[20, 40, 60, 80], ax=axes[1]) axes[1].set_title('年龄对生存概率的影响') axes[1].set_ylabel('') plt.tight_layout() plt.show() # 这个图显示了在固定其他变量为均值(或指定值)时,改变某个变量,生存曲线如何变化。 # 3. 预测个体风险 # 我们可以预测新个体的风险评分或生存概率 # 风险评分(越高风险越大) partial_hazard = cph.predict_partial_hazard(df_survival.iloc[0:5]) # 前5个样本 print("部分风险(线性预测值):\n", partial_hazard) # 预测在特定时间点的生存概率 survival_prob = cph.predict_survival_function(df_survival.iloc[0:5], times=[30, 90, 180]) print("\n在30, 90, 180天的生存概率预测:\n", survival_prob)

5. 参数生存模型:当我们需要一个明确的分布形式

Cox模型灵活,但有时我们需要对生存时间的分布本身进行建模和预测,或者PH假设严重不成立。这时,参数模型就派上用场了。参数模型假设生存时间服从某个特定的参数分布,如指数分布、威布尔分布、对数正态分布等。

5.1 常用参数模型简介

  • 指数模型(ExponentialFitter):最简单,假设风险函数是常数(不随时间变化)。适用于事件发生率为恒定的场景,现实中较少见。
  • 威布尔模型(WeibullFitter):非常灵活,风险函数可以随时间递增、递减或恒定。是工程可靠性分析中最常用的模型之一。
  • 对数正态模型(LogNormalFitter):假设生存时间的对数服从正态分布。适用于生存时间分布右偏(长尾)明显的场景,比如某些疾病的复发时间。
  • 对数逻辑斯蒂模型(LogLogisticFitter):风险函数呈单峰形状(先增后减),适用于有“早期风险高,然后降低”特征的数据。

5.2 拟合与模型选择

from lifelines import WeibullFitter, LogNormalFitter, LogLogisticFitter # 拟合多个参数模型 wf = WeibullFitter() lnf = LogNormalFitter() llf = LogLogisticFitter() wf.fit(df_survival['T'], event_observed=df_survival['E'], label='Weibull') lnf.fit(df_survival['T'], event_observed=df_survival['E'], label='LogNormal') llf.fit(df_survival['T'], event_observed=df_survival['E'], label='LogLogistic') # 绘制生存函数对比 fig, ax = plt.subplots(figsize=(10, 6)) wf.plot_survival_function(ax=ax) lnf.plot_survival_function(ax=ax) llf.plot_survival_function(ax=ax) kmf.plot_survival_function(ax=ax, label='Kaplan-Meier (Non-Parametric)') # 加入KM曲线作为基准 ax.set_title('不同参数模型与KM估计的生存曲线对比') ax.set_xlabel('生存时间(天)') ax.set_ylabel('生存概率 S(t)') ax.legend() plt.tight_layout() plt.show()

通过图形可以直观地看哪个参数模型的曲线最贴近非参数的KM曲线,越贴近通常说明拟合越好。

5.3 使用AIC进行模型选择

对于参数模型,我们可以使用赤池信息准则(AIC)来客观比较。AIC值越小,模型在拟合优度和复杂度之间的权衡越好。

print(f"Weibull AIC: {wf.AIC_:.2f}") print(f"LogNormal AIC: {lnf.AIC_:.2f}") print(f"LogLogistic AIC: {llf.AIC_:.2f}")

选择AIC最小的模型。注意,参数模型通常用于描述整体生存时间分布或进行简单的预测。当引入多个协变量时,需要使用对应的参数回归模型(如WeibullAFTFitter),它们类似于Cox模型,但直接对生存时间(或其对数)进行建模,解释的是时间加速/减速因子,而非风险比。

6. 模型评估与验证:你的模型可靠吗?

拟合模型只是第一步,评估其性能和泛化能力至关重要。生存模型的评估有其特殊性,因为数据包含删失。

6.1 区分度评估:Concordance Index (C-index)

C-index是生存模型最常用的区分度指标,可以理解为模型预测排序与实际结果一致的概率。对于任意两个可比较的个体,如果模型预测风险高的那个个体实际事件发生得更早,则认为预测是一致的。C-index的范围是0到1,0.5表示随机猜测,1表示完美预测。通常大于0.7认为模型有一定区分能力。

from lifelines.utils import concordance_index # 计算C-index (在训练集上,可能过于乐观) c_index_train = cph.score(df_survival, scoring_method='concordance_index') print(f"Cox模型在训练集上的C-index: {c_index_train:.4f}") # 更可靠的做法是使用交叉验证 from lifelines.utils import k_fold_cross_validation from lifelines import CoxPHFitter def cph_fitter(df, training_timeline, testing_timeline): cph = CoxPHFitter(penalizer=0.1) # 可以加入一点L2正则化防止过拟合 cph.fit(df, duration_col='T', event_col='E') return cph scores = k_fold_cross_validation(cph_fitter, df_survival, 'T', event_col='E', k=5, scoring_method='concordance_index') print(f"5折交叉验证C-index得分: {np.mean(scores):.4f} (+/- {np.std(scores):.4f})")

6.2 校准度评估:预测 vs 实际

校准度衡量模型预测的生存概率与实际观测到的生存概率是否一致。一个常用的方法是绘制校准曲线。

from lifelines.calibration import survival_probability_calibration # 将数据分为训练集和测试集(时间序列数据需注意避免未来信息泄露) from sklearn.model_selection import train_test_split df_train, df_test = train_test_split(df_survival, test_size=0.3, random_state=42) # 在训练集上拟合模型 cph_train = CoxPHFitter() cph_train.fit(df_train, duration_col='T', event_col='E') # 在测试集上评估校准度 fig, ax = plt.subplots(figsize=(8, 6)) survival_probability_calibration(cph_train, df_test, t0=180, ax=ax) # 评估在t=180天时的校准度 ax.set_title('模型在180天生存概率上的校准曲线 (测试集)') ax.set_xlabel('预测的生存概率') ax.set_ylabel('实际观测的生存概率') ax.plot([0, 1], [0, 1], 'k:', label='理想校准线') ax.legend() ax.grid(True, linestyle='--', alpha=0.7) plt.tight_layout() plt.show()

理想情况下,点应该沿着对角线分布。如果点系统地落在对角线之上,说明模型过于悲观(预测的生存概率低于实际);反之则过于乐观。

6.3 综合评估:时间相关的ROC曲线

对于生存数据,我们还可以计算不同时间点的ROC曲线和AUC,这比单一的C-index提供了更动态的评估。

from lifelines.utils import survival_events_from_table from sklearn.metrics import roc_auc_score import warnings warnings.filterwarnings('ignore') # 计算在特定时间点的风险评分(作为预测值) df_test['risk_score'] = cph_train.predict_partial_hazard(df_test) # 我们需要一个函数来计算在时间t的AUC def auc_at_time(t, df, risk_scores): # 在时间t,事件发生=1,删失或事件发生在t之后=0 y_true = (df['T'] <= t) & (df['E'] == 1) # 注意:对于在时间t之前删失的个体,我们不知道他们是否会在t之前发生事件,通常将其排除或视为阴性 # 这里简化处理,仅使用在时间t之前事件发生或确定存活(T>t)的样本 mask = (df['T'] > t) | y_true if mask.sum() < 10: # 样本太少则返回NaN return np.nan return roc_auc_score(y_true[mask], risk_scores[mask]) # 计算多个时间点的AUC times_to_evaluate = [30, 90, 180, 365] aucs = {} for t in times_to_evaluate: auc_val = auc_at_time(t, df_test, df_test['risk_score']) aucs[t] = auc_val print(f"在时间点 t={t} 天, AUC = {auc_val:.3f}")

7. 进阶技巧与避坑指南

掌握了基础建模和评估后,我们来看看一些实战中提升模型效果和效率的技巧,以及那些容易踩的坑。

7.1 处理连续变量与分类变量

  • 连续变量:在放入Cox模型前,考虑是否需要标准化或缩放。这不会改变模型的拟合优度,但可以使系数更容易比较。特别是当变量量纲差异很大时(如年龄和月收入),使用sklearnStandardScaler是个好习惯。
    from sklearn.preprocessing import StandardScaler scaler = StandardScaler() df_survival[['age_scaled', 'spend_scaled']] = scaler.fit_transform(df_survival[['age', 'avg_monthly_spend']]) # 然后用缩放后的变量建模
  • 分类变量lifelinespandas可以很好地处理通过pd.get_dummies生成的虚拟变量。但要注意避免虚拟变量陷阱(完全多重共线性)。通常我们设置drop_first=True
    df_with_dummies = pd.get_dummies(df_survival, columns=['city_tier'], drop_first=True, prefix='city')

7.2 特征工程与变量选择

  • 时间依赖型协变量:有些变量会随时间变化(如病人的血压、客户的累计消费)。标准的Cox模型需要这些变量在时间起点处的值。处理时变协变量需要更复杂的扩展Cox模型或使用lifelinesCoxTimeVaryingFitter,这对数据格式有特殊要求(每行代表一个个体在某个时间区间内的状态)。
  • 变量选择:与线性回归类似,可以使用前向选择、后向消除或LASSO(L1正则化)。CoxPHFitter支持penalizer参数进行L2正则化。对于更激进的变量选择,可以结合scikit-learnRFE(递归特征消除)或使用专门的生存分析特征选择包。

7.3 处理极端值与缺失值

  • 生存时间异常值:检查生存时间的最大值是否合理。有时由于数据错误,会出现极大的生存时间(如99999天)。需要根据业务逻辑设定一个合理的上限进行截断(Winsorizing)或视为删失。
  • 缺失值lifelines的模型不支持数据中有NaN。必须进行缺失值处理。对于协变量,可以使用中位数/众数填充、基于模型的插补,或者简单删除缺失率过高的变量/样本。需要谨慎评估删除样本是否引入偏差。

7.4 模型不稳定或收敛警告

如果拟合Cox模型时出现收敛警告或系数非常大/非常小,可能的原因和解决办法:

  1. 完全分离:某个预测变量完美区分了事件发生与否。需要检查数据或合并类别。
  2. 多重共线性:预测变量之间高度相关。检查方差膨胀因子(VIF),考虑删除或合并相关变量。
  3. 样本量不足:特别是事件数太少。经验法则是,每个待估计的参数(协变量)至少需要10-20个事件。如果事件数少,就只放入最重要的几个变量。
  4. 数值问题:尝试对连续变量进行标准化。

7.5 一个完整的建模示例流程

最后,我们用一个简化的伪代码串起整个流程:

# 1. 导入与数据准备 import pandas as pd, numpy as np, matplotlib.pyplot as plt from lifelines import CoxPHFitter, KaplanMeierFitter from lifelines.statistics import logrank_test from sklearn.model_selection import train_test_split from sklearn.preprocessing import StandardScaler # 加载和清洗数据 df = pd.read_csv('your_survival_data.csv') df['T'] = df['end_date'] - df['start_date'] # 计算生存时间 df['E'] = df['event_flag'] # 事件标志 df = df.dropna(subset=['T', 'E']) # 删除关键列缺失的样本 # 2. 探索性分析 print(df[['T', 'E']].describe()) print(f"事件率: {df['E'].mean():.2%}") # 绘制KM曲线 kmf = KaplanMeierFitter().fit(df['T'], df['E']) kmf.plot_survival_function() # 3. 特征工程 # 处理分类变量 df = pd.get_dummies(df, columns=['category_var'], drop_first=True) # 缩放连续变量 scaler = StandardScaler() cont_vars = ['age', 'income'] df[cont_vars] = scaler.fit_transform(df[cont_vars]) # 选择用于建模的特征 features = ['age', 'income', 'category_var_B', 'category_var_C'] modeling_df = df[['T', 'E'] + features].copy() # 4. 划分训练测试集 train_df, test_df = train_test_split(modeling_df, test_size=0.3, random_state=42) # 5. 拟合Cox模型 cph = CoxPHFitter(penalizer=0.01) # 加入轻微正则化 cph.fit(train_df, duration_col='T', event_col='E') cph.print_summary() # 6. 模型诊断 cph.check_assumptions(train_df, show_plots=True, p_value_threshold=0.05) # 7. 模型评估 # C-index c_index_train = cph.score(train_df) c_index_test = cph.score(test_df) print(f"训练集C-index: {c_index_train:.3f}, 测试集C-index: {c_index_test:.3f}") # 8. 预测与应用 # 预测新客户的180天留存概率 new_customer = pd.DataFrame({ 'age': [0.5], # 标准化后的值 'income': [1.2], 'category_var_B': [1], 'category_var_C': [0] }) survival_prob = cph.predict_survival_function(new_customer, times=[180]) print(f"预测的新客户180天生存概率: {survival_prob.iloc[0, 0]:.3f}")

生存分析是一个强大而深邃的领域,lifelines库让它变得触手可及。从一张简单的Kaplan-Meier曲线开始,到构建一个可以量化风险因素的Cox模型,你已经在用数据解读“时间”的秘密。记住,模型是工具,对业务问题的深刻理解才是核心。多看看你的数据分布,多想想每个系数背后的业务含义,不断用交叉验证和诊断图来审视你的模型,你就能越来越得心应手地运用生存分析来解决实际的预测和归因问题。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询