Time-to-event (TTE) endpoints are fundamental in drug development and biomedical research. These endpoints are used to measure the duration until a specific event, such as disease progression, treatment failure, or death, which makes them critical to understanding treatment efficacy and patient prognosis. Analyzing TTE data requires specialized statistical and computational methods to account for censoring, where the event of interest may not be observed at the time of analysis. Traditional statistical models, such as the Cox proportional hazards model, have long been the cornerstone of the analysis of TTE outcomes. The Cox model [14] provides a semiparametric framework that estimates the effect of covariates on the hazard function while making minimal assumptions about the baseline hazard. Its interpretability and flexibility have made it a preferred choice in clinical and biomedical research. However, the Cox model assumes proportional hazards over time, which may not hold in real-world datasets, and it may face challenges when handling high-dimensional data or complex relationships, e.g., interactions among predictors.

In recent years, machine learning (ML) methods have gained attention for predicting TTE outcomes due to their flexibility and capacity to model nonlinear relationships and interactions among variables. Techniques such as tree-based models [31, 20] and deep learning approaches [34, 36] have been shown to achieve improved predictive performance in scenarios where traditional models fall short. ML methods can accommodate large-scale data, automatically select relevant features, and handle complex, high-dimensional predictors. These strengths have positioned ML as a promising alternative to classical statistical models in survival analysis.

Moreover, recent studies indicate the improved performance after incorporating post-baseline time-varying predictors into predictive models [20, 56, 48]. Post-baseline predictors, which evolve over time during the course of observation, can capture dynamic changes in patient status or in responses to treatment that may significantly influence TTE outcomes. Integrating these predictors into ML models has been reported to further enhance prediction accuracy, offering a more comprehensive and adaptive representation of survival dynamics.

Given these advances, a critical question remains: Which method performs better, traditional statistical models such as the Cox model, or flexible ML approaches, with or without post-baseline predictors in typical clinical trials with moderate sample size? Cuthbert et. al. [15] compared the prediction performance of survival analysis methods with only baseline predictors present. This paper aims to address this question through a systematic comparison of these approaches using multiple evaluation metrics, and to provide insights into their respective strengths and limitations, particularly in the context of incorporating post-baseline time-varying predictors, and to guide researchers in selecting the most appropriate methodology for their TTE analyses.

In this study, we first review a range of statistical and ML models for TTE outcomes, including the Cox proportional hazards model and tree-based models, with or without post-baseline predictors. Then, we examine several widely used model evaluation metrics for TTE outcomes, including the time-dependent area under the receiver operating characteristic (ROC) curve [AUC, 29], the concordance index [C-index, 26], and Brier score and integrated Brier score [9, 23, 58], which are commonly employed to assess predictive performance in survival analysis. In addition, we evaluate and compare the predictive performance of traditional statistical models, such as the Cox proportional hazards model and ML models for TTE outcomes. This evaluation considers scenarios with only baseline covariates and those incorporating post-baseline covariates. Lastly, we investigate the performance of various model evaluation metrics under different true data generation mechanisms. By comparing these metrics across different models, we aim to elucidate their strengths and limitations, offering a perspective on their capability to handle complex survival data.

While machine learning methods offer substantial modeling flexibility, their suitability depends critically on the data characteristics typical of clinical trials. Advanced approaches such as gradient boosting machines [18, 13] and deep learning–based survival models (e.g., DeepSurv [34]) have demonstrated strong performance in large-scale observational datasets. However, clinical trials—particularly Phase II and III oncology studies— usually have moderate sample sizes, where such methods may face challenges related to optimization stability, overfitting, and sensitivity to hyperparameter tuning. Empirical evidence suggests that deep learning and boosting-based methods often require substantially larger sample sizes to achieve stable and competitive performance in clinical tabular data [5, 44]. In contrast, random forest–based approaches have been shown to exhibit greater robustness under moderate sample sizes, with relatively low sensitivity to hyperparameter choices [42]. In addition, boosting algorithms may overemphasize noisy or mislabeled observations in small samples, increasing the risk of overfitting [16, 10, 19, 52], and tuning complexity [42, 6, 4]. Given these considerations, this study focuses on methods that balance nonlinear flexibility, robustness under moderate sample sizes, and feasibility for modeling post-baseline time-varying covariates—namely, the Cox proportional hazards model, tree-based methods for LTRC data [20], and survival random forest [31]. A more detailed discussion is provided in the Discussion section.

This paper is structured as follows. In Section 2, we describe the time-to-event (TTE) models, including the traditional Cox proportional hazards model, and modern ML methods tailored for TTE outcomes. These models incorporate baseline predictors and, where applicable, post-baseline predictors to enhance predictive performance. In Section 3, we explore model evaluation metrics specific to TTE variables. In Section 4, simulation studies are conducted in various settings to demonstrate the strengths and limitations of the TTE prediction models and their associated evaluation metrics. In Section 5, we conclude with a summary and a discussion.

Let X i = ( X i 1 , … , X i p ) ′ denote the p time-invariant covariates for subject i. The hazard function for the Cox proportional hazards model [14, 47] is given by: (2.1) λ ( t ∣ X i ) = λ 0 ( t ) exp X i ′ β , where λ 0 ( t ) is the baseline hazard function, β represents the regression coefficients for time-invariant covariates. The model assumes that the hazard ratios between different levels of covariates are constant over time, which means that the effect of the covariates on the hazard rate is multiplicative and does not change as time progresses. The model does not assume a specific form for the baseline hazard function, which allows flexibility in capturing time-dependent risk. The baseline hazard function can be estimated by the Breslow estimator [8]. The regression coefficients can be estimated by maximizing the partial likelihood independent of the baseline hazard [14, 47]. We use the R package survival to implement the fitting and prediction of the Cox proportional hazards model [49].

Ishwaran et al. [31] developed survival random forest models (SRFs), which offer a more flexible nonparametric alternative that can capture nonlinear relationships and interactions between covariates, often achieving superior predictive accuracy in complex data settings. SRF is an ensemble of multiple survival trees. Individual trees are grown using bootstrap samples and, at each node, a subset of variables is randomly selected to determine the optimal split. The nodes are divided according to survival differences using statistical measures such as the logrank test statistic [38]. Each tree produces a cumulative hazard function, which is averaged across all trees in the forest to produce ensemble predictions. In this study, we used random forest to model survival outcomes using only baseline covariates. The SRF is implemented in the R package randomForestSRC [32]. To optimize the performance of the model, we tune key parameters, including the number of trees in the forest (ntree), the number of variables considered for splitting at each node (mtry), and the minimum size of terminal nodes (nodesize). A five-fold cross-validation approach is used to determine the combination of tuning parameters that produces the locally optimal concordance index (C-index), a widely used evaluation metric for survival models.

When both baseline and post-baseline covariates are incorporated into the analysis, the modeling framework becomes more challenging due to the time-varying nature of some predictors. In the presence of time-varying predictors, let Z i j = ( Z i 1 j , … , Z i q j ) ′ denote the vector of q time-varying covariates for subject i observed at the j-th time interval. The hazard function for the Cox proportional hazards model can be conceptually extended to: (2.2) λ ( t ∣ X i , Z i j ) = λ 0 ( t ) exp X i ′ β + Z i j ′ γ , where γ denotes the coefficients for the q time-varying covariates. However, a key rule for time-dependent covariates is that they must not “look into the future”, or the hazard depends on the covariate values just prior to the event time [48]. An approach to handling time-dependent covariates is to use the methods for time intervals that account for left truncation and right censoring [LTRC 50, 48]. Thernea et al. [48] provide one example of coding time-dependent covariates in intervals of time. Imagine a subject whose follow-up period extends from time 0 to death at 185 days. Suppose a time-dependent covariate, such as the repeating laboratory test of creatinine in a clinical trial, is measured on day 0, day 90, and day 120, with recorded values of 0.9, 1.5 and 1.2 mg / dL, respectively. To represent these data in a suitable format for analysis, the follow-up time can be divided into three intervals: 0–90, 90–120, and 120–185 days. Each interval is represented as a separate row of data. The structured data appear in Table 1. The predictors of TTE in each interval are observed at the beginning of the interval, preventing the use of future data that would bias the results. Each row is also referred to as a pseudo-subject [20, 56], as each row represents a “subject” according to the interval time-to-event (TTE) models. However, the rows do not correspond to distinct true subjects in reality. Of note, the number of pseudo-subjects can be substantially larger than that of the true subjects, depending on the number of post-baseline observations for each true subject.

Additionally, ML methods such as survival trees for LTRC data [20] provide a nonparametric or semiparametric approach to handling TTE predictive modeling with post-baseline predictors. The survival tree model effectively segments pseudo-subjects into homogeneous subgroups according to time-invariant and time-dependent covariates. Two survival tree models were proposed including the LTRC tree based on conditional inference tree (LTRCIT) and LTRC tree based on classification and regression tree (LTRCART). The construction of LTRCIT is based on the logrank test score specifically adjusted for LTRC data, ensuring unbiased selection of splitting predictors. At each terminal node of the tree, the Kaplan-Meier estimate of the survival function is used to summarize survival times within that subgroup. In comparison, LTRCART is a likelihood-based method and uses deviation reduction and proportional hazards assumptions to select splits. The survival distribution at each terminal node of the tree is estimated through the baseline cumulative hazard function by the Nelson-Aalen estimator and the estimated relative risk of each node. Consequently, LTRCART is well suited for small datasets. In addition, practical experience suggests that the computational speed of modeling fitting and prediction for LTRCART is considerably faster than that of LTRCIT. Hence, in this paper, we use LTRCART as the ML model to predict TTE with predictors varying over time.

The LTRCART method is implemented in the R package LTRCtrees [21].The Key tuning parameters include the minimum sum of weights required in a terminal node (minbucket) and the maximum allowable depth of the tree (maxdepth), both of which influence the complexity and performance of the tree. To optimize these parameters, a five-fold cross-validation approach was implemented, with the C-index used as the evaluation metric to determine the local optimal settings. This methodology enables the survival tree to provide accurate and interpretable survival predictions while effectively addressing the challenges posed by the left-truncated and right-censored data.

In this paper, three metrics for evaluating survival models were investigated, including the C-index, integrated Brier score, and time-dependent AUC. They are briefly described in the following.

In this study, we systematically evaluated and compared the predictive performance of the Cox proportional hazards (PH) model and ML methods, including tree-based models, for analyzing time-to-event (TTE) outcomes. We investigated four distinct simulation settings that varied in complexity, incorporating both baseline and post-baseline time-varying covariates, as well as scenarios where no model was correctly specified. To comprehensively assess model performance, we utilized evaluation metrics such as the concordance index (C index), the integrated Brier score (IBS) and the time-dependent area under the curve (AUC).

The simulation studies highlight the importance of correct model specification, the strengths and limitations of ML approaches, and the role of sample size in model performance. In settings where the true model is correctly specified (e.g., Cox-post-baseline in Setting 1), it consistently outperforms alternatives, emphasizing the value of including both baseline and post-baseline covariates when relevant. ML models, such as SRF, excel in capturing complex nonlinear relationships and interactions, as seen in Setting 3, but require larger sample sizes to stabilize performance, particularly for tree-based models like LTRCART. When only baseline covariates are relevant (Setting 2), post-baseline predictors add no value, and simpler models perform comparably. However, evaluation metrics such as the time-dependent AUCs can fail to distinguish over-fitted models. In misspecified scenarios (Setting 4), the Cox model demonstrates robustness, performing similarly to ML models despite violating assumptions. The Cox model may serve as a valuable option, offering performance comparable to that of machine learning (ML) methods in practical applications, particularly when the sample size is moderate and/or the model assumptions, such as proportional hazards and linear effects of predictors on the log hazard ratio, are reasonably satisfied. In general, the studies underscore the utility of ML methods for complex data structures, the robustness of traditional models with mild misspecification, and the need for sufficient sample sizes for reliable performance in ML approaches.

Several studies in the literature have investigated model comparison criteria that include penalties for model complexity in the context of time-to-event (TTE) outcomes with baseline predictors only. Karabey and Tutkun [33] used the Akaike Information Criterion [AIC; 2] and the Bayesian Information Criterion [BIC; 43] to compare nested survival models. Habibi et al. [24] employed AIC to compare various survival models, including Exponential, Weibull, Gompertz, Log-normal, Log-logistic, and Generalized Gamma models. Ozaki and Ninomiy [39] utilized AIC to identify change-points in the Cox proportional hazards model. Similarly, Fagbamigb et al. [17] compared parametric and semi-parametric survival models using both AIC and BIC. However, these model evaluation methods require the specification of a likelihood function and may not be directly applicable to many machine learning models. More research is needed to develop model evaluation metrics that are more effective in detecting and addressing over-fitting in machine learning contexts and with post-baseline predictors.

For ML methods with post-baseline predictors, we only considered the tree model LTRCART and did not include ensemble methods. In the literature, Yao et al. [56] proposed ensemble approaches to estimate survival functions with time-varying covariates, based on conditional inference [51] and relative risk forests [30]. These methods are implemented in the R package LTRCforests [57]. However, as noted in this study, the number of pseudo-subjects can be substantially larger than the number of true subjects, leading to computational challenges. The computational demands of LTRCforests are significantly higher than those of SRF and LTRCART, exceeding our computational capacity for simulation studies. Furthermore, the large sample size requirement observed for LTRCART in simulation studies is likely applicable to LTRCforests as well. Further research in this direction is warranted to better understand the performance of ensemble tree methods for time-to-event outcomes with post-baseline predictors.

It is important to emphasize that our comparison was not intended to serve as an exhaustive benchmark of all machine learning approaches for survival prediction. Rather, we deliberately focused on methods that are most appropriate for clinical trial settings characterized by moderate sample sizes and post-baseline time-varying covariates. We did not include gradient boosting or deep learning–based survival models for three interrelated reasons. First, sample size considerations play a critical role. Deep learning models generally require substantially larger datasets to achieve stable optimization and outperform classical approaches. As demonstrated by Billichová et al. [5], DeepSurv required a sample size of approximately 6,000 to match the performance of the Cox model. Likewise, Silvey and Liu [44] showed that gradient boosting methods required considerably larger sample sizes than random forests to achieve stable AUC estimates in clinical tabular data. These requirements can exceed the typical sample sizes available in Phase II and III clinical trials. Second, hyperparameter robustness was a key consideration. Survival random forests were selected because random forest–based methods have been shown to exhibit low hyperparameter sensitivity, with default parameter settings often yielding near-optimal performance [42]. In contrast, gradient boosting methods rely on extensive tuning of multiple interdependent hyperparameters, including learning rate, tree depth, number of boosting iterations, and regularization parameters, to achieve optimal performance [42, 6]. In simulation studies, this sensitivity can introduce variability that reflects the efficiency of the tuning strategy rather than the intrinsic predictive capability of the method itself. Third, overfitting concerns are particularly relevant in smaller or noisy datasets. Boosting algorithms such as AdaBoost are known to overemphasize misclassified observations, which may include mislabeled or noisy data points, leading to excessive fitting of noise rather than underlying signal [16, 19]. Although regularization and early stopping can mitigate these effects, their successful application further depends on careful hyperparameter tuning, compounding the robustness issues discussed above [10, 52, 4]. For these reasons, we restricted our comparison to the Cox model, tree-based method for LTRC data and survival random forest, which are computationally feasible, robust to moderate sample sizes, and well suited for time-varying covariates. Future work leveraging large-scale real-world evidence databases may enable more comprehensive comparisons that include gradient boosting and deep learning methods under conditions where their advantages can be more fully realized.

Handling missing data is critical in data analysis and modeling, particularly in contexts like survival analysis with longitudinal post-baseline predictors. In the current paper, we only consider statistical and ML methods for datasets without missingness. Specifically, only pseudo-subjects with complete observations of baseline and post-baseline covariates consider are kept in the model. This approach allows us to focus on the modeling fitting and evaluation of TTE models. However, in practice, such an approach could lead to a significant reduction in sample size and potential bias when the missing mechanism is not trivial, for example, it is not random [37]. Methods for fitting Cox model with missing data have been proposed in the literature [12, 54]. SRF introduces a novel adaptive tree imputation algorithm to manage missing covariates and outcomes during the tree growth and prediction phases [31]. Evaluating TTE models with missing data will be considered in future studies. The R code supporting the computation in this paper is available at https://github.com/zhaohualu/SurvivalPredictiveModelEvaluation.