In survival analysis, independent variables can either be time, categorical, or continuous. Meanwhile, the dependent variable is always binary. Data is non-normal and skewed due to a short or long interval t. In order to analyze data in survival analysis, data must be rearranged in terms of the “time in the study”, rather than the actual year of entry and exit, since dates are not important but the duration is.
There are three main terms in survival analysis, namely: –
- Survival time and event.
- Censoring.
- Survival function and hazard function.
Censoring is permitted in survival analysis. This happens when information about survival time is incomplete. Right censoring occurs when subjects leave the study before an event occurs or the study ends before an event has occurred. Censored data is a type of missing data that helps avoid bias. Other terms are discussed in the next sections.
Import Python libraries
Underneath we imported key Python libraries (a chuck of code that will used in the project) . These libraries helps with creating arrays and profile tables, manipulating dataframes, computing graphs and making estimates.
import numpy as np
import pandas as pd
import warnings
warnings.filterwarnings("ignore")
import seaborn as sns
sns.set("talk","ticks",font_scale=1,font="sans-serif",color_codes=True)
import matplotlib.pyplot as plt
%matplotlib inline
from pylab import rcParams
plt.rcParams["figure.figsize"] = [10,10]
from lifelines.datasets import load_waltons
from lifelines import *
from lifelines.utils import survival_table_from_events
df = load_waltons()
df.head()| T | E | group | |
| 0 | 6.0 | 1 | miR-137 |
| 1 | 13.0 | 1 | miR-137 |
| 2 | 13.0 | 1 | miR-137 |
| 3 | 13.0 | 1 | miR-137 |
| 4 | 16.0 | 1 | miR-137 |
Detect missing values
Underneath we checked whether there were any missing values in the data or not. It is important to note that survival analysis has no strict rules on missing data.
sns.heatmap(df.isnull(),cmap="Blues")
plt.title("Detect Missing Values")
plt.show()
Creating T and E array
Underneath we define object types
T = df["T"]
E = df["E"]T represents time.
E represents events.
Descriptive statistics
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(16,7))
sns.boxplot(df["T"],ax=ax1)
ax1.set_title("Time Box Plot")
sns.distplot(df["E"],ax=ax2)
ax2.set_title("Time Histogram")
ax2.set_ylabel("Related Time Frequency")
plt.show()
The boxplot above shows that there are outliers in the data. Both boxplot and histogram diagrams indicate that the data is skewed to the right. Observations are not spread symmetrically around the true value. The knowledge of mean and variance is not sufficient to summarize data from this distribution.
Survival analysis makes no assumptions about normal distribution.
Life table
In the underneath section, we constructed life tables. To create life tables in Python, we used the survival_from_events() function. Life tables are often used in life insurance to set insurance premiums.
- Removed – the number of subjects who die (or suffer the event of interest) during interval t.
- Observed – the number of subjects who are event-free and well-thought-out as being at risk during interval t.
- Censored – the number of participants who are censored during interval t.
- Entrance – the number of subjects who entered during interval t.
- At-risk – the aggregate number of subjects at risk during interval t.
| removed | observed | censored | entrance | at_risk | |
| event_at | |||||
| 0.0 | 0 | 0 | 0 | 163 | 163 |
| 6.0 | 1 | 1 | 0 | 0 | 163 |
| 7.0 | 2 | 1 | 1 | 0 | 162 |
| 9.0 | 3 | 3 | 0 | 0 | 160 |
| 13.0 | 3 | 3 | 0 | 0 | 157 |
| 15.0 | 2 | 2 | 0 | 0 | 154 |
| 17.0 | 1 | 1 | 0 | 0 | 152 |
| 19.0 | 3 | 3 | 0 | 0 | 151 |
| 22.0 | 4 | 4 | 0 | 0 | 148 |
| 26.0 | 5 | 5 | 0 | 0 | 144 |
| 29.0 | 5 | 5 | 0 | 0 | 139 |
| 32.0 | 1 | 1 | 0 | 0 | 134 |
| 33.0 | 3 | 3 | 0 | 0 | 133 |
| 36.0 | 2 | 2 | 0 | 0 | 130 |
| 38.0 | 2 | 2 | 0 | 0 | 128 |
| 41.0 | 7 | 7 | 0 | 0 | 126 |
| 43.0 | 1 | 1 | 0 | 0 | 119 |
| 45.0 | 10 | 9 | 1 | 0 | 118 |
| 47.0 | 1 | 1 | 0 | 0 | 108 |
| 48.0 | 8 | 8 | 0 | 0 | 107 |
| 51.0 | 3 | 3 | 0 | 0 | 99 |
| 53.0 | 7 | 7 | 0 | 0 | 96 |
| 54.0 | 2 | 2 | 0 | 0 | 89 |
| 56.0 | 18 | 18 | 0 | 0 | 87 |
| 58.0 | 4 | 4 | 0 | 0 | 69 |
| 60.0 | 16 | 15 | 1 | 0 | 65 |
| 61.0 | 11 | 9 | 2 | 0 | 49 |
| 62.0 | 2 | 2 | 0 | 0 | 38 |
| 63.0 | 9 | 9 | 0 | 0 | 36 |
| 66.0 | 3 | 3 | 0 | 0 | 27 |
| 68.0 | 10 | 9 | 1 | 0 | 24 |
| 69.0 | 13 | 12 | 1 | 0 | 14 |
| 75.0 | 1 | 1 | 0 | 0 | 1 |
Kaplan-Meier
The Kaplan-Meier (KM) method is a non-parametric method used to estimate the survival probability from observed survival times. It also referred to as a product limit estimator. This estimator is upright for dealing with ungrouped data. When conducting survival analysis, one should look at the number of patients dying in a specified period, considering that they have already survived.
Underneath , the KaplanMeierFitter().fit() function is used to calculate the Kaplan-Meier survival estimate. The function calculates all distinct times of death occurring in the study.
Main arguments are as follows: –
- a survival object constructed by means of the KaplanMeierFitter().fit() function.
- data containing variables (T and E).
kmf = KaplanMeierFitter().fit(T,E)Cumulative density
Density function refers to the probability of the failure time occurring at exactly time t (out of the whole range of possible t’s.) This function is concerned with finding the
kmf.cumulative_density_| KM_estimate | |
| timeline | |
| 0.0 | 0.000000 |
| 2.0 | 0.000000 |
| 4.0 | 0.000000 |
| 6.0 | 0.006135 |
| 8.0 | 0.012270 |
| 10.0 | 0.030790 |
| 12.0 | 0.030790 |
| 14.0 | 0.049310 |
| 16.0 | 0.061656 |
| 18.0 | 0.067830 |
| 20.0 | 0.086350 |
| 22.0 | 0.111043 |
| 24.0 | 0.111043 |
| 26.0 | 0.141910 |
| 28.0 | 0.141910 |
| 30.0 | 0.172776 |
| 32.0 | 0.178949 |
| 34.0 | 0.197469 |
| 36.0 | 0.209816 |
| 38.0 | 0.222163 |
| 40.0 | 0.222163 |
| 42.0 | 0.265376 |
| 44.0 | 0.271549 |
| 46.0 | 0.327109 |
| 48.0 | 0.383183 |
| 50.0 | 0.383183 |
| 52.0 | 0.401875 |
| 54.0 | 0.457949 |
| 56.0 | 0.570097 |
| 58.0 | 0.595019 |
| 60.0 | 0.688476 |
| 62.0 | 0.759079 |
| 64.0 | 0.819310 |
| 66.0 | 0.839386 |
| 68.0 | 0.899616 |
| 70.0 | 0.985659 |
| 72.0 | 0.985659 |
| 74.0 | 0.985659 |
| 76.0 | 1.000000 |
| 78.0 | 1.000000 |
| 80.0 | 1.000000 |
| 82.0 | 1.000000 |
| 84.0 | 1.000000 |
| 86.0 | 1.000000 |
| 88.0 | 1.000000 |
| 90.0 | 1.000000 |
| 92.0 | 1.000000 |
| 94.0 | 1.000000 |
| 96.0 | 1.000000 |
| 98.0 | 1.000000 |
kmf.survival_function_In survival analysis, event outcomes and follow up times are used to estimate the survival function.
| KM_estimate | |
| timeline | |
| 0.0 | 1.000000 |
| 2.0 | 1.000000 |
| 4.0 | 1.000000 |
| 6.0 | 0.993865 |
| 8.0 | 0.987730 |
| 10.0 | 0.969210 |
| 12.0 | 0.969210 |
| 14.0 | 0.950690 |
| 16.0 | 0.938344 |
| 18.0 | 0.932170 |
| 20.0 | 0.913650 |
| 22.0 | 0.888957 |
| 24.0 | 0.888957 |
| 26.0 | 0.858090 |
| 28.0 | 0.858090 |
| 30.0 | 0.827224 |
| 32.0 | 0.821051 |
| 34.0 | 0.802531 |
| 36.0 | 0.790184 |
| 38.0 | 0.777837 |
| 40.0 | 0.777837 |
| 42.0 | 0.734624 |
| 44.0 | 0.728451 |
| 46.0 | 0.672891 |
| 48.0 | 0.616817 |
| 50.0 | 0.616817 |
| 52.0 | 0.598125 |
| 54.0 | 0.542051 |
| 56.0 | 0.429903 |
| 58.0 | 0.404981 |
| 60.0 | 0.311524 |
| 62.0 | 0.240921 |
| 64.0 | 0.180690 |
| 66.0 | 0.160614 |
| 68.0 | 0.100384 |
| 70.0 | 0.014341 |
| 72.0 | 0.014341 |
| 74.0 | 0.014341 |
| 76.0 | 0.000000 |
| 78.0 | 0.000000 |
| 80.0 | 0.000000 |
| 82.0 | 0.000000 |
| 84.0 | 0.000000 |
| 86.0 | 0.000000 |
| 88.0 | 0.000000 |
| 90.0 | 0.000000 |
| 92.0 | 0.000000 |
| 94.0 | 0.000000 |
| 96.0 | 0.000000 |
| 98.0 | 0.000000 |
Survival curve and cumulative density curve
In a survival curve, the x-axis represents time in years, and the y-axis represents the probability of surviving or the proportion of people surviving. As we can see, the survival curve shows a step function rather than a smooth curve – it does not drop sharply towards zero. We can state that we have decent survival. We plotted the curve using the plot_survival_function() function. Furthermore, we plotted the Kaplan-Meier curve for the cumulative risk.
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(16,7))
kmf.plot_survival_function(ax=ax1)
kmf.plot_cumulative_density(ax=ax2)
plt.show()
Notice that the survival probability is 100% for 5 years and then drops to 98%. The median survival is 45 years (i.e., 50% of the population survive 45 years).