Regarding the training cohort, We selected 4517 medical cases from Database: Incidence - SEER Research Plus Data, 18 Registries, Nov 2019 Sub (2000–2017) - Linked To County Attributes - Total U.S., 1969-2018 Counties, National Cancer Institute,DCCPS, Surveillance Research Program, released April 2021, based on the November 2019 submission. We included Data records if they meet the criterion (1), patients pathologically diagnosed between January 2012 and December 2017 with primary stage-teriii non-small cell lung cancer (NSCLC) and (2) the existence of one malignant lesion. On the contrary, We excluded clinical cases according to the standard (1), patients whose regional lymph nodes performed during the initial work-up or first course of therapy are unknown or missing. Then we choose the features relevant to the OS (overall survival) of the NSCLC, including demographic information (Age and Sex) and NSCLC-cancer-related characteristics (TNM stage, histology type, primary site, tumor size, regional node number examined, regional node positive number and laterality of the tumor), and treatment details(surgery of primary site, radiation, and chemotherapy), The outcome is the patient survival time and death indicator. As for the cohort for external validation of the model, the inclusion criteria and exclusion criteria are consistent with the training group, So we randomly collected 100 stage-iii non-small cell lung cancer patients who underwent surgery (Lobectomy WITH mediastinal lymph node dissection and Pneumonectomy) from January 2012 to December 2017 in Shaanxi Provincial People’s Hospital, China.

The training data and the testing data are stored in CSV files. Both datasets contain two types of variables in the covariates, numerical variable and categorical variable. In the dataset, we have 3 numerical variable fields, including regional node positive number, regional node number examined and tumor size as well as other 10 categorical variable features. In order to avoid the evaluation problems via using label encoding conversion to categorical, we converse the 10 categorical features by utilizing one hot encoding to identify the different categorical values in the feature in a binary fashion. To illustrate, Regarding feature surgery on the primary site, before conversion, this field contains two values encoded for two surgery types (Lobectomy WITH mediastinal lymph node dissection, Pneumonectomy WITH mediastinal lymph node dissection). After transformation, the very field will be replaced by two surgery types, the value of the two features could only be 0 or 1 to identify the specific surgery type. In addition, as for feature tumor size, in the training set the unit is millimeter while in the testing dataset, the unit is centimeter. So we divide the value in the training set by 10 to make the unit the same. Finally, we perform normalization in order to accelerate the training process.

In this section, we created three machine learning models to perform the survival analysis to select the optimal one.

We developed a deep learning model based on DeepSurv to predict personal hazard rate according to the patient’s current clinical condition. From the input to output, the patient’s baseline data is the input to the neural network, followed by the fully-connected hidden layers of nodes as well as a drop layer after each hidden layer. The output of the network is the hazard rate. Regarding the activation function of each node, in order to overcome the problem of vanishing gradients, we select ReLU to add nonlinearity to the model which could help the model learn the complex relationship between covariates and the hazard rate. As for the loss function, we train the model to minimize the average negative log partial likelihood with regularization:

where θ is the weight of every node in the network, 1 N E = 1 is the number of dead patients and λ is the l ₂ regularization parameter, h ^ θ ( x ) is the predicated hazard rate. we use Adam for the gradient descent algorithm to update the parameter of the model for lots of epochs, because Adam is more efficient when working with problems involving high dimensional data and requiring less memory for optimization process compared with SGD method (18). We utilize random Search to optimize the hyper-parameters because compared to Grid Search, Random Search could try more cases for important hyper-parameters. In the experiment, we perform this on the log space of the learning rate in [0.00001, 0.1], the dropout rate in [0.2-0.5], the number of hidden layers in [1, 7] and the number of nodes in each hidden layer in [5,90].

We also trained a random forest model, this model is reliable because it forces each split to consider only a subset of the predictors. In this study, Random Search is still used to tune the number of the comprising trees in [100,300], the minimum number of samples required to split an internal node in [2,50] as well as the minimum number of samples required to be at a leaf node in [1,20].

Lastly, we trained the Penalized Cox Proportional hazard model with the same loss function as the deep learning model. we tuned the hyperparameter by using Random Search Method, specifically, the penalizer in [0.001,1] and the learning rate in [0.001,1].

The concordance index(C-index) is used to measure the performance of the model. The C-index is the ratio of pairs of patients ordered correctly to all pairs. Thus the higher C-index, the better performance of the model. In the study, The 4517 SEER data records were divided into two groups, 3534(80%) records were used for training while 883(20%) records were treated as the validation set. The five-fold cross-validation was performed to tune the hyper-parameters of each model and select the best model for survival prediction. Additionally, external validation was performed on the selected optimal model and TNM staging system and compare the generalizability of the two models. Eventually, we performed the attribution analysis for the deep learning survival model by the integrated gradients (19) method based on the testing dataset to rank the clinical feature importance.

The deep learning algorithm could recommend treatment for patients according to their current clinical conditions (20). we could load the model and set the input according to the patient’s demographic feature(age and gender), Surgery Type(lobectomy and pneumonectomy), Type(histology type and laterality) and the stage information of NSCLC(TNM, the number of the examined regional node, the number of the positive regional node and the tumor size). As for Adjuvant therapy, we predict the hazard rate under four adjuvant therapy treatments (with radiation and chemotherapy, with radiation and without chemotherapy, without radiation and with chemotherapy, and without radiation and chemotherapy). Then we could get the four cumulative hazard functions under each adjuvant therapy treatment and finally derive the four 5-year survival functions after negating and exponentiating the cumulative hazard function. In this application, we develop the backend code to calculate the four 5-year adjuvant therapy survival functions and implement the UI code to display the predicated survival functions in the line race chart.

The three models are trained with Python v 3.9, PyTorch v 1.11.0 is used to train the deep learning algorithm and PySurvival v 0.1.2 is utilized to train the random survival forest and penalized cox proportional hazard model. The Front UI of the adjuvant therapy recommender system is developed with Vue.js javascript framework and a Material Design component framework called Vuetify. The backend code of the web application is implemented by the Django REST framework. The recommender system is deployed on Tencent Cloud, which could be accessed through a web browser.

Based on the inclusion criteria, we include 4617 stage-iii NSCLC patients who received Surgeries (Lobectomy and Pneumonectomy with mediastinal lymph node dissection) in this study. The 4517 patients out of 4617 are extracted from the SEER database and used as a training set while the other 100 patients are from China Database for model testing. The baseline medical characteristics of the two cohorts are shown in Table 1 . From the AJCC TNM staging system’s perspective, all patients in the training set and the testing set are stage-iii NSCLC patients. In the SEER cohort, most patients’ histology type is Adenocarcinoma, which takes 44.28%. The next one is Squamous cell carcinoma, which takes 23.27%. Regarding the Received surgeries, 85.51% of patients received Lobectomy WITH mediastinal lymph node dissection while the rest (14.48%) accepted Pneumonectomy WITH mediastinal lymph node dissection for treatment. Concerning Adjuvant treatment, 74.12% of patients accepted chemotherapy and about 41.88% received beam radiation. On the contrary, in the test cohort, most patients received Lobectomy WITH mediastinal lymph node dissection, the two leading histology types are Squamous cell carcinoma and Adenocarcinoma, respectively 46% and 43% of the population. As for Adjuvant treatment, 1/3 received beam radiation and almost everyone received chemotherapy.

After the process of random search, we finally settled down on the hyperparameter of the deep learning model, the model consists of 2 hidden layers, from input to output, including 60, 43 neurons in each layer with a dropping out unit between each layer. we improve neural network generalization by setting the learning rate to 0.001 and 0.5 as the dropout rate to avoid overfitting. Figure 2 shows the training loss curves of the survival network. At the beginning of the training process, the loss of the validation and training set decreases continually. After 331 epochs of parameter optimization, the loss of the validation set begins at 3.6936 and stops decreasing at 3.1753 while the training loss continues to decrease from 3.3844 started at 3.8446. Then we terminate the optimization to avoid overfitting and save the model for test.

In the random survival forest, We set the number of the estimating trees to 959, the minimum number of samples required to split an internal node to 10 and the minimum number of samples required to be at a leaf node to 15. In the Penalized Cox Proportional hazard model, we configure the penalizer to 0.005 and the learning rate to 0.01

Then we perform 5-fold cross-validation to select the optimal model for survival prediction. Figure 3 displays the exact value and the line chart of each model in every fold validation, the deep learning model shows a more stable and exceptional performance on the concordance index compared to the other two models. The mean of the concordance index of the deep learning algorithm is 0.843, which is much higher than the random forest (0.678) and cox proportional hazard model (0.678) ( Table 2 ). Based on the result of cross-validation, deep learning is selected to compare the TNM staging system on external validation. The performance of the deep learning model is better (0.82 vs 0.65)

As for the feature importance for the network, from the Figure 4 we can observe four of the top important features: regional positive nodes (0.6634), regional examied nodes (-0.7648), tumor size (-0.5633) and Age(-0.4633). In terms of least important features, we observe that the surgery on the primary site (0.0632) is voted to be least significant based on attribution algorithm. The absolute value for attribution scores of other features is greater than 0.1 and less than 0.5.

Then we perform 5-fold cross-validation to select the optimal model for survival prediction. Figure 3 displays the exact value and the line chart of each model in every fold validation, the deep learning model shows a more stable and exceptional performance on the concordance index compared to the other two models. The mean of the concordance index of the deep learning algorithm is 0.843, which is much higher than the random forest (0.678) and cox proportional hazard model (0.678) ( Table 2 ). Based on the result of cross-validation, deep learning is selected to compare the TNM staging system on external validation. The performance of the deep learning model is better (0.82 vs 0.65)

Since the deep learning model has better performance than the TNM staging system, we could not only predict the survival function of the current patient but also offer an adjuvant therapy reference to the oncology doctor based on prediction over different therapy treatment plans. Thus we deployed the recommender system to the Internet, which could be accessed with a browser in [http://1.15.80.136/nsclc/], input the current clinical status, including Demographic, surgery type, cancer type and stage information, of one patient, and click the submit button ( Figure 5 ).

Then the browser will redirect to the result page ( Figure 6 ), and we could see four 5-year predicted survival curves for each treatment plan. Based on the plot, the predicted optimal treatment plan is only receiving beam radiation for adjuvant treatment, whose survival probability is highest in the next 60 months.

Thus, the specialist could get the reference for adjuvant treatment plan decision-making. Code related to this application can be found at https://github.com/snowflake-Zhao/nsclc.