Python による Kaplan-Meier 法、ログランク検定、Cox 比例ハザードモデル

lifelines を使い、ログランク検定および Cox 比例ハザードモデル分析を行う。

必要なライブラリを import する。

!pip install lifelines #必要な場合のみ
from lifelines import KaplanMeierFitter, CoxPHFitter
from lifelines.statistics import logrank_test
import pandas as pd

通常は Excel 等で作成した CSV ファイルを R に読み込んで分析する。しかし、ここでは便宜的に以下の簡単なデータフレームを作成し、data という変数に格納しておく。打ち切り = 0 なら打ち切りありとしている。

data = pd.DataFrame({"患者番号":[1,2,3,4,5,6,7,8,9,10], "抗がん剤":["A","A","A","A","A","B","B","B","B","B"], "生存期間":[20,30,40,50,60,55,75,105,115,120], "打ち切り":[1,1,1,1,0,1,1,1,1,0], "全身状態":[2,3,2,3,2,0,1,0,0,1]})

data の内容を確認する。

data

data から A 群、B 群をそれぞれ抽出し、Kaplan-Meier モデルに適合させる。

#A 群
A = data.query('抗がん剤 == "A"')
kmf_A = KaplanMeierFitter()
kmf_A.fit(A.生存期間, A.打ち切り)

#B 群
B = data.query('抗がん剤 == "B"')
kmf_B = KaplanMeierFitter()
kmf_B.fit(B.生存期間, B.打ち切り)

Kaplan-Meier 法により、抗がん剤 A 群と抗がん剤 B 群を比較する。

#生存期間中央値
kmf_A.median_survival_time_
kmf_B.median_survival_time_

#生存確率、信頼区間
kmf_A.confidence_interval_survival_function_
kmf_B.confidence_interval_survival_function_

#5年 (60カ月) 生存率
kmf_A.predict(60)
kmf_B.predict(60)

#Kaplan-Meier 曲線の図示
kmf_A.plot() and kmf_B.plot()

ログランク検定を行う。

result = logrank_test(A.生存期間, B.生存期間, A.打ち切り, B.打ち切り)
result.print_summary()

p 値をもっと正確に知りたい場合は以下のようにする。

result.p_value

次に、Cox 比例ハザードモデル分析を行う。lifelines で Cox 比例ハザードモデル分析を行う場合、data に “A”、”B” などの文字が含まれているとうまく解析できない。また、変量名に “生存期間” などの日本語(2バイト文字)を使うとグラフのラベルに表示されない。

以上を踏まえ、あらためて以下の簡単なデータフレームを使う。

data_cox = pd.DataFrame({"Group":[0,0,0,0,0,1,1,1,1,1], "Duration":[20,30,40,50,60,55,75,105,115,120], "Censor":[1,1,1,1,0,1,1,1,1,0], "PS":[2,3,2,3,2,0,1,0,0,1]})

data_cox の内容を確認する。

data_cox

Cox 比例ハザードモデル分析を行う。

cph = CoxPHFitter()
cph.fit(data_cox, "Duration", "Censor")
cph.print_summary()

各変量の信頼区間をグラフに描出することもできる。

cph.plot()