최근에 python으로도 생존분석이 가능하다는 사실을 알게 되어 여기에 간단히 적어봅니다. 사실은 lifelines 라는 패키지의 설명서(https://lifelines.readthedocs.io/en/latest/index.html# )를 그냥 번역한 수준인데 나중에 공부하고 나서 내용을 수정하거나 추가할 수 있습니다.
생존분석 패키지가 몇 가지 있는데 여기서는 lifelines를 가지고 하겠습니다. 그 외에 scikit-survival 도 참고하시기 바랍니다.(https://scikit-survival.readthedocs.io/en/stable/user_guide/00-introduction.html )
Google colab에서 lifeline을 설치하고 import로 가져옵니다.
lifelines.datasets에서 load_waltons라는 데이터를 가져옵니다. T는 survival time, E는 event (1=event, 0=censoring)를 의미합니다. group의 값은 miR-137과 control 두 가지입니다.
Kaplan-Meier (아래 Meyer는 오타입니다) 그림은 lifelines에 있는 KaplanMeierFitter의 KaplanMeierFitter() 함수로 추정합니다. kmf.fit(T, E)로 지정하고 나서 kmf.plot_survival_function()으로 그림을 얻을 수 있습니다.
누적 밀도함수는 plot_cumulative_density()를 이용합니다.
kmf.fit(T, E, timeline=range(0, 100, 2)) 로 시간 구간을 바꿔 보았습니다.
중앙 생존 시간( 개체의 50%에서 event가 관측된 시점)은 media_survival_times에 있는 media_survival_time_ 으로 구할 수 있습니다. 여기서는 56이 중앙 생존 시간이고 신뢰구간은 (54, 58)입니다.
비모수적 추정법인 Kaplan-Meier와 모수적 생존시간인 Weibull, Exponential Lognormal 등으로 적합시간 생존곡선을 구해서 비교합니다.
Control과 miR-137 집단 간의 Kaplan-Meier 곡선을 각각 그려서 비교했습니다. miR-137 집단의 위험도가 더 높은 듯 합니다.
survival_table_from_event를 불러오면 생존구간 표를 얻을 수 있습니다.
lifelines에서 Cox proportional 회귀분석은 CoxPHFitter로 추정합니다. 여기서 cph로 이름을 두고 ,cph.fit( 데이터, '생존시간', event_col='이벤트')로 구한 후, cph.print_summary()로 추정값을 확인합니다.
추정된 계수에 대한 plot은 다음과 같습니다.
모수적 생존회귀분석 기법인 Weibull 회귀모형은 WeibullAFTFitter()를 이용합니다.
Weibull 모수적 생존회귀모형에서 추정된 계수와 신뢰구간 그림입니다.
또다는 비모수적 추정법인 Aalen의 additive model을 이용하여 추정한 후 Cox, Weibull, Aalen 세 가지 추정 생존곡선을 그려 비교합니다.
아래쪽 내용은 (https://medium.com/the-researchers-guide/survival-analysis-in-python-km-estimate-cox-ph-and-aft-model-5533843c5d5d )의 내용을 구현하려다 일단 데이터만 열어본 겁니다.
colab에서 외부 데이터를 여는 방법은 Google의 gdrive에 데이터를 올리고 불러오는 방법이 있고, 아래에는 로컬 컴퓨터 상에 있는 데이터를 불러온 후 io.BytesIO() 를 이용하여 불러왔습니다.
일단 소개삼아 대강의 내용을 적어보았고, 공부하는 대로 업데이트할 예정입니다. 코드는 Surv_lifelines.ipynb 라는 이름으로 첨부했습니다. Jupiter notebook이나 Google colab에서 열 수 있습니다.
첫댓글 좋은 튜토리얼 감사합니다. 저는 요즘 scikit survival을 사용하고 있어요.
봐주셔서 감사합니다.