Survival analysis analyzes data where the outcome variable is the time until the occurrence of an event of interest. The distribution of the event times is typically described by a survival function.
Non-parametric Distribution Fitting (NPDF) provides an estimate of the survival function without making any assumptions concerning the distribution of the data. NPDF includes Kaplan-Meier estimation, life tables, and specialized extension algorithms to support left censored, interval censored, and recurrent event data.
Survival analysis analyzes data where the outcome variable is the time until the occurrence of an event of interest.The distribution of the event times is typically described by a survival function.
Parametric Distribution Fitting (PDF) provides an estimate of the survival function by comparing the functions for several known distributions (exponential, Weibull, log-normal, and log-logistic) to determine which, if any, describes the data best. In addition, the distributions for two or more groups of cases can be compared.
Parametric regression modeling (PRM) is a survival analysis technique that incorporates the effects of covariates on the survival times. PRM includes two model types: accelerated failure time and frailty. Accelerated failure time models assume that the relationship of the logarithm of survival time and the covariates is linear. Frailty, or random effects, models are useful for analyzing recurrent events, correlated survival data, or when observations are clustered into groups.
PRM automatically selects the survival time distribution (exponential, Weibull, log-normal, or log-logistic) that best describes the survival times.
Example code:
import com.ibm.spss.ml.survivalanalysis.ParametricRegression
var prm = ParametricRegression().
setBeginField("startTime").
setEndField("endTime").
setStatusField("status").
setPredictorFields(Array("age", "surgery", "transplant")).
setDefinedStatus(
DefinedStatus(
failure = Some(StatusItem(points = Some(Points("0.0")))),
intervalCensored = Some(StatusItem(points = Some(Points("1.0")))))).
val prmModel = prm.fit(data)
val PMML = prmModel.toPMML()
val statXML = prmModel.statXML()
val predictions = prmModel.transform(data)
predictions.show()
I was a post-graduate in Mechanical Engineering when I joined the analytics industry as a fresher.
The only background I had in analytics industry was based on a few courses in Operations Research. I was scared to join the industry, where I knew lesser than a student of statistics in 12th standard. In less than a month, I realized that Analytics industry does not require you to have a masters in statistics or economics, but requires structured thinking with sharp mathematical reflexes. It took me no more than 3 months to build and implement my first logistic regression model. Most of the people in India still use basic analytics tools such as CART, regressions and time series. We still are scared of using complex statistical techniques, such as neural network and survival analysis. Last year, I used survival analysis in one of the analytics projects and realized the power of the tool without getting into the Limbo of statistics behind the tool. This article will help you find if the survival analysis is the right tool for your next project. The article will end with a case study, which we will solve using survival analysis in the next article.
How is survival analysis different from regression models?
Regression models have a single output function. In case of a logistic regression, the output is the response function which can only take two values. In any other model, we can define the output function in a single objective function. For instance, if we are building a customer attrition model, which predicts whether a customer will attrite in next 3 months, following is the objective of the logistic model :
f(x) = 0 if no attrition in next 3 months
= 1 if customer attrited in next 3 months
Say, we want to profile the customers, who are likely to attrite early and restrict the acquisition of such profile customers. Let’s assume for simplicity that there are only 2 variables : Gender and Tenure. We build a logistic model in Jan ’13 and found that out of 100 Males 30 attrite till Jan’13, whereas out of 100 females only 10 attrite till Jan’13. My model profiles females as better profile, but because of some reason we did not implement this factor into our acquisition strategy. Now, we stand in Jul’13 and if we look at the same population as considered in the logistic model, out of 100 Males 35 attrite and out of 100 Females 55 attrite.
The results seem to have swapped. In last 6 months, females saw high attrition, whereas male population seems to be very stable over this period. Now you observe that females had a lower tenure in Jan’13 as compared to the male population. The possible solution in this case is to take same month tranches/acquisitions to build a model. You now take only customer who was acquired before Jan’12. You get a population of 50 Males and 10 Females. However, you have reduced the noise coming from new tranches, you have also reduced the population you are building a model on. How do you address this issue?
Such data for whom the results (attrition in this case) are unknown, are called censored data. We can include this data without compromising on the model accuracy in a survival analysis model. This is because the output or target variable of a survival analysis is a combination of death (attrition in this case) and time on books (tenure of customer).
Applicationsof survival analysis
There are four major applications of survival analysis into analytics:
1. Business Planning : Profiling customers who has a higher survival rate and make strategy accordingly.
2. Lifetime Value Prediction : Engage with customers according to their lifetime value
3. Active customers : Predict when the customer will be active for the next time and take interventions accordingly.
4. Campaign evaluation : Monitor effect of campaign on the survival rate of customers.
Following are some industrial specific applications of survival analysis :
• Banking – customer lifetime and LTV
• Insurance – time to lapsing on policy
• Mortgages – time to mortgage redemption
• Mail Order Catalogue – time to next purchase
• Retail – time till food customer starts purchasing non-food
• Manufacturing – lifetime of a machine component
• Public Sector – time intervals to critical events
Case study
You are the head of the analytics team with an online Retail chain Mazon. You have received a limited number of offers which costs you $200/customer targeted . The offer breaks even if a customer makes a purchase of minimum $20,000 in his entire lifetime. You want to target customers who are likely to make a purchase of $20,000 as early as possible. You have their relationship card segment (Platinum cards are expensive and Gold card is cheaper) and their response rate on last campaign they were targeted with. You have been targeting customers with this offer for past 1 year. Now you want to learn from the past response data and target accordingly. You need to find which customer base should you target for this offer. You have data for a similar test campaign in the past, based on which you will have to build the strategy.
You will have to use survival analysis in this case because the dependent variable is the time to respond the campaign. This again contains censored data which are people who did not respond till date. We will solve this case study in the next article where we will lay out a step by step process doing a survival analysis to find profile of customers who respond early.
End notes
In the last three years, I have realized that analytics project are not always created by the business. The business itself is constrained by what value they have seen in the past coming from the analytics team. Learning new modelling techniques and learning from nextgen analytics project make us capable of seeing the value beyond what business can see. We will like to know the new techniques you have learnt in recent past and their applications. The objective to start a discussion on survival analysis here is not restricted to only this technique.
Did you find any opportunity in your line of business where you can implement survival analysis? Did you find this article helpful? Have you worked on other cutting edge modelling techniques in the recent past? Were the results encouraging? Do let us know your thoughts in the comments below.
Survival analysis techniques are useful for studying the longevity of groups of elements or individuals, taking into account time considerations and right censorship. This episode explores how survival analysis can describe marriages, in particular, using the non-parametric Cox proportional hazard model.
This episode discusses some good summaries of survey data on marriage and divorce which can be found here.
The python lifelines library is a good place to get started for people that want to do some hands on work.
Bayesian survival analysis in A Song of Ice and Fire
This post originally appeared on Allen Downey's personal blog.
Last fall I taught an introduction to Bayesian statistics at Olin College. My students worked on some excellent projects, and I invited them to write up their results as guest articles for my blog.
One of the teams applied Bayesian survival analysis to the characters in A Song of Ice and Fire, the book series by George R. R. Martin. Using data from the first 5 books, they generate predictions for which characters are likely to survive and which might die in the forthcoming books.
With Season 5 of the Game of Thrones television series starting on April 12, we thought this would be a good time to publish their report.
Bayesian survival analysis in A Song of Ice and Fire
By Erin Pierce and Ben Kahle
The Song of Ice and Fire series has a reputation for being quite deadly. No character, good or bad, major or minor is safe from Martin's pen. The reputation is not unwarranted; of the 916 named characters that populate Martin's world, a third have died, alongside uncounted nameless ones.
In this report, we take a closer look at the patterns of death in the novels and create a Bayesian model that predicts the probability that characters will survive the next two books.
Using data from A Wiki of Ice and Fire, we created a dataset of all 916 characters that appeared in the books so far. For every character, we know what chapter and book they first appeared, if they are male or female, if they are part of the nobility or not, what major house they are loyal to, and, if applicable, the chapter and book of their death. We used this data to predict which characters will survive the next couple books.
Methodology
We extrapolated the survival probabilities of the characters through the seventh book using Weibull distributions. A Weibull distribution provides a way to model the hazard function, which measures the probability of death at a specific age. The Weibull distribution depends on two parameters, k and lambda, which control its shape.
To estimate these parameters, we start with a uniform prior. For each alive character, we check how well that value of k or lambda predicted the fact that the character was still alive by comparing the calculated Weibull distribution with the character's hazard function. For each dead character, we check how well the parameters predicted the time of their death by comparing the Weibull distribution with the character's survival function.
The main code used to update these distributions is:
class GOT(thinkbayes2.Suite, thinkbayes2.Joint):
def Likelihood(self, data, hypo):
"""Determines how well a given k and lam predict the life/death of a character """
age, alive = data k, lam = hypo
if alive:
prob = 1-exponweib.cdf(age, k, lam)
else:
prob = exponweib.pdf(age, k, lam)
return prob
def Update(k, lam, age, alive):
"""Performs the Baysian Update and returns the PMFs of k and lam"""
joint = thinkbayes2.MakeJoint(k, lam)
suite = GOT(joint)
suite.Update((age, alive))
k = suite.Marginal(0, label=k.label), lam = suite.Marginal(1, label=lam.label)
return k, lam
def MakeDistr(introductions, lifetimes,k,lam):
"""Iterates through all the characters for a given k and lambda. It then updates the k and lambda distributions"""
k.label = 'K'
lam.label = 'Lam'
print("Updating deaths")
for age in lifetimes:
k, lam = Update(k, lam, age, False)
print('Updating alives')
for age in introductions:
k, lam = Update(k, lam, age, True)
return k,lam
For the Night’s Watch, this lead to the posterior distribution in Figure 3:
The distribution for lambda is quite tight, around 0.27, but the distribution for k is broader.
To translate this back to a survival curve, we took the mean of k and lambda, as well as the 90 percent credible interval for each parameter. We then plot the original data, the credible interval, and the survival curve based on the posterior means.
Jon Snow
Using this analysis, we can can begin to make a prediction for an individual character like Jon Snow. At the end of A Dance with Dragons, the credible interval for the Night's Watch survival (Figure 4) stretches from 36 percent to 56 percent. The odds are not exactly rosy that Jon snow is still alive. Even if Jon is still alive at the end of book 5, the odds that he will survive the next two books drop to between 30 percent and 51 percent.
The credible interval closely encases the data, and the mean-value curve appears to be a reasonable approximation.
However, it is worth considering that Jon is not an average member of the Night's Watch. He had a noble upbringing and is well trained at arms. We repeated the same analysis with only members of the Night's Watch considered noble due to their family, rank, or upbringing. There have only been 11 nobles in the Night's Watch, so the credible interval as seen in Figure 5 is understandably much wider, however, the best approximation of the survival curve suggests that a noble background does not increase the survival rate for brothers of the Night's Watch.
When only noble members of the Night’s Watch are included, the credible interval widens significantly and the lower bound gets quite close to zero.
The houses of ASOIF
The 90 percent credible intervals for all of the major houses. This includes the 9 major houses, the Night’s Watch, the Wildlings, and a "None" category which includes non-allied characters.
90 percent credible interval for Arryn (Blue), Lannister (Gold), None (Green), and Stark (Grey)
90 percent credible interval for Tyrell(Green), Tully(Blue), Baratheon(Orange), and Night’s Watch (Grey)
90 percent credible interval for Martell(Orange), Targaryen (Maroon), Greyjoy (Yellow), and Wildling (Purple)
These intervals, shown in Figures 6, 7, and 8, demonstrate a much higher survival probability for the houses Arryn, Tyrell, and Martell. Supporting these results, these houses have stayed out of most of the major conflicts in the books, however this also means there is less information on them. We have 5 or fewer examples of dead members for those houses, so the survival curves don’t have very many points. This uncertainty is reflected in the wide credible intervals.
In contrast, our friends in the north, the Starks, Night’s Watch, and Wildlings have the lowest projected survival rates and smaller credible intervals given their warring positions in the story and the many important characters included amongst their ranks. This analysis considers entire houses, but there are also additional ways to sort the characters.
Men and women
While A Song of Ice and Fire has been lauded for portraying women as complex characters who take an a variety of roles, there are still many more male characters (769) than female ones (157). Despite a wider credible interval, the women tend to fare better than their male counterparts, out-surviving them by a wide margin as seen in Figure 9.
The women of Westeros appear to have a better chance of surviving then the men.
Class
The ratio between noble characters(429) and smallfolk characters (487) is much more even than gender and provides an interesting comparison for analysis. Figure 10 suggests that while more smallfolk tend to die quickly after being introduced, those that survive their introductions tend to live for a longer period of time and may in fact outpace the nobles.
The nobility might have a slight advantage when introduced, but their survival probability continues to fall while the smallfolk’s levels much more quickly.
Selected Characters
The same analysis can be extended to combine traits, sorting by gender, house, and class to provide a rough model for individual characters. One of the most popular characters in the books is Arya and many readers are curious about her fate in the books to come. The category of noblewomen loyal to the Starks also includes other noteworthy characters like Sansa and Brienne of Tarth (though she was introduced later). Other intriguing characters to investigate are the Lannister noblewomen Cersei and poor Myrcella. As it turns out, not a lot noble women die. In order to get more precise credible intervals for the specific female characters we included the data of both noble and smallfolk women.
While both groups have very wide ranges of survival probabilities, the Lannister noblewomen may be a bit more likely to die than the Starks.
The data presented in Figure 11 is inconclusive, but it looks like Arya has a slightly better chance of survival than Cersei.
Two minor characters we are curious about are Val, the wildling princess, and the mysterious Quaithe.
Representing the survival curves of more minor characters, Quaithe and Val have dramatically different odds of surviving the series.
They both had more data than the Starks and Lannisters, but they have the complication that they were not introduced at the beginning of the series. Val is introduced at 2.1 books, and so her chances of surviving the whole series are between 10 percent and 53 percent, which are not the most inspiring of chances.
Quaithe is introduced at 1.2 books, and her chances are between 58 percent and 85 percent, which are significantly better than Val’s. These curves are shown in Figure 12.
For most of the male characters (with the exception of Mance), there was enough data to narrow to house, gender and class.
The survival curves of different classes and alliances of men shown through various characters.
Figure 13 shows the Lannister brothers with middling survival chances ranging from 35 percent to 79 percent. The data for Daario is less conclusive, but seems hopeful, especially considering he was introduced at 2.5 books. Mance seems to have to worst chance of surviving until the end. He was introduced at 2.2 books, giving him a chance of survival between 19 percent and 56 percent.
The survival curves of different classes and alliances of men shown through various characters.
Some characters who many wouldn’t mind seeing kick the bucket include Lord Walder Frey and Theon Greyjoy. However, Figure 14 suggests that neither are likely meet untimely (or in Walder Frey’s case, very timely) deaths. Theon seems likely to survive to the bitter end. Walder Frey was introduced at 0.4 books, putting his chances at 44 percent to 72 percent. As it is now, Hoster Tully may be the only character to die of old age, so perhaps Frey will hold out until the end.
Conclusion
Of course who lives and who dies in the next two books has more to do with plot and storyline than with statistics. Nonetheless, using our data we were able we were able to see patterns of life and death among groups of characters. For some characters, especially males, we are able to make specific predictions of how they will fare in the next novels. For females and characters from the less central houses, the verdict is still out.
Most characters were fairly easy to classify, but there are always edge cases.
Gender - This was the most straight forward. There are not really any gender-ambigous characters.
Nobility - Members of major and minor Westeros houses were counted as noble, but hedge knights were not. For characters from Essos, I used by best judgement based on money and power, and it was usually an easy call. For the wildlings, I named military leaders as noble, though that was often a blurry line. For members of the Night’s Watch, I looked at their status before joining in the same way I looked at other Westeros characters. For bastards, we decided on a case by case basis. Bastards who were raised in a noble family and who received the education and training of nobles were counted as noble. Thus Jon Snow was counted as noble, but someone like Gendry was not.
Death - Characters that have come back alive-ish (like Beric Dondarrion) were judged dead at the time of their first death. Wights are not considered alive, but others are. For major characters whose deaths are uncertain, we argued and made a case by case decision.
Houses - This was the trickiest one because some people have allegiances to multiple houses or have switched loyalties. We decided on a case by case basis. The people with no allegiance were of three main groups:
— People in Essos who are not loyal to the Targaryens.
— People in the Riverlands, either smallfolk who’s loyalty is not known, or groups like the Brotherhood Without Banners or the Brave Companions with ambiguous loyalty.
— Nobility that are mostly looking out for their own interests, like the Freys, Ramsay Bolton, or Petyr Baelish.
Time-to-event outcomes are common in medical research as they offer more information than simply whether or not an event occurred. To handle these outcomes, as well as censored observations where the event was not observed during follow-up, survival analysis methods should be used. Kaplan-Meier estimation can be used to create graphs of the observed survival curves, while the log-rank test can be used to compare curves from different groups. If it is desired to test continuous predictors or to test multiple covariates at once, survival regression models such as the Cox model or the accelerated failure time model (AFT) should be used. The choice of model should depend on whether or not the assumption of the model (proportional hazards for the Cox model, a parametric distribution of the event times for the AFT model) is met. The goal of this paper is to review basic concepts of survival analysis. Discussions relating the Cox model and the AFT model will be provided. The use and interpretation of the survival methods model are illustrated using an artificially simulated dataset.
Keywords: Survival analysis, Cox proportional hazards model, accelerated failure time model
Clinical outcomes come in a variety of statistical forms. Some are continuous, such as systolic blood pressure, and can be easily analyzed with linear regression. Others, such as mortality or myocardial infarction (MI), are distinct events and have forms that are slightly more complex to analyze statistically. If the clinical outcome observed is “either-or,” such as if a patient suffers an MI or not, logistic regression can be used. However, if the information on the time to MI is the observed outcome, data are analyzed using statistical methods for survival analysis. It should be noted that despite the name “survival analysis,” methods can be used in any time-to-event outcome, such as the time until a patient experiences an MI or the time to hospitalization. Such studies will discuss event-free survival which is the proportion of subjects who have not yet experienced an event.
There are important clinical and statistical reasons for investigating a time-to-event outcome using survival analysis. For example, consider a study that found that the final observed proportion of events between two treatment groups is identical. However, if one group had all events occur shortly after randomization, while the other had no events until just before the end of follow-up then the two treatments would logically be considered to have different clinical effects despite the identical proportions at the end of follow-up. Similarly, if all-cause mortality is the outcome, then a sufficiently long follow-up would reveal equal survival proportions of 0% between any groups. In such cases, the time to an event contains much more clinical information than whether or not the event occurred. There is much more statistical information as well, as survival analyses tend to have greater statistical power to detect a significant treatment or exposure effect than methods for binary outcomes such as logistic regression.
It is typical in these types of studies to have subjects who did not experience the event before the end of a study or dropped out before the event of interest occurs. These subjects are said to be right-censored. Although these may seem to be cases of missing data as the time-to-event is not actually observed, these subjects are highly valuable as the observation that they went a certain amount of time without experiencing an event is itself informative. One of the most important properties of survival methods is their ability to handle such censored observations which are ignored by methods such as a t-test (or analysis of variance) for comparing survival times of two (or more) groups and linear regression. It is worth noting that survival methods of analyses can handle other types of censoring such as left-censoring, where a subject had already experienced the event at the time they enrolled in the study, and interval- censoring, where it is known only that the event happened in a particular interval, say between the last and current visit. This article will only discuss right-censored data.
OVERVIEW OF SURVIVAL METHODS AND THEIR USE IN NUCLEAR CARDIOLOGY
In comparing the survival distributions of two or more groups (for example, new therapy vs standard of care), Kaplan-Meier estimation1 and the log-rank test2 are the basic statistical methods of analyses. These are non-parametric methods in that no mathematical form of the survival distributions is assumed. If an investigator is interested in quantifying or investigating the effects of known covariates (e.g., age or race) or predictor variables (e.g., blood pressure), regression models are utilized. As in the conventional linear regression models, survival regression models allow for the quantification of the effect on survival of a set of predictors, the interaction of two predictors, or the effect of a new predictor above and beyond other covariates.
Among the available survival regression models, the Cox proportional hazards model developed by Sir David Cox3 has seen great use in epidemiological and medical studies, and the field of nuclear cardiology is no exception. What follows are some examples of Cox models being used in nuclear cardiology. Xu et al4looked at how myocardial scarring (assessed with positron emission tomography [PET] or single photon emission computed tomography [SPECT]) and other demographic and medical history factors predicted mortality in patients with advanced heart failure who received cardiac resynchronization therapy. Bourque et al5 looked at how left ventricular ejection fraction (LVEF, assessed with angiography) and nuclear summed rest score (SRS, assessed with SPECT) interacted to change the risk of mortality. Hachamovitch and Berman6 looked at the incremental prognostic value of myocardial perfusion SPECT (MPS) parameters in the prediction of sudden cardiac death. Nakata et al7 looked at how the heart-to-mediastinum ratio (assessed with metaiodobenzylguanidine [MIBG] imaging) predicted cardiac death.
Survival models other than the Cox model have been used in nuclear cardiology as well. For example, in a study of diagnosis strategies for quantifying myocardial perfusion with SPECT, Duvall et al8 utilized a log-normal survival model, a member of the parametric family of regression survival models, since initial data exploration revealed that the proportional hazards assumption of the Cox model was invalid. While this is an excellent example of when to utilize other survival models, it has been more common to see such data presented in conjunction with a Cox model analysis. In earlier studies of MPS-derived predictors of cardiac events, Hachamovitch et al9 used Cox models to identify significant predictors and parametric models, specifically the accelerated failure time (AFT) model, to make estimates of the time to certain percentiles of survival. An identical analysis strategy was used by the research group comprised of Cuocolo, Acampa, Petretta, Daniele et al10–13 in their research of the impact of various SPECT-derived predictors on the occurrence of cardiac events.
In analyzing survival or time-to-event data, there are several important quantities of interest to define. One of the most important quantities is the survival function, denoted by S(t), which provides the probability of survival at a given time. To illustrate, suppose that death is the event of interest, and time is measured in years from study enrollment. Examples of survival functions for two groups are displayed in Figure 1A. Group 1 has a higher risk of experiencing death than Group 2, because its survival curve decreases faster than the curve for Group 2. It is expected that about 61% of Group 1 and about 76% in Group 2 will survive past 5 years of study enrollment; while about 25% in Group 1 and 47% in Group 2 will survive past 10 years. The median time for each group is found by looking at the intersection between the horizontal line associated with the probability of survival equal to 0.5 and the particular survival curve. The median time is about 6.2 years for Group 1 and about 9.4 years for Group 2 using Figure 1A.
Plot of sample Weibull survival functions (A) and the corresponding hazard functions (B) The solid red curve represents the hazard function of Group 1, and the blue dashed curve represents the hazard function of Group 2.
Another important quantity in the analysis of survival data is the rate at which a person who is event-free at a given point in time will instantaneously experience the event. This rate is quantified by the hazard function, denoted by h(t). The value of the hazard function is not a probability, but it is an indicator of the risk of experiencing the event. The higher the value of the hazard function, the higher the risk of event. Mathematically, the hazard function is related to how fast the survival function decreases over time. Therefore, the faster the survival function decreases, the higher the hazard. Returning to Figure 1A, since the survival curve for Group 1 decreases faster than Group 2, it is expected that the hazard function for Group 1 is higher than that of Group 2 as seen in Figure 1B which displays the hazard function corresponding to the survival curves in Figure 1A. It is also worth noting that both hazard functions increase over time.
The hazard function plays a very important role in survival analysis. The Cox proportional hazards model, the most popularly used survival regression model, investigates the relationship of predictors and the time-to-event through the hazard function. It assumes that the predictors have a multiplicative effect on the hazard and that this effect is constant over time, i.e.,
[Math Processing Error]
where h(t|x) is the hazard at time t for a subject with a set of predictors x1,…,xp, h0(t) is the baseline hazard function, and β1,…,βp are the model parameters describing the effect of the predictors on the overall hazard. Therefore, the interpretation of the Cox model is done using hazard ratios (HR), defined as the ratio of the predicted hazard function under two different values of a predictor variable. A hazard ratio greater than 1 means the event is more likely to occur, and a ratio less than one means an event is less likely to occur. A hazard ratio of 1 means the predictor has no effect on the hazard of the event. Also, due to the regression framework of the model, one can get hazard ratio estimates that are controlled for other covariates in the model such as age, sex, and race.
Cox models have achieved great popularity, because they do not require the investigator to assume a particular survival distribution for the data. Instead, these models use a hazard function. In estimating the baseline hazard function, a Cox model uses the so-called Aalen-Breslow estimator, which is a generalization of the non-parametric Nelson-Aalen estimator of the cumulative hazard function.14 The lack of a parametric form of the survival distribution gives the Cox model its other name, the semiparametric model, since the only parameters to estimate in the model are those describing how the predictors affect the hazard.
Despite the seeming ease at which the Cox model can be implemented, there is one significant assumption which must be checked. The Cox model assumes proportional hazards between the values of the predictors regardless of how the underlying hazard may change over time (admittedly, the addition of time-varying covariates affects this assumption). Since the proportional hazards model is built entirely around this assumption, if it happens to be invalid for a set of predictors in a given dataset, then the Cox model should not be used on that dataset, and any results would be questionable. A way around this issue is to fit a stratified Cox model for which the baseline hazard can be different from stratum to stratum or fit a model that includes time-varying covariates.14 The latter case can address broader cases but the interpretation of the resulting model will not be straight forward.
It is important to examine the proportional hazards assumption for all predictors tested in a Cox model. A way to assess proportional hazards for a continuous predictor is to plot the Schoenfeld residuals vs time. If the proportional hazards assumption is valid, then the Schoenfeld residuals should look like a random scatter around zero.15 When examining a categorical predictor like a medical treatment or disease status, it is easiest to compare a log-log transformation of the Kaplan-Meier survival curves for the different categories. Under proportional hazards, the curves should be approximately parallel and should not intersect after time apart. Note that a bit of crossing at early time points may be a product of noise in the survival estimates and may not constitute a violation of the proportional hazards assumption. Both Schoenfeld residuals and Kaplan-Meier estimates can be easily obtained from statistical software such as SAS.
An alternative to the Cox model is a parametric survival model wherein a particular form of the survival distribution is assumed. There are several classes of parametric models: (1) parametric proportional hazards model which takes the form of the Cox model but assumes a parametric form on the baseline hazard; (2) the additive hazards model where the predictors affect the hazard function in an additive manner instead of multiplicative; and (3) the AFT model which is most similar to conventional linear regression.
Unlike the proportional hazards model and the additive hazards model, which examine how predictors affect the hazard function, the AFT model postulates a direct relationship between the predictors and the survival time, making its interpretation easier. Suppose there is only one predictor denoting the presence or absence a disease. The AFT model assumes that the disease either accelerates or decelerates the rate of decrease of the survival function. In other words, if S1(t) and S2(t) denote the survival functions of the presence and absence of the disease, respectively, then the AFT model assumes the relationship S1(t) = S2(ηt), where η is the acceleration factor. If median time-to-event is of interest, then the AFT model implies that the median time for those with the disease is η times the median time of those without the disease. If η > 1 (or η < 1), then the group with the disease has a longer (or shorter) median time-to-event relative to those without the disease. This interpretation is true for other percentiles of the survival distribution.
In general, the AFT model can be expressed two ways:
[Math Processing Error]
or
[Math Processing Error]
where T is the time-to-event (the failure time); x1,…,xp, and β1,…,βp are predictor variables and their corresponding coefficients, respectively; ε is the error term assumed to have a particular parametric distribution; and ln(ε) is the natural log of the error term. Traditional regression and the AFT model differ in the following points: (1) predictor variables in the AFT model affect the event time multipli-catively; (2) the AFT model accommodates censored observations; and (3) the error terms in the AFT model, although still independent and identically distributed, no longer follow a normal distribution. Some commonly assumed parametric distributions in survival models include (see for instance Klein and Moeschberger15 for discussions regarding properties of these distributions): exponential, Weibull, generalized gamma, log-normal, and log-logistic. These are used in place of a normal distribution since the event times are positively valued and generally have a skewed distribution, making the symmetric normal distribution a poor choice for fitting the data closely.
In practice, the choice of which parametric distribution to use is done by comparing the model fit for a variety of different distributions. The choice of candidate distributions to be considered is driven by prior assumptions or scientific insight about the data. The comparison may be done graphically using probability plots which will show how the observed data follow an assumed parametric model. Distributions with multiple parameters defining their shape may have a better fit, but if parsimony is desired, it would be better to rely on a penalized metric provided by model selection indices such as the Akaike information criterion (AIC)16 or Bayesian information criterion (BIC)17 to choose which distribution gives the best fit with the fewest parameters among candidate distributions. These indices allow for numeric comparison which may be less subjective than comparing graphs. However, there is no formal statistical test associated with these indices. The choice of the distribution should not be based on which distribution gives a favorable P value. It should be noted that there is no distribution that provides a perfect fit, and it is possible that more than one distribution may fit the data well.
Once the distribution of the outcome has been decided, an investigator can focus on the effects of variables of interest on the time to an event. As previously noted, the effects of individual predictors in the AFT model are interpreted using time ratios (TR) where the ratio denotes the acceleration factor. Contrary to HR, a time ratio greater than one means that an event is less likely to occur as it means that an investigator must wait longer for the event to happen. Similarly, a time ratio less than one implies that the event is more likely to happen.
An important point to note is that when the survival distribution of the event of interest follows a Weibull distribution, the AFT model and the Cox proportional hazard model coincide.15 In other words, the AFT model assumes proportional hazard if the distribution is Weibull and vice versa. For all other parametric distributions, the AFT model assumes non-proportional hazards. This underlines the important distinction between the two models: for a given set of data, the AFT model and the Cox model (without covariates that vary with time) cannot both be correct, unless the survival distribution is Weibull.
A fictitious study enrolled a selected cohort of 200 patients with New York Heart Association (NYHA) Class II-III diastolic heart failure who were followed over time. Suppose that 100 of these patients have diabetes mellitus (DM), while the other 100 patients are non-diabetic (non-DM). Let the goal of the study be the comparison of cardiovascular-related mortality between diabetics and non-diabetics who all have NYHA Class II-III diastolic heart failure. Using the statistical package SAS version 9.3, data on time to death for both groups were artificially generated through simulation based on the mathematical model that generated the survival curves in Figure 1A, which are Weibull survival functions. Censoring was also artificially generated by assuming a maximum length of follow-up of 10 years and allowing for dropouts and loss to follow-up (Details of the simulation and analyses are in Appendix A.4).
In the simulated data, there were 38 deaths in the DM group and 22 deaths in the non-DM group. Figure 2 displays the Kaplan-Meier curve for the two groups and the P value of the log-rank test. Patients with diabetes have significantly lower survival than those without diabetes (P = .002). Kaplan-Meier curves do not go all the way down to zero when the largest observed time (which is around 9.5 years for both groups for this example) is censored.
Fitting a Cox model with only one predictor variable (i.e., presence of DM), a significant group difference (P = .003) was found just as in the log-rank test. However, in the Cox model, one can estimate the hazard ratio. In this case, the estimate of the model coefficient (β) for diabetes is 0.809 where non-DM is the reference group. Thus, the hazard ratio estimate is HR = e0.809 = 2.24 (95% confidence interval (CI): 1.3–3.8). This means that the DM group is estimated to have a hazard rate about twice that of the non-DM group. Figure 3 displays a graphical check for proportionality of hazards showing the transformed Kaplan-Meier (also known as the product limit estimator) curves. Since the curves in Figure 3 are approximately parallel to each other, there is no evidence of violation of the proportional hazard assumption from Figure 3.
Plot to check proportional hazards assumption of the Cox model.
In fitting an AFT model, it is worth noting that since the data were simulated, the survival distribution (Weibull in this case) is known. However, in practice, the true distribution of the event times is unknown. Therefore, in the analysis of data collected, it is recommended to fit several parametric distributions. In this case, the Weibull, log-normal, log-logistic, and Gamma distributions were fitted. Weibull was chosen to be the best fitting model using the AIC and BIC criteria (see Table 1). The fact that Weibull was the best fit model also supports the appropriateness of the proportional hazards assumption of the Cox model.
Model selection indices using several parametric distributions
As in the log-rank and Cox models, the Weibull AFT model with only DM as a predictor variable found significant group differences (P = .0034). The estimate for the model parameter associated with DM is −0.58 (95% CI: −0.96, −0.19) where non-DM is the reference group. Consequently, the estimate of the TR ratio is exp{−0.58} = 0.56 (95% CI: 0.38, 0.83), meaning DM shortened the survival time.
One may also prefer to provide estimates of the median time to death for each group. Table 2 displays the estimates of the median and its 95% CI for each group. Using this information, it is estimated that a patient from this artificially generated population with DM has a median time to death of 5.76 years (95% CI: 54.59–7.23). Note that the ratio of the median time-to-death estimates for patients with DM and non-DM is 0.56 (=5.76 ÷ 10.26), which was the TR obtained previously.
This paper reviews some basic concepts of survival analyses including discussions and comparisons between the semiparametric Cox proportional hazards model and the parametric AFT model. The appeal of the AFT model lies in the ease of interpreting the results, because the AFT models the effect of predictors and covariates directly on the survival time instead of through the hazard function. If the assumption of proportional hazards of the Cox model is met, the AFT model can be used with the Weibull distribution, while if proportional hazard is violated, the AFT model can be used with distributions other than Weibull.
It is essential to consider the model assumptions and recognize that if the assumptions are not met, the results may be erroneous or misleading. The AFT model assumes a certain parametric distribution for the failure times and that the effect of the covariates on the failure time is multiplicative. Several different distributions should be considered before choosing one. The Cox model assumes proportional hazards of the predictors over time. Model diagnostic tools and goodness of fit tests should be utilized to assess the model assumptions before statistical inferences are made.
In conclusion, although the Cox proportional hazards model tends to be more popular in the literature, the AFT model should also be considered when planning a survival analysis. It should go without saying that the choice should be driven by the desired outcome or the fit to the data, and never by which gives a significant P value for the predictor of interest. The choice should be dictated only by the research hypothesis and by which assumptions of the model are valid for the data being analyzed.
Glossary of statistical quantities used in this paper
Term
Symbol
Description
Survival function
S(t)
Proportion of subjects who are event-free at time t
Hazard function
h(t)
Instantaneous rate of experiencing an event, given the subject is event-free at time t
Predictor variable/covariates
x1, x2,…,xp
The P independent variables in a regression model such as age, sex, race, and treatment group
Predictor parameters
β1, β2,…,βp
The values corresponding to the P predictors that quantify how the predictor affects the outcome
Table 4
Glossary of statistical abbreviations used in this paper
Term
Definition
Explanation
AIC
BIC
Akaike’s information criteria
Bayesian information criteria
Indices that quantify how well the statistical model fits the data, with a penalty for added complexity (model parameters)
AFT
Accelerated failure time
Statistical model that can test the effects of multiple predictors on survival, controlling for the others. Assumes the event times follow a parametric distribution defined by the analyst
HR
Hazard ratio
The proportion of the hazard changes in the presence of a categorical predictor variable or from a one-unit increase in a continuous predictor. Calculated from the Cox model
TR
Time ratio
The proportion of the time-to-event changes in the presence of a categorical predictor variable or from a one-unit increase in a continuous predictor. Calculated from the AFT model
SAS
Statistical analysis system
Commercial software for statistical analysis. Developed by SAS Institute Inc., Cary, NC
Table 5
List of commonly used survival distributions
Name
# of parameters
Notes
Exponential
1
Assumes constant hazard
Special case of Weibull
Weibull
2
Assumes proportional hazards
Hazard is monotonically increasing or decreasing
Special case of gamma
Gamma
3
Hazard is monotonically increasing or decreasing
Log-normal
2
Hazard increases and later decreases
Log-logistic
2
Hazard is monotonically increasing or decreasing
Note that the number of parameters relates to the model complexity; if a simpler model is desired, a survival distribution with fewer parameters is preferred. Furthermore to a limited extent, additional parameters results in a greater fit to the data
1. Kaplan EK, Meier P. Nonparametric estimation from incomplete observations. J Am Stat Assoc. 1958;53(282):457–481.
2. Savage IR. Contributions to the theory of rank order statistics: The two sample case. Ann Math Stat. 1956;27(3):590–615.
3. Cox DR. Regression models and life-tables. J R Stat Soc Ser B (Methodol) 1972;34(2):187–220.
4. Xu Y-Z, Cha Y-M, Feng D, Powell BD, Wise HJ, Hua W, et al. Impact of myocardial scarring on outcomes of cardiac resyn-chronization therapy: Extent or location? J Nucl Med. 2012;53(1):47–54.[PubMed]
5. Bourque JM, Velazquez EJ, Tuttle RH, Shaw LK, O’Connor CM, Borges-Neto S. Mortality risk associated with ejection fraction differs among resting nuclear perfusion findings. J Nucl Cardiol. 2007;14(2):165–173. [PubMed]
6. Hachamovitch R, Berman DS. The use of nuclear cardiology in clinical decision making. Semin Nucl Med. 2005;35(1):62–72. [PubMed]
7. Nakata T, Miyamoto K, Doi A, Sasao H, Wakabayashi T, Ko-bayashi H, et al. Cardiac death prediction and impaired cardiac sympathetic innervation assessed by MIBG in patients with failing and nonfailing hearts. J Nucl Cardiol. 1998;5(6):579–590. [PubMed]
8. Duvall WL, Wijetunga MN, Klein TM, Razzouk L, Godbold J, Croft LB, et al. The prognosis of a normal stress-only Tc-99m myocardial perfusion imaging study. J Nucl Cardiol. 2010;17(3):370–377. [PubMed]
9. Hachamovitch R, Hayes S, Friedman JD, Cohen I, Shaw LJ, Germano G, et al. Determinants of risk and its temporal variation in patients with normal stress myocardial perfusion scans: What is the warranty period of a normal scan? J Am Coll Cardiol. 2003;41(8):1329–1340. [PubMed]
10. Acampa W, Evangelista L, Petretta M, Liuzzi R, Cuocolo A. Usefulness of stress cardiac single-photon emission computed tomographic imaging late after percutaneous coronary intervention for assessing cardiac events and time to such events. Am J Cardiol. 2007;100(3):436–441. [PubMed]
11. Petretta M, Acampa W, Evangelista L, Daniele S, Ferro A, Cuo-colo A. Impact of inducible ischemia by stress SPECT in cardiac risk assessment in diabetic patients: Rationale and design of a prospective multicenter trial. J Nucl Cardiol. 2008;15(1):100–104. [PubMed]
12. Daniele S, Nappi C, Acampa W, Storto G, Pellegrino T, Ricci F, et al. Incremental prognostic value of coronary flow reserve with single-photon emission computed tomography. J Nucl Cardiol. 2011;18(4):612–619. [PubMed]
13. Acampa W, Petretta M, Cuocolo R, Daniele S, Cantoni V, Cuo-colo A. Warranty period of normal stress myocardial perfusion imaging in diabetic patients: A propensity score analysis. J Nucl Cardiol. 2014;21(1):50–56. [PubMed]
14. Kalbfleisch JD, Prentice RL. The statistical analysis of failure time data. 2nd ed. Hoboken: John Wiley and Sons; 2002.
15. Klein JP, Moeschberger ML. Survival analysis: Techniques for censored and truncated data. New York: Springer; 2003.
16. Akaike H. A new look at the statistical model identification. IEEE Trans Autom Control. 1974;19(6):716–723.
17. Schwarz GE. Estimating the dimension of a model. Ann Stat. 1978;6(2):461–464
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570048/NWK_lam.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570068/NWExp.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570076/NWN_Exp.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570080/House1.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570082/Houses2.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570088/Houses3.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570108/Gender.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570110/Class.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570124/Fem1.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570132/Fem2.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570148/Male1.0.png)
/cdn0.vox-cdn.com/uploads/chorus_asset/file/3570154/Male2.0.png)
