Estimating cumulative treatment effects in the presence of nonproportional hazards

Summary .   Often in medical studies of time to an event, the treatment effect is not constant over time. In the context of Cox regression modeling, the most frequent solution is to apply a model that assumes the treatment effect is either piecewise constant or varies smoothly over time, i.e., the Cox nonproportional hazards model. This approach has at least two major limitations. First, it is generally difficult to assess whether the parametric form chosen for the treatment effect is correct. Second, in the presence of nonproportional hazards, investigators are usually more interested in the cumulative than the instantaneous treatment effect (e.g., determining if and when the survival functions cross). Therefore, we propose an estimator for the aggregate treatment effect in the presence of nonproportional hazards. Our estimator is based on the treatment-specific baseline cumulative hazards estimated under a stratified Cox model. No functional form for the nonproportionality need be assumed. Asymptotic properties of the proposed estimators are derived, and the finite-sample properties are assessed in simulation studies. Pointwise and simultaneous confidence bands of the estimator can be computed. The proposed method is applied to data from a national organ failure registry.     

Bayesian cure rate frailty models with application to a root canal therapy study

Due to natural or artificial clustering, multivariate survival data often arise in biomedical studies, for example, a dental study involving multiple teeth from each subject. A certain proportion of subjects in the population who are not expected to experience the event of interest are considered to be "cured" or insusceptible. To model correlated or clustered failure time data incorporating a surviving fraction, we propose two forms of cure rate frailty models. One model naturally introduces frailty based on biological considerations while the other is motivated from the Cox proportional hazards frailty model. We formulate the likelihood functions based on piecewise constant hazards and derive the full conditional distributions for Gibbs sampling in the Bayesian paradigm. As opposed to the Cox frailty model, the proposed methods demonstrate great potential in modeling multivariate survival data with a cure fraction. We illustrate the cure rate frailty models with a root canal therapy data set.     

Hierarchical proportional hazards regression models for highly stratified data

In clinical trials conducted over several data collection centers, the most common statistically defensible analytic method, a stratified Cox model analysis, suffers from two important defects. First, identification of units that are outlying with respect to the baseline hazard is awkward since this hazard is implicit (rather than explicit) in the Cox partial likelihood. Second (and more seriously), identification of modest treatment effects is often difficult since the model fails to acknowledge any similarity across the strata. We consider a number of hierarchical modeling approaches that preserve the integrity of the stratified design while offering a middle ground between traditional stratified and unstratified analyses. We investigate both fully parametric (Weibull) and semiparametric models, the latter based not on the Cox model but on an extension of an idea by Gelfand and Mallick (1995, Biometrics 51, 843-852), which models the integrated baseline hazard as a mixture of monotone functions. We illustrate the methods using data from a recent multicenter AIDS clinical trial, comparing their ease of use, interpretation, and degree of robustness with respect to estimates of both the unit-specific baseline hazards and the treatment effect.     

Analysis of the survival rate in a combined radiation lesion. The evaluation of synergism (or antagonism) in the case of 2 exposures

The method of the estimation of interaction between two harmful effects based on the statistical survival rate data is discussed. In terms of the competitive risk model, an independent effect was defined as a product of survival functions corresponding to an isolated effect of each injurious factor. The proposed method is based on the comparison between the animals' lifetimes in the case of the isolated effects and the actual survival rates after the combined radiation effect. The comparison of the two alternatives was made by the regression model of Cox adapted for the purposes of this study.     

Correlated versus uncorrelated frailty Cox models: A comparison of different estimation procedures

In many studies in medicine, including clinical trials and epidemiological investigations, data are clustered into groups such as health centers or herds in veterinary medicine. Such data are usually analyzed by hierarchical regression models to account for possible variation between groups. When such variation is large, it is of potential interest to explore whether additionally the effect of a within‐group predictor varies between groups. In survival analysis, this may be investigated by including two frailty terms at group level in a Cox proportional hazards model. Several estimation methods have been proposed to estimate this type of frailty Cox models. We review four of these methods, apply them to real data from veterinary medicine, and compare them using a simulation study.     

Semiparametric Contrasts of Cumulative Pre-Treatment Mortality in the Presence of Dependent Censoring

In clinical settings, the necessity of treatment is often measured in terms of the patient’s prognosis in the absence of treatment. Along these lines, it is often of interest to compare subgroups of patients (e.g., based on underlying diagnosis) with respect to pre-treatment survival. Such comparisons may be complicated by at least two important issues. First, mortality contrasts by subgroup may differ over follow-up time, as opposed to being constant, and may follow a form that is difficult to model parametrically. Moreover, in settings where the proportional hazards assumption fails, investigators tend to be more interested in cumulative (as opposed to instantaneous) effects on mortality. Second, pre-treatment death is censored by the receipt of treatment and in settings where treatment assignment depends on time-dependent factors that also affect mortality, such censoring is likely to be informative. We propose semiparametric methods for contrasting subgroup-specific cumulative mortality in the presence of dependent censoring. The proposed estimators are based on the cumulative hazard function, with pre-treatment mortality assumed to follow a stratified Cox model. No functional form is assumed for the nature of the non-proportionality. Asymptotic properties of the proposed estimators are derived, and simulation studies show that the proposed methods are applicable to practical sample sizes. The methods are then applied to contrast pre-transplant mortality for acute versus chronic End-Stage Liver Disease patients.     

Joint Analysis of Survival Time and Longitudinal Categorical Outcomes

In biomedical or public health research, it is common for both survival time and longitudinal categorical outcomes to be collected for a subject, along with the subject’s characteristics or risk factors. Investigators are often interested in finding important variables for predicting both survival time and longitudinal outcomes which could be correlated within the same subject. Existing approaches for such joint analyses deal with continuous longitudinal outcomes. New statistical methods need to be developed for categorical longitudinal outcomes. We propose to simultaneously model the survival time with a stratified Cox proportional hazards model and the longitudinal categorical outcomes with a generalized linear mixed model. Random effects are introduced to account for the dependence between survival time and longitudinal outcomes due to unobserved factors. The Expectation–Maximization (EM) algorithm is used to derive the point estimates for the model parameters, and the observed information matrix is adopted to estimate their asymptotic variances. Asymptotic properties for our proposed maximum likelihood estimators are established using the theory of empirical processes. The method is demonstrated to perform well in finite samples via simulation studies. We illustrate our approach with data from the Carolina Head and Neck Cancer Study (CHANCE) and compare the results based on our simultaneous analysis and the separately conducted analyses using the generalized linear mixed model and the Cox proportional hazards model. Our proposed method identifies more predictors than by separate analyses.     

The Mizon-Richard encompassing test for the Cox and Aalen additive hazards models.

The Cox hazards model (Cox, 1972, Journal of the Royal Statistical Society, Series B34, 187-220) for survival data is routinely used in many applied fields, sometimes, however, with too little emphasis on the fit of the model. A useful alternative to the Cox model is the Aalen additive hazards model (Aalen, 1980, in Lecture Notes in Statistics-2, 1-25) that can easily accommodate time changing covariate effects. It is of interest to decide which of the two models that are most appropriate to apply in a given application. This is a nontrivial problem as these two classes of models are nonnested except only for special cases. In this article we explore the Mizon-Richard encompassing test for this particular problem. It turns out that it corresponds to fitting of the Aalen model to the martingale residuals obtained from the Cox regression analysis. We also consider a variant of this method, which relates to the proportional excess model (Martinussen and Scheike, 2002, Biometrika 89, 283-298). Large sample properties of the suggested methods under the two rival models are derived. The finite-sample properties of the proposed procedures are assessed through a simulation study. The methods are further applied to the well-known primary biliary cirrhosis data set.     

The Mizon–Richard Encompassing Test for the Cox and Aalen Additive Hazards Models

Summary .   The Cox hazards model ( Cox, 1972 , Journal of the Royal Statistical Society, Series B 34, 187–220) for survival data is routinely used in many applied fields, sometimes, however, with too little emphasis on the fit of the model. A useful alternative to the Cox model is the Aalen additive hazards model ( Aalen, 1980 , in Lecture Notes in Statistics-2 , 1–25) that can easily accommodate time changing covariate effects. It is of interest to decide which of the two models that are most appropriate to apply in a given application. This is a nontrivial problem as these two classes of models are nonnested except only for special cases. In this article we explore the Mizon–Richard encompassing test for this particular problem. It turns out that it corresponds to fitting of the Aalen model to the martingale residuals obtained from the Cox regression analysis. We also consider a variant of this method, which relates to the proportional excess model ( Martinussen and Scheike, 2002 , Biometrika 89, 283–298). Large sample properties of the suggested methods under the two rival models are derived. The finite-sample properties of the proposed procedures are assessed through a simulation study. The methods are further applied to the well-known primary biliary cirrhosis data set.     

Empirical likelihood for cumulative hazard ratio estimation with covariate adjustment

In medical studies, it is often of scientific interest to evaluate the treatment effect via the ratio of cumulative hazards, especially when those hazards may be nonproportional. To deal with nonproportionality in the Cox regression model, investigators usually assume that the treatment effect has some functional form. However, to do so may create a model misspecification problem because it is generally difficult to justify the specific parametric form chosen for the treatment effect. In this article, we employ empirical likelihood (EL) to develop a nonparametric estimator of the cumulative hazard ratio with covariate adjustment under two nonproportional hazard models, one that is stratified, as well as a less restrictive framework involving group-specific treatment adjustment. The asymptotic properties of the EL ratio statistic are derived in each situation and the finite-sample properties of EL-based estimators are assessed via simulation studies. Simultaneous confidence bands for all values of the adjusted cumulative hazard ratio in a fixed interval of interest are also developed. The proposed methods are illustrated using two different datasets concerning the survival experience of patients with non-Hodgkin's lymphoma or ovarian cancer.     

Statistical Methods for Analyzing Right‐Censored Length‐Biased Data under Cox Model

Summary Length‐biased time‐to‐event data are commonly encountered in applications ranging from epidemiological cohort studies or cancer prevention trials to studies of labor economy. A longstanding statistical problem is how to assess the association of risk factors with survival in the target population given the observed length‐biased data. In this article, we demonstrate how to estimate these effects under the semiparametric Cox proportional hazards model. The structure of the Cox model is changed under length‐biased sampling in general. Although the existing partial likelihood approach for left‐truncated data can be used to estimate covariate effects, it may not be efficient for analyzing length‐biased data. We propose two estimating equation approaches for estimating the covariate coefficients under the Cox model. We use the modern stochastic process and martingale theory to develop the asymptotic properties of the estimators. We evaluate the empirical performance and efficiency of the two methods through extensive simulation studies. We use data from a dementia study to illustrate the proposed methodology, and demonstrate the computational algorithms for point estimates, which can be directly linked to the existing functions in S‐PLUS or R .     

Heterogeneous treatment effects in stratified clinical trials with time‐to‐event endpoints

When analyzing clinical trials with a stratified population, homogeneity of treatment effects is a common assumption in survival analysis. However, in the context of recent developments in clinical trial design, which aim to test multiple targeted therapies in corresponding subpopulations simultaneously, the assumption that there is no treatment‐by‐stratum interaction seems inappropriate. It becomes an issue if the expected sample size of the strata makes it unfeasible to analyze the trial arms individually. Alternatively, one might choose as primary aim to prove efficacy of the overall (targeted) treatment strategy. When testing for the overall treatment effect, a violation of the no‐interaction assumption renders it necessary to deviate from standard methods that rely on this assumption. We investigate the performance of different methods for sample size calculation and data analysis under heterogeneous treatment effects. The commonly used sample size formula by Schoenfeld is compared to another formula by Lachin and Foulkes, and to an extension of Schoenfeld's formula allowing for stratification. Beyond the widely used (stratified) Cox model, we explore the lognormal shared frailty model, and a two‐step analysis approach as potential alternatives that attempt to adjust for interstrata heterogeneity. We carry out a simulation study for a trial with three strata and violations of the no‐interaction assumption. The extension of Schoenfeld's formula to heterogeneous strata effects provides the most reliable sample size with respect to desired versus actual power. The two‐step analysis and frailty model prove to be more robust against loss of power caused by heterogeneous treatment effects than the stratified Cox model and should be preferred in such situations.     

Clustered restricted mean survival time regression

For multicenter randomized trials or multilevel observational studies, the Cox regression model has long been the primary approach to study the effects of covariates on time-to-event outcomes. A critical assumption of the Cox model is the proportionality of the hazard functions for modeled covariates, violations of which can result in ambiguous interpretations of the hazard ratio estimates. To address this issue, the restricted mean survival time (RMST), defined as the mean survival time up to a fixed time in a target population, has been recommended as a model-free target parameter. In this article, we generalize the RMST regression model to clustered data by directly modeling the RMST as a continuous function of restriction times with covariates while properly accounting for within-cluster correlations to achieve valid inference. The proposed method estimates regression coefficients via weighted generalized estimating equations, coupled with a cluster-robust sandwich variance estimator to achieve asymptotically valid inference with a sufficient number of clusters. In small-sample scenarios where a limited number of clusters are available, however, the proposed sandwich variance estimator can exhibit negative bias in capturing the variability of regression coefficient estimates. To overcome this limitation, we further propose and examine bias-corrected sandwich variance estimators to reduce the negative bias of the cluster-robust sandwich variance estimator. We study the finite-sample operating characteristics of proposed methods through simulations and reanalyze two multicenter randomized trials.     

Survival of Patients with Primary Brain Tumors: Comparison of Two Statistical Approaches

Purpose

We reviewed the survival time for patients with primary brain tumors undergoing treatment with stereotactic radiation methods at the Masaryk Memorial Cancer Institute Brno. We also identified risk factors and characteristics, and described their influence on survival time.

Methods

In summarizing survival data, there are two functions of principal interest, namely, the survival function and the hazard function. In practice, both of them can depend on some characteristics. We focused on nonparametric methods, propose a method based on kernel smoothing, and compared our estimates with the results of the Cox regression model. The hazard function is conditional to age and gross tumor volume and visualized as a color-coded surface. A multivariate Cox model was also designed.

Results

There were 88 patients with primary brain cancer, treated with stereotactic radiation. The median survival of our patient cohort was 47.8 months. The estimate of the hazard function has two peaks (about 10 months and about 40 months). The survival time of patients was significantly different for various diagnoses (p≪0.001), KI (p = 0.047) and stereotactic methods (p = 0.033). Patients with a greater GTV had higher risk of death. The suitable threshold for GTV is 20 cm³. Younger patients with a survival time of about 50 months had a higher risk of death. In the multivariate Cox regression model, the selected variables were age, GTV, sex, diagnosis, KI, location, and some of their interactions.

Conclusion

Kernel methods give us the possibility to evaluate continuous risk variables and based on the results offer risk-prone patients a different treatment, and can be useful for verifying assumptions of the Cox model or for finding thresholds of continuous variables.     

Confidence bands for multiplicative hazards models: Flexible resampling approaches

We propose new resampling‐based approaches to construct asymptotically valid time‐simultaneous confidence bands for cumulative hazard functions in multistate Cox models. In particular, we exemplify the methodology in detail for the simple Cox model with time‐dependent covariates, where the data may be subject to independent right‐censoring or left‐truncation. We use simulations to investigate their finite sample behavior. Finally, the methods are utilized to analyze two empirical examples with survival and competing risks data.     

Finding Risk Groups by Optimizing Artificial Neural Networks on the Area under the Survival Curve Using Genetic Algorithms

We investigate a new method to place patients into risk groups in censored survival data. Properties such as median survival time, and end survival rate, are implicitly improved by optimizing the area under the survival curve. Artificial neural networks (ANN) are trained to either maximize or minimize this area using a genetic algorithm, and combined into an ensemble to predict one of low, intermediate, or high risk groups. Estimated patient risk can influence treatment choices, and is important for study stratification. A common approach is to sort the patients according to a prognostic index and then group them along the quartile limits. The Cox proportional hazards model (Cox) is one example of this approach. Another method of doing risk grouping is recursive partitioning (Rpart), which constructs a decision tree where each branch point maximizes the statistical separation between the groups. ANN, Cox, and Rpart are compared on five publicly available data sets with varying properties. Cross-validation, as well as separate test sets, are used to validate the models. Results on the test sets show comparable performance, except for the smallest data set where Rpart’s predicted risk groups turn out to be inverted, an example of crossing survival curves. Cross-validation shows that all three models exhibit crossing of some survival curves on this small data set but that the ANN model manages the best separation of groups in terms of median survival time before such crossings. The conclusion is that optimizing the area under the survival curve is a viable approach to identify risk groups. Training ANNs to optimize this area combines two key strengths from both prognostic indices and Rpart. First, a desired minimum group size can be specified, as for a prognostic index. Second, the ability to utilize non-linear effects among the covariates, which Rpart is also able to do.     

A Flexible Alternative to the Cox Proportional Hazards Model for Assessing the Prognostic Accuracy of Hospice Patient Survival

Prognostic models are often used to estimate the length of patient survival. The Cox proportional hazards model has traditionally been applied to assess the accuracy of prognostic models. However, it may be suboptimal due to the inflexibility to model the baseline survival function and when the proportional hazards assumption is violated. The aim of this study was to use internal validation to compare the predictive power of a flexible Royston-Parmar family of survival functions with the Cox proportional hazards model. We applied the Palliative Performance Scale on a dataset of 590 hospice patients at the time of hospice admission. The retrospective data were obtained from the Lifepath Hospice and Palliative Care center in Hillsborough County, Florida, USA. The criteria used to evaluate and compare the models'' predictive performance were the explained variation statistic R², scaled Brier score, and the discrimination slope. The explained variation statistic demonstrated that overall the Royston-Parmar family of survival functions provided a better fit (R² = 0.298; 95% CI: 0.236–0.358) than the Cox model (R² = 0.156; 95% CI: 0.111–0.203). The scaled Brier scores and discrimination slopes were consistently higher under the Royston-Parmar model. Researchers involved in prognosticating patient survival are encouraged to consider the Royston-Parmar model as an alternative to Cox.     

Analysis of failure time data with multilevel clustering, with application to the child vitamin a intervention trial in Nepal

We consider the problem of estimating covariate effects in the marginal Cox proportional hazard model and multilevel associations for child mortality data collected from a vitamin A supplementation trial in Nepal, where the data are clustered within households and villages. For this purpose, a class of multivariate survival models that can be represented by a functional of marginal survival functions and accounts for hierarchical structure of clustering is exploited. Based on this class of models, an estimation strategy involving a within-cluster resampling procedure is proposed, and a model assessment approach is presented. The asymptotic theory for the proposed estimators and lack-of-fit test is established. The simulation study shows that the estimates are approximately unbiased, and the proposed test statistic is conservative under extremely heavy censoring but approaches the size otherwise. The analysis of the Nepal study data shows that the association of mortality is much greater within households than within villages.     

Targeted maximum likelihood estimation of effect modification parameters in survival analysis

The Cox proportional hazards model or its discrete time analogue, the logistic failure time model, posit highly restrictive parametric models and attempt to estimate parameters which are specific to the model proposed. These methods are typically implemented when assessing effect modification in survival analyses despite their flaws. The targeted maximum likelihood estimation (TMLE) methodology is more robust than the methods typically implemented and allows practitioners to estimate parameters that directly answer the question of interest. TMLE will be used in this paper to estimate two newly proposed parameters of interest that quantify effect modification in the time to event setting. These methods are then applied to the Tshepo study to assess if either gender or baseline CD4 level modify the effect of two cART therapies of interest, efavirenz (EFV) and nevirapine (NVP), on the progression of HIV. The results show that women tend to have more favorable outcomes using EFV while males tend to have more favorable outcomes with NVP. Furthermore, EFV tends to be favorable compared to NVP for individuals at high CD4 levels.     

Predicting overall survival of patients with hepatocellular carcinoma using a three‐category method based on DNA methylation and machine learning

Hepatocellular carcinoma (HCC) is closely associated with abnormal DNA methylation. In this study, we analyzed 450K methylation chip data from 377 HCC samples and 50 adjacent normal samples in the TCGA database. We screened 47,099 differentially methylated sites using Cox regression as well as SVM‐RFE and FW‐SVM algorithms, and constructed a model using three risk categories to predict the overall survival based on 134 methylation sites. The model showed a 10‐fold cross‐validation score of 0.95 and satisfactory predictive power, and correctly classified 26 of 33 samples in testing set obtained by stratified sampling from high, intermediate and low risk groups.