Radiomics and deep learning predict neoadjuvant immunochemotherapy response in locally advanced esophageal squamous cell carcinoma

Abstract

BACKGROUND

Neoadjuvant immunochemotherapy (NICT) offers a promising treatment for locally advanced esophageal squamous cell carcinoma (ESCC), yet predicting which patients will achieve a pathological response remains challenging caused by the lack of reliable noninvasive biomarkers.

AIM

To construct and assess multimodal models integrating radiomics, deep learning (DL), and clinical data to predict NICT response in locally advanced ESCC.

METHODS

188 patients who received NICT and radical surgery were enrolled in this retrospective multicenter study. Clinical variables, radiomic features from contrast-enhanced computed tomography, and DL features from three pre-trained architectures (Densenet201, Resnet101, Vgg11) were collected. These features were used to build clinical, radiomics, radiomics-DL, and clinical-radiomics-DL models using four machine learning algorithms. End-to-end DL models were also developed. Model performance was evaluated using area under the curve (AUC), calibration, and decision curve analysis.

RESULTS

In the external validation cohort, single-modality models (clinical, radiomics, and DL) showed limited predictive performance. The radiomics-DL model constructed with RandomForest (RF) with Vgg11 features achieved the best overall performance, with AUCs of 0.998, 0.855, and 0.779 in the training, internal validation, and external validation cohorts, respectively, and demonstrated good calibration and consistent clinical net benefit. The clinical-radiomics-DL model did not further improve predictive performance.

CONCLUSION

The radiomics-DL model (Vgg11-RF) predicts NICT response in locally advanced ESCC with promising generalizability, offering a noninvasive tool for individualized treatment decisions.

Key Words: Esophageal squamous cell carcinoma; Neoadjuvant immunochemotherapy; Radiomics; Deep learning; Treatment response prediction

Core Tip: Accurate prediction of response to neoadjuvant immunochemotherapy in locally advanced esophageal squamous cell carcinoma remains challenging. This multicenter study developed and compared clinical, radiomics, deep learning (DL), and multimodal models using pretreatment contrast-enhanced computed tomography images from 188 patients. The radiomics-DL model combining handcrafted radiomics features with Vgg11-derived DL features using RandomForest achieved the best performance, with areas under the curve of 0.998, 0.855, and 0.779 in training, internal validation, and external validation cohorts, respectively. This model demonstrated favorable generalizability across centers and offers a noninvasive tool for individualized treatment decision-making.

INTRODUCTION

Globally, esophageal cancer ranks as the sixth most frequent cause of cancer-related death[1]. In China, > 90% are esophageal squamous cell carcinoma (ESCC)[2,3]. ESCC is frequently diagnosed at locally advanced stage because of absence of specific early symptoms and its aggressive biology.

Neoadjuvant immunochemotherapy (NICT) is increasingly used for locally advanced resectable ESCC due to the widespread adoption of immune checkpoint inhibitors[4-6]. However, clinical responses to NICT are heterogeneous: Only about 30% of patients achieve a pathological complete response[7], while the non-responders risk disease progression that may preclude curative surgery. Although various biomarkers have been evaluated to predict NICT response, none has yet demonstrated sufficient reliability for routine clinical application. Tissue-based biomarkers, including microsatellite instability status, tumor mutational burden, and programmed cell death protein ligand 1 expression, have been the most extensively investigated; however, their predictive value remains controversial, with conflicting results across studies[8]. Moreover, these biomarkers rely on invasive tissue acquisition and are affected by both spatial and temporal heterogeneity, limiting their practical utility. On the other hand, conventional computed tomography (CT)/magnetic resonance imaging (MRI) features (e.g., tumor shrinkage, necrosis) are highly subjective, heavily dependent on physician expertise, and difficult to standardize. Thus, accurate predictive instruments are urgently required to identify which patients are most likely to respond to NICT among those with locally advanced resectable ESCC.

Radiomics employs computational methods to extract high-throughput quantitative features from medical images, revealing imaging information not readily appreciable by the human eye[9-11]. Deep learning (DL), unlike conventional radiomics, uses nonlinear neural networks to automatically learn image features related to texture and morphology. Prior work in other malignancies has indicated that multimodal fusion models that integrate conventional radiomics, DL features, and clinical data often outperform single-modality models[12,13]. The use of artificial intelligence (AI) in esophageal cancer imaging has evolved along multiple directions, broadly encompassing computer-aided detection/diagnosis and treatment response prediction. For instance, Weng et al[14] and Tsai et al[15] recently developed DL-based frameworks to enhance endoscopic detection of esophageal cancer, demonstrating the potential of AI to improve diagnostic sensitivity and reduce missed lesions. However, beyond diagnostic assistance, a more therapeutically pertinent application of AI lies in predicting tumor response to neoadjuvant therapy, which could carry direct implications for treatment selection and patient management. Recently, radiomics studies for NICT response prediction have emerged[16-20]. Despite these advances, several challenges remain. First, most existing studies remain limited by relatively small sample sizes and a singlecenter design. Second, the optimal strategy for integrating radiomics, DL, and clinical features remains underexplored, with inconsistent findings across studies regarding whether adding clinical data improves or compromises model generalizability.

This study aimed to construct and validate models to predict the NICT response in resectable, locally advanced ESCC patients using pretreatment contrast-enhanced chest CT. We constructed clinical, conventional radiomics, DL, and several fusion models, with the goal of generating evidence to advance personalized treatment strategies and support clinical decision-making.

MATERIALS AND METHODS

Patients

Approval for this retrospective, two-center study was obtained from the Ethics Committee of the Affiliated Hospital of Qingdao University (Approval No. QYFY WZLL 30959) and the Research Ethics Committee of the Affiliated Hospital of Jining Medical University (Approval No. 2026-01-C021), and the study was carried out in accordance with the Declaration of Helsinki. Written informed consent was not required, given the retrospective nature of the research.

A total of 154 patients presenting with resectable, locally advanced ESCC and undergoing radical surgery after NICT between March 2021 and June 2025 were included in this study. They radical surgery following NICT at the Affiliated Hospital of Qingdao University between March 2021 and June 2025 were enrolled. They were randomly assigned at a 7:3 ratio to a training set (n = 107) and an internal validation set (n = 47). A separate external validation set consisted of 34 patients who had received identical regimen at the Affiliated Hospital of Jining Medical University between April 2023 and June 2025.

Inclusion criteria were: (1) Pathologically confirmed ESCC by endoscopic biopsy; (2) Clinical stage cT1b-4bN+M0 or cT3-4bN0M0 in accordance with the American Joint Committee on Cancer (AJCC) staging system; (3) Treatment consisting of NICT followed by radical esophagectomy; and (4) Contrast-enhanced chest CT performed within two weeks before initiation of NICT.

Exclusion criteria were: (1) Incomplete clinical or pathological data; (2) Significant CT artifacts impairing analysis; (3) Prior anti-tumor therapy before NICT; and (4) Concurrent other primary malignancies. Figure 1 illustrates the patient recruitment flowchart.

Figure 1 Patient enrollment process. ESCC: Esophageal squamous cell carcinoma; AJCC: American Joint Committee on Cancer; NICT: Neoadjuvant immunochemotherapy; CT: Computed tomography.

NICT regimen and clinicopathological data collection

Chemotherapy regimens were platinum-based, consisting of a two-drug combination of albumin-bound paclitaxel and a platinum agent, administered at standard doses and schedules in accordance with the National Comprehensive Cancer Network and the Chinese Society of Clinical Oncology guidelines. The immunotherapy involved standard doses of programmed cell death protein 1 (PD-1) monoclonal antibodies (tislelizumab, nivolumab, sintilimab, pembrolizumab, toripalimab, or camrelizumab). All patients were treated with 2-4 cycles of NICT, and each treatment cycle lasted 3 weeks. Tumor response was assessed by contrast-enhanced chest CT after every two cycles. Radical esophagectomy in the McKeown or Ivor Lewis procedure with lymphadenectomy was carried out between 4 weeks and 8 weeks following completion of the final NICT cycle. Details of the NICT regimens for patients in the three cohorts are provided in Supplementary Table 1.

Baseline characteristics, defined as data recorded in patient medical records prior to the initiation of NICT, included: Age, gender, smoking and alcohol consumption history, cT stage, cN stage, TNM stage (according to the AJCC 8^th edition), and Ki-67 expression status in the initial endoscopic biopsy specimen at diagnosis. The total number of NICT cycles completed was also recorded.

Pathological evaluation

The therapeutic response to NICT was assessed using the AJCC tumor regression grade (TRG) system. Treatment response was classified into four grades: TRG 0 indicated no histologically identifiable cancer cells; TRG 1 represented single cells or rare small groups of cancer cells; TRG 2 represented residual cancer with evident tumor regression but more than single cells or rare small groups of cancer cells; and TRG 3 represented extensive residual cancer with no evident tumor regression. Based on postoperative pathological findings, patients attaining TRG 0-1 were assigned to the responder group, while those with TRG 2-3 were designated as non-responders.

Clinical feature selection and model construction

Within the training cohort, clinical characteristics were compared between the groups with different treatment response. Clinical features with statistically significant differences between responders and nonresponders (P < 0.05) were selected for model construction. Using these selected clinical features, four distinct clinical models were developed employing four machine learning algorithms: Support Vector Machine (SVM), RandomForest (RF), K-Nearest Neighbors (KNN), and Multilayer Perceptron (MLP). These models were subsequently applied to predict treatment response in both the internal and external validation cohorts.

CT acquisition protocol

All enrolled patients underwent contrast-enhanced chest CT prior to their initial NICT treatment within 2 weeks. A 64-detector row configuration was employed for the examinations, utilizing either the GE MEDICAL SYSTEMS Optima CT620/670 or the Siemens SOMATOM Definition Flash scanners. The scanning parameters were established as follows: Tube voltage of 120 kVp and slice thickness of 5 mm. Patients were positioned supine, with scan coverage extending from the level of the hyoid bone to 2 cm below the inferior margin of the liver. A non-ionic contrast agent (iohexol, 300 mg/mL, Omnipaque, GE Healthcare) was administered via an antecubital vein using an automatic injector, at an injection rate of 3.5 mL/second and a dose of 1.5 mL/kg. Following contrast injection, arterial phase imaging was conducted at a delay of 40 seconds.

Tumor segmentation and radiomics feature extraction

Primary esophageal tumors were manually delineated layer-by-layer as regions of interest (ROIs) on arterial-phase axial contrast-enhanced chest CT images by two radiologists, each with over five years of experience, using ITK-SNAP software (v.3.8.0, http://www.itksnap.org). Any disagreements during ROI delineation were resolved by consensus, with review by another radiologist with ten years of clinical experience. Tumor boundaries were determined with reference to multiple diagnostic methods, such as endoscopic findings and positron emission tomography/CT. Air and fluid within the esophageal lumen were carefully excluded from the ROIs during delineation.

The CT images and corresponding ROIs were resampled to achieve an isotropic spatial resolution of 1 mm × 1 mm × 1 mm using RIAS software (version 0.2.1). Image feature data were then standardized using the Z-score method. Radiomic features were extracted from the images and ROIs of patients in the training cohort using the PyRadiomics package (version 3.0.1).

Radiomics feature extraction and model development

From ROIs we extracted 1197 radiomic features, which encompassing first-order statistics, shape-based features, texture features, and wavelet-based features. The intra- and inter-class correlation coefficients (ICC) were used to assess the stability and reproducibility of features, and those with an ICC value of 0.8 or higher were retained. An independent samples t test was then applied, identifying 302 features with P < 0.05. Subsequently, the least absolute shrinkage and selection operator (LASSO) regression method was employed to determine the most informative predictive features, with the optimal penalty parameter λ determined via 5-fold cross-validation.

Using the selected radiomic features, four radiomics models were constructed by applying four machine learning algorithms: SVM, RF, KNN, and MLP.

DL model development

CT datasets from different centers were harmonized using a mediastinum window setting. Voxel spacing was resampled to 1 mm × 1 mm × 1 mm using cubic spline interpolation. For each patient, the single axial CT slice showing the maximum cross-sectional area of the primary esophageal tumor was selected as the model input. From this selected slice, the ROI was cropped using the minimum rectangular bounding box that encompassed the tumor, resampled using linear interpolation, and uniformly resized to 448 × 448 pixels. To expand the dataset and enhance model robustness, data augmentation strategies, including random horizontal and vertical flipping, along with random cropping to dimensions of 448 × 448 pixels, were applied. The standardized size for the final input images was set to 448 × 448 pixels. During model training, each image was processed independently. The detailed description of the model architectures used in our study was shown in Supplementary Table 2.

In this study, Densenet201, Resnet101, and Vgg11 networks pre-trained on the ImageNet dataset, which comprises approximately 14 million annotated natural images, were utilized. Optimization of the models was carried out using the Adam optimizer with a batch size of 64 and an initial learning rate of 1 × 10^-2. L2 regularization alongside an early stopping strategy were implemented to prevent overfitting. After training, the DL models were then applied to both the internal and external validation cohorts, generating a probability score of NICT response for each patient.

To interpret the DL models’ predictions and visualize the image regions contributing to the prediction of NICT response, gradient-weighted class activation mapping (Grad-CAM) was applied for network visualization. The resulting heatmap values were normalized to a range of 0-1.

DL feature extraction

After training DL models, we extracted DL features from the images for subsequent construction of the combined models. These features were derived from the final fully connected layer and reduced to 32 dimensions, yielding 32 DL features for each ROI.

CT images from the internal and external validation cohorts underwent the same cropping and standardization procedures applied to the training cohort. The processed images were then input into the trained DL models, and DL features were extracted from the fully connected layer to validate the predictive performance of the models.

Construction of the fusion models

We first constructed the combined radiomics-DL models. A total of 1197 radiomic features were separately combined with the 32 DL features extracted from Densenet201, Resnet101, and Vgg11, generating three radiomics-DL feature sets for feature selection. After Z-score normalization was applied to all radiomics-DL features, an independent samples t test was employed for each feature set to identify features with P < 0.05. LASSO regression was subsequently used to further filter the most informative features. For each of the three filtered radiomics-DL feature sets, four machine learning algorithms—SVM, RF, KNN, and MLP—were applied to construct the models, resulting in a total of 12 radiomics-DL fusion models.

Subsequently, combined clinical-radiomics-DL models were developed. Based on the predictive performance of the different algorithms in the single-modality models and the radiomics-DL models, four clinical features and the 1197 radiomic features were combined with the DL features extracted from the best-performing DL architecture, generated a clinical-radiomics-DL feature set, which was then proceeded to the feature selection stage. After undergoing the same Z-score normalization, independent samples t test, and LASSO regression procedures applied to the radiomics-DL features, the filtered clinical-radiomics-DL feature set was used to construct clinical-radiomics-DL fusion models using the four machine learning algorithms (SVM, RF, KNN, and MLP).

Model evaluation

All feature selection steps (including univariable filtering, correlation screening, and LASSO regression), algorithm selection (SVM, RF, KNN, and MLP), DL architecture selection (Densenet201, Resnet101, and Vgg11), and hyperparameter tuning were performed strictly within the training cohort using 5-fold cross-validation. The internal and external validation cohorts were not accessed at any stage during model development, feature selection, or parameter optimization. The final models were then locked and applied once to the validation cohorts to obtain unbiased performance estimates. The predictive performance of the different models was comprehensively evaluated using receiver operating characteristic (ROC) curve analysis. The evaluation metrics included the area under the curve (AUC), sensitivity, specificity, negative predictive value (NPV), positive predictive value (PPV), and accuracy. For the final selected model with the best performance, the optimal classification threshold was determined using the Youden index within the training cohort. All other candidate models were evaluated at the default 0.5 threshold, as the primary model selection criterion was the AUC, which is threshold-independent. Model calibration, which reflects the goodness-of-fit, was assessed using Brier score. Decision curve analysis (DCA) was used to determine the clinical utility of the models across various threshold probabilities by quantifying the net benefit. Figure 2 shows a schematic diagram of the study.

Figure 2 Study workflow. SVM: Support Vector Machine; KNN: K-Nearest Neighbors; MLP: Multilayer Perceptron; AUC: Area under the curve; TPR: True positive rate; FPR: False positive rate.

Statistical analysis

Statistical analyses were performed using R software (version 4.2.1) and Python (version 3.12.2). For baseline clinical characteristics, continuous variables were expressed as means with standard deviations, while categorical variables were presented as frequencies. Comparisons between groups were conducted using appropriate statistical tests: For continuous variables, either the Student's t test or the Mann-Whitney U test was applied based on the data distribution, while categorical variables were analyzed using the χ² test or Fisher's exact test. Specifically, for 2 × 2 contingency tables, the χ² test optionally with Yates' continuity correction was used, whereas Fisher's exact test was applied in other relevant scenarios. A two-sided approach was adopted for all statistical tests, while P value < 0.05 was deemed to indicate statistical significance.

RESULTS

Patient characteristics

This study enrolled 188 patients with resectable locally advanced ESCC, who were divided into three cohorts: A training cohort (n = 107), an internal validation cohort (n = 47), and an external validation cohort (n = 34). The response rates to NICT were 45.8% for the training cohort, 53.2% for the internal validation cohort, and 58.8% for the external validation cohort. Baseline patient characteristics of each cohort are presented in Table 1. Within the training cohort, significant intergroup differences between responders and nonresponders were noted for smoking history, alcohol consumption history, cT stage, and Ki-67 expression. In the internal validation cohort, there were significant differences between groups for alcohol consumption history and Ki-67 expression, while in the external validation cohort, smoking history was the sole variable showing a significant difference. No notable differences in baseline characteristics were observed across the three cohorts.

Table 1 Baseline characteristics of the patients.

Characteristics	Overall			Training cohort			Test cohort			Validation cohort			P value
	Non-responder	Responder	P value	Non-responder	Responder	P value	Non-responder	Responder	P value	Non-responder	Responder	P value
Number	94	94		58	49		22	25		14	20		0.367
Age	61 (57, 66)	63 (57.25, 69)	0.165	61 (57, 66)	64 (57, 69)	0.294	59 (55, 66.5)	61 (58.5, 68.5)	0.352	64 (59, 68)	66.5 (60.75, 69)	0.958	0.153
Gender			0.266			1.000			1.000			0.126	0.083
Male	89 (94.7)	85 (90.4)		56 (96.6)	47 (95.9)		19 (86.4)	22 (88)		14 (100)	16 (80)
Female	5 (5.3)	9 (9.6)		2 (3.4)	2 (4.1)		3 (13.6)	3 (12)		0 (0)	4 (20)
Smoking			< 0.001b			0.003b			0.058			0.011a	0.623
Former or current	82 (87.2)	56 (59.6)		50 (86.2)	30 (61.2)		18 (81.8)	14 (56)		14 (100)	12 (60)
Never	12 (12.8)	38 (40.4)		8 (13.8)	19 (38.8)		4 (18.2)	11 (44)		0 (0)	8 (40)
Alcohol			< 0.001b			0.012a			0.004b			0.422	0.961
Former or current	82 (87.2)	60 (63.8)		49 (84.5)	31 (63.3)		21 (95.5)	15 (60)		12 (85.7)	14 (70)
Never	12 (12.8)	34 (36.2)		9 (15.5)	18 (36.7)		1 (4.5)	10 (40)		2 (14.3)	6 (30)
cT stage			0.012a			0.007b			0.303			0.781	0.685
T2	1 (1.1)	2 (2.1)		0 (0)	1 (2)		0 (0)	1 (4)		1 (7.1)	0 (0)
T3	63 (67)	79 (84)		37 (63.8)	42 (85.7)		15 (68.2)	20 (80)		11 (78.6)	17 (85)
T4	30 (31.9)	13 (13.8)		21 (36.2)	6 (12.2)		7 (31.8)	4 (16)		2 (14.3)	3 (15)
cN stage			0.285			0.261			0.391			0.126	0.298
N0	31 (33)	27 (28.7)		16 (27.6)	15 (30.6)		9 (40.9)	10 (40)		6 (42.9)	2 (10)
N1	42 (44.7)	34 (36.2)		28 (48.3)	15 (30.6)		10 (45.5)	7 (28)		4 (28.6)	12 (60)
N2	17 (18.1)	27 (28.7)		12 (20.7)	16 (32.7)		3 (13.6)	7 (28)		2 (14.3)	4 (20)
N3	4 (4.3)	6 (6.4)		2 (3.4)	3 (6.1)		0 (0)	1 (4)		2 (14.3)	2 (10)
cTNM stage			0.100			0.057			0.202			0.427	0.865
II	21 (22.3)	27 (28.7)		12 (20.7)	14 (28.6)		3 (13.6)	9 (36)		6 (42.9)	4 (20)
III	41 (43.6)	48 (51.1)		23 (39.7)	26 (53.1)		12 (54.5)	11 (44)		6 (42.9)	11 (55)
IV	32 (34)	19 (20.2)		23 (39.7)	9 (18.4)		7 (31.8)	5 (20)		2 (14.3)	5 (25)
Tumor location			0.807			0.416			0.993			0.797	0.714
Upper	19 (20.2)	17 (18.1)		12 (20.7)	7 (14.3)		5 (22.7)	6 (24)		2 (14.3)	4 (20)
Middle	48 (51.1)	46 (48.9)		30 (51.7)	23 (46.9)		10 (45.5)	11 (44)		8 (57.1)	12 (60)
Lower	27 (28.7)	31 (33)		16 (27.6)	19 (38.8)		7 (31.8)	8 (32)		4 (28.6)	4 (20)
Tumor length (cm)	8.5 (6, 12)	9 (6.5, 13.5)	0.779	9.4(5.5, 12.375)	10 (6.5, 13.5)	0.341	8 (5.125, 12)	7 (5.5, 13)	0.898	8.5 (8, 9.625)	7.75 (7, 10)	0.370	0.287
Tumor diameter (mm)	36 (30, 41)	38 (34, 44)	0.069	36.5 (30, 41)	37 (34, 44)	0.231	36 (30.5, 44.25)	36 (31, 44)	0.522	36 (31, 39.5)	39 (34, 40)	0.211	0.796
Cycles	2.5 (2, 4)	2 (2, 3)	0.756	3 (2, 4)	3 (2, 3)	0.699	3 (2, 3.75)	2 (2, 3)	0.322	2 (2, 3)	2 (2, 3)	1.000	0.111
Ki-67	60 (40, 70)	35 (10, 60)	< 0.001b	60 (50, 70)	40 (30, 70)	0.003b	50 (40, 60)	30 (10, 32.5)	< 0.001b	30 (25, 50)	60 (10, 70)	0.311	0.121

^aP < 0.05.

^bP < 0.01.

Data were showed as median (IQR) or n (%).

Single-modal model construction and evaluation

Based on the baseline characteristic comparison between responders and nonresponders in the training cohort (Table 1), smoking history, alcohol consumption history, cT stage, and Ki-67 expression were selected as clinical features. Clinical prediction models were subsequently constructed using the SVM, RF, KNN, and MLP algorithms. The ROC curves for the clinical models are presented in Supplementary Figure 1, and detailed evaluation metrics regarding their predictive performance are supplied in Supplementary Table 3.

A total of 1197 radiomic features were derived from the ROIs. An independent samples t test identified 302 features with P < 0.05. Subsequently, Pearson correlation analysis and LASSO regression were applied. Ultimately, 24 features were selected for model construction. Radiomics models were developed using the SVM, RF, KNN, and MLP algorithms. The ROC curves for the radiomics models are displayed in Supplementary Figure 2, and detailed evaluation metrics for their predictive performance are provided in Supplementary Table 4.

Three DL models were developed using the Densenet201, Resnet101, and Vgg11 architectures. Within the training cohort, the AUC values amounted to 0.985 for Densenet201, 0.954 for Resnet101, and 0.889 for Vgg11, respectively. In the internal validation cohort, the corresponding AUC values were 0.820, 0.785, and 0.702. For the external validation cohort, the AUC values were 0.664, 0.664, and 0.682, respectively. The ROC curves for the DL models are presented in Supplementary Figure 3, and detailed evaluation metrics for their predictive performance are provided in Supplementary Table 5.

To interpret the DL models, we used Grad-CAM to visualize key regions contributing to the model predictions. Figure 3 presents representative Grad-CAM images for patients with different treatment responses, illustrating the spatial distribution of the decision-making regions identified by each DL model on CT images. In these visualizations, red and yellow regions indicate areas predicted by the model to be responsive to NICT, whereas green and blue regions denote areas predicted to be non-responsive.

Figure 3  Representative gradient-weighted class activation mapping heatmaps of deep learning models for predicting neoadjuvant immunochemotherapy response.

The activation regions generated by the Densenet201 and Vgg11 models were relatively well-defined and largely overlapped with the primary tumor sites. In contrast, the activation maps produced by the Resnet101 model appeared more diffuse and less clearly localized. These visualization results were consistent with the predictive performance of the models as reflected by their ROC curves (Supplementary Figure 3).

Combined model development and performance evaluation

The 1197 radiomic features were combined with 32 DL features extracted by each of the three architectures (Densenet201, Resnet101, and Vgg11), resulting in three sets of combined radiomics-DL features. For each radiomics-DL feature set, an independent samples t test was performed, identifying 307, 304, and 308 features with P < 0.05, respectively. Subsequently, LASSO regression was applied to each feature set, selecting 24, 22, and 20 features, respectively, for model construction.

A total of 12 radiomics-DL combined models were then developed by applying four machine learning algorithms—SVM, RF, KNN, and MLP—to each of the three selected radiomics-DL feature sets. The ROC curves for these 12 radiomics-DL models are presented in Supplementary Figure 4. Overall, the radiomics-DL combined models demonstrated improved predictive performance compared with the previously established single-modality models. Notably, the best-performing model, Vgg11-RF, yielded AUC values of 0.998, 0.855, and 0.779 across the training, internal validation, and external validation cohorts, respectively (Figure 4). Detailed evaluation metrics for the radiomics-DL models are provided in Supplementary Table 6.

Figure 4 Receiver operating characteristic curves of the best performing models across different feature types in training, internal validation, and external validation cohort. A: Training; B: Internal validation; C: External validation cohort. RD: Radiomics-deep learning; CRD: Clinical-radiomics-deep learning; AUC: Area under the curve; TPR: True positive rate; FPR: False positive rate.

Given that, among the radiomics-DL combined models, the model incorporating DL features extracted by Vgg11 demonstrated the best predictive performance, these features were subsequently used to construct the clinical-radiomics-DL combined models. The four clinical features, 1197 radiomic features, and 32 DL features extracted by Vgg11 were integrated to form the clinical-radiomics-DL feature set. An independent samples t test was performed on this feature set, identifying 311 features with P < 0.05. Subsequently, LASSO regression was applied, ultimately selecting 26 features for model construction. Using these filtered clinical-radiomics-DL features, clinical-radiomics-DL models were developed with four machine learning algorithms: SVM, RF, KNN, and MLP. Supplementary Figure 5 shows the ROC curves corresponding to these clinical-radiomics-DL models.

In the training cohort, the clinical-radiomics-DL models performed at a level similar to the radiomics-DL models, and in the internal validation cohort, the clinical-radiomics-DL models showed enhanced predictive performance compared with the radiomics-DL models. However, in the external validation cohort, the clinical-radiomics-DL models underperformed compared with the radiomics-DL models. The AUC values for the SVM, RF, KNN, and MLP models were 0.643, 0.600, 0.736, and 0.657, respectively (Figure 4). Detailed evaluation metrics for the clinical-radiomics-DL models are provided in Supplementary Table 7.

Model comparison

From the models constructed using each type of feature set, the performance was comprehensively compared primarily based on the AUC, and the model with the best overall performance across all three cohorts was chosen as the representative model. These included the clinical model using RF architecture, the radiomics model using SVM architecture, the DL model using Densenet201 achitecture, the radiomics-DL model (Vgg11-RF), and the clinical-radiomics-DL model (Vgg11-SVM) model. A comprehensive comparison was subsequently conducted in terms of predictive performance, goodness-of-fit, and clinical utility.

Figure 4 presents the ROC curves of the best-performing models developed using different feature types. In the training cohort, all models exhibited excellent predictive performance. In the internal validation cohort, all models maintained good predictive capability. However, there was a marked decline in predictive performance across all models in the external validation cohort, except the radiomics-DL model, whose performance remained comparatively stable. Table 2 summarizes predictive performance of best-performing models across different feature types. The radiomics-DL model yielded AUCs of 0.998, 0.855, and 0.779 across the three cohorts. After 500-iteration bootstrap optimism correction, the training-set AUC was 0.921 (95%CI: 0.853-0.987). Moreover, with the Youden-derived threshold of 0.4173, the radiomics-DL model achieved a sensitivity of 0.900, specificity of 0.500, PPV of 0.720, NPV of 0.778, and accuracy of 0.735 in the external validation cohort, indicating a more balanced trade-off between sensitivity and specificity compared with the default 0.5 threshold (Table 3). Figure 5 is the confusion matrices of the radiomics-DL model in training, internal validation, and external validation cohort. Supplementary Figure 6 shows the ROC curves of the radiomics-DL model for predicting NICT response across different treatment regimen subgroups (stratified by platinum agent type, period from treatment to surgery, and number of treatment cycles) within the training and internal validation cohorts.

Figure 5 Confusion matrices of the radiomics-deep learning model in the training, internal validation, and external validation cohorts at the Youden-derived threshold of 0.4173. A: Training; B: Internal validation; C: External validation. Values are presented as counts (percentages of the cohort).

Table 2 Predictive performances of the best performing models developing on different kind of features.

Model	Cohort	AUC	95%CI	Sensitivity	Specificity	NPV	PPV	Accuracy
Clinical model	Train	0.998	0.995-1.000	0.927	1.000	0.957	1.000	0.972
	Internal validation	0.899	0.813-0.985	0.720	0.909	0.741	0.900	0.809
	External validation	0.636	0.441-0.830	0.700	0.571	0.700	0.571	0.647
Radiomics model	Train	0.996	0.989-1.000	0.951	0.985	0.970	0.975	0.972
	Internal validation	0.904	0.816-0.991	0.600	0.955	0.677	0.938	0.766
	External validation	0.439	0.228-0.650	0.450	0.429	0.529	0.353	0.441
Deep learning model	Train	0.985	0.968-1.000	0.951	0.924	0.968	0.886	0.935
	Internal validation	0.820	0.699-0.941	0.760	0.864	0.760	0.864	0.809
	External validation	0.664	0.470-0.859	0.900	0.429	0.750	0.692	0.706
Radiomics-deep learning model	Train	0.998	0.994-1.000	0.927	1.000	0.957	1.000	0.972
	Internal validation	0.855	0.742-0.967	0.600	0.864	0.655	0.833	0.723
	External validation	0.779	0.609-0.948	1.000	0.357	1.000	0.690	0.735
Clinical-radiomics-deep learning model	Train	1.000	1.000-1.000	1.000	1.000	1.000	1.000	1.000
	Internal validation	0.887	0.794-0.980	0.680	0.909	0.714	0.895	0.787
	External validation	0.643	0.450-0.836	0.700	0.286	0.400	0.583	0.529

AUC: Area under the curve; NPV: Negative predictive value; PPV: Positive predictive value.

Table 3 Predictive performances of the radiomics-deep learning model With the Youden-derived threshold of 0.4173.

Radiomics-deep learning model (Youden-derived threshold)	Sensitivity	Specificity	NPV	PPV	Accuracy
Train	1.000	0.970	1.000	0.953	0.981
Internal validation	0.760	0.682	0.714	0.731	0.723
External validation	0.900	0.500	0.778	0.720	0.735

NPV: Negative predictive value; PPV: Positive predictive value.

DCA was additionally employed to further evaluate the clinical utility of the different models (Figure 6). In both the training and internal validation cohorts, all models provided a certain degree of clinical net benefit. However, in the external validation cohort, only the radiomics-DL model consistently demonstrated stable clinical benefit across a broad spectrum of threshold probabilities. The DCA of the radiomics-DL model in the three cohorts with confidence ranges are provided in Supplementary Figure 7.

Figure 6 Decision curve analysis of the best-performing models across different feature types in the training, internal validation, and external validation cohorts. A: Training; B: Internal validation; C: External validation. RD: Radiomics-deep learning; CRD: Clinical-radiomics-deep learning.

Table 4 presents the Brier scores of the different models. All models demonstrated good calibration in the training cohort, as indicated by their low Brier scores. In the internal validation cohort, there was a decline in calibration performance across all models. However, the calibration performance of most models further deteriorated in the external validation cohort, whereas the radiomics-DL model maintained a calibration level comparable to that observed in the internal validation cohort. Supplementary Figure 8 shows the calibration curves of the radiomics-DL model in the three cohorts.

Table 4 Brier score of the best performing models developing on different kind of features.

	Train	Internal validation	External validation
Clinical model	0.026	0.131	0.279
Radiomics model	0.030	0.164	0.323
Deep learning model	0.058	0.179	0.294
Radiomics-deep learning models	0.052	0.182	0.189
Clinical-radiomics-deep learning models	< 0.001	0.168	0.296

DISCUSSION

Multiple clinical trials and meta-analyses have demonstrated that, for locally advanced ESCC patients, NICT achieves therapeutic efficacy comparable with neoadjuvant chemoradiotherapy (NCRT), which is currently considered the standard of care[4,21]. However, treatment response to NICT varies substantially among individual patients, and no reliable method is currently available to accurately predict this response. Therefore, the development of an effective predictive model for NICT response is of considerable clinical importance. In this research, we constructed clinical models, conventional radiomics models, DL models, radiomics-DL models, and clinical-radiomics-DL models for the purpose of predicting treatment response in patients with resectable locally advanced ESCC. Among these, the radiomics-DL model (Vgg11-RF) demonstrated the best overall predictive performance, yielding AUC values of 0.998, 0.855, and 0.779 for the training, internal validation, and external validation cohorts, respectively. Furthermore, the Brier score indicated that the radiomics-DL model exhibited good calibration, and DCA also indicated favorable clinical utility.

Conventional radiomics techniques enable us to extract numerous quantitative features from medical images that remain imperceptible to the human eye. By leveraging these features to construct models, radiomics can help clinicians in predicting treatment outcomes. For NCRT, which remains the standard preoperative treatment for esophageal cancer, multiple studies have developed radiomics models to predict therapeutic response, demonstrating promising performance[22-28]. In recent years, the clinical use of NICT for resectable locally advanced ESCC has increased substantially. Correspondingly, radiomics studies aimed at predicting response to NICT have also emerged[16-20,29-33]. However, there are limitations in most of these studies, such as small sample sizes, single-center designs and reliance on single-modality data, which may restrict the generalizability and clinical applicability of the proposed models.

Validation of model performance in multicenter cohorts is essential to ensure their clinical applicability across different institutions and patient populations[18]. In the present study, the single-modality models and the clinical-radiomics-DL model demonstrated strong performance in both the training and internal validation cohorts; however, their predictive performance declined substantially in the external validation cohort. This finding highlights the significance of multicenter assessment in validating model robustness and generalizability. Currently, many studies developing predictive models for NICT response in ESCC are conducted within single-center cohorts[16,29-32]. Although these models may exhibit strong performance in the populations from which they are derived, their generalizability to independent cohorts is often not evaluated, thereby limiting their reliability and clinical applicability. In contrast, the radiomics-DL model in this research maintained its predictive performance in the external validation cohort, suggesting improved robustness and generalizability across different centers. Nevertheless, owing to the limited sample size of the external validation cohort (n = 34), these results should be viewed with appropriate caution, and further validation using data from additional institutions is warranted to facilitate its translation into routine clinical practice.

Unlike conventional radiomics, DL approaches leverage the nonlinear representational capacity of neural networks to automatically extract complex imaging features and construct predictive models. However, DL methods may be more sensitive to variations in image acquisition and preprocessing, whereas handcrafted radiomics features tend to exhibit greater robustness. In addition, when training data are limited, DL models are more susceptible to overfitting. Consequently, multimodal approaches that integrate complementary features from different sources have been increasingly explored and have demonstrated superior performance compared to single-modality models[12,13,34,35]. Notably, in our DL pipeline, we used only the single axial slice with the maximal crosssectional tumor area as the model input. While this single-slice strategy is consistent with prior radiomics and DL studies, it may not fully represent wholetumor heterogeneity. Volumetric analysis using multiple slices, which may further enhance predictive performance, is a valuable direction for future investigation[19,36].

In the context of predicting NICT response in patients with resectable locally advanced ESCC, several multicenter studies have adopted such multimodal strategies and reported improved predictive performance. A fusion model incorporating habitat-based radiomics features with multi-channel 2D DL features was developed by Zhu et al[19], achieving AUC values of 0.910, 0.829, and 0.793 across the training, internal validation, and external validation cohorts, respectively. Similarly, Zhang et al[18] constructed a fusion model using a Vision-Mamba-based DL framework combining voxel-level radiomics features with CT images, achieving an AUC of 0.91 in the training cohort and AUCs ranging from 0.83 to 0.92 across multiple validation cohorts. In both studies, the combined models consistently outperformed single-modality models, particularly in external validation cohorts where predictive performance remained relatively stable.

Consistent with these findings, the radiomics-DL model constructed in our investigation exhibited enhanced predictive performance relative to all single-modality models and maintained stable performance in the external validation cohort. This improvement could be attributed to complementary nature of radiomics and DL features, which capture both handcrafted and high-level abstract imaging information. Notably, although the radiomics-DL model demonstrated favorable predictive performance across different treatment subgroups in both the training and internal validation cohorts, effective subgroup analyses could not be performed in the external validation cohort due to its limited sample size. In addition, given the wide variety of PD-1 monoclonal antibodies used in this study and the consequently small patient numbers within each immunotherapy subgroup, meaningful subgroup analyses stratified by specific immunotherapy agents were also not feasible. Future studies should address this limitation by performing subgroup analyses according to immunotherapy type, or by conducting prospective studies with standardized immunotherapy regimens.

In some radiomics-based predictive modeling studies, integrating clinical features with imaging features has been shown to further improve predictive performance[27,30,31,37-39]. However, many of these studies lack external validation of the combined models. Qi et al[40] investigated a fusion strategy combining conventional radiomics and pathomics to develop a multimodal model for predicting NICT response in patients with ESCC. Their results suggested that directly constructing a model by integrating heterogeneous features yielded better performance than approaches that fuse outputs from independently developed single-modality models. Accordingly, we adopted this direct integration strategy in the present study. Despite this, although in training and internal validation cohorts the clinical-radiomics-DL models demonstrated favorable performance, in the external validation cohort their predictive performance, calibration, and clinical utility declined substantially. In contrast, the radiomics-DL models maintained stable performance across all cohorts. This discrepancy may be attributable to the relatively lower consistency of clinical features across centers, which could adversely affect model generalizability after feature integration. While no significant differences in baseline characteristics were detected among the training, internal validation, and external validation cohorts, unmeasured heterogeneity likely remained. Such heterogeneity may include differences in socioeconomic conditions across regions and variations in diagnostic and treatment practices between institutions, which are inherently difficult to quantify and control. Therefore, future radiomics studies incorporating clinical features into multimodal models should emphasize rigorous multicenter validation to ensure robustness and generalizability.

It is worth noting that delta radiomics has shown substantial progress in predicting treatment response in oncology. This approach characterizes changes in radiomic features between pre- and post-treatment images and incorporates these variations as novel features into predictive models[27]. Compared with conventional radiomics, delta radiomics places greater emphasis on capturing treatment-induced alterations in the tumor microenvironment, thereby potentially improving the assessment of therapeutic response[33,41]. However, a key limitation of delta radiomics is that it relies on post-treatment imaging, thereby restricting its utility for early clinical decision-making. Consequently, patients with poor response to NICT may experience delays in surgical intervention, potentially increasing the risk of disease progression. In contrast, the strategy adopted in this research—constructing predictive models solely using pre-treatment imaging—may sacrifice some temporal information but enables earlier prediction of treatment response. This characteristic may offer greater clinical value by facilitating timely decision-making and supporting the implementation of precision medicine.

Beyond predictive performance, the clinical utility of our model lies in its potential to inform multidisciplinary treatment decisions before the initiation of NICT. For patients predicted to respond, the model provides reassurance to proceed with immunotherapy-based regimens; for those predicted to be non-responsive, clinicians may consider alternative approaches such as NCRT or upfront surgery, thereby avoiding unnecessary toxicity and potential disease progression. Nevertheless, for the model to be seamlessly adopted into complex and busy clinical workflows, simplicity and user-friendliness are paramount. In this regard, recent work in the non-small cell lung cancer field has provided an encouraging precedent: Researchers have developed and openly published a software for predicting neoadjuvant chemoimmunotherapy response[42], setting an important example for the translational application of radiomics-based predictive models in clinical practice. Prospective studies and user-friendly implementation strategies for the prediction of NICT response of ESCC are also required before routine clinical adoption. In addition, the calibration curve suggested a tendency of the radiomics-DL model to underestimate response probabilities in the internal validation cohort. This is consistent with the observation that the Youdenderived optimal threshold (0.4173) was lower than the default 0.5. Nevertheless, we acknowledge that threshold adjustment does not correct the underlying probability estimates; future work may consider Platt scaling or isotonic calibration to improve the absolute accuracy of the predicted probabilities.

Several limitations of this study should be acknowledged. First, despite its multi-center design, the sample size remained relatively modest, and model performance was evaluated using only a single external validation cohort. While our results demonstrate the potential utility of the radiomics-DL model, widespread clinical adoption should be preceded by rigorous prospective validation in diverse patient populations and across multiple institutions. In particular, the limited size of the external validation cohort (n = 34) warrants caution when interpreting the model's generalizability. The model should currently be viewed as a research tool to support clinical decision-making, rather than a standalone diagnostic or therapeutic guide. Second, the present study was retrospective in nature. Future prospective work is warranted to validate the clinical value of the proposed model for guiding treatment decisions. Third, due to incomplete follow-up data, postoperative survival outcomes were unavailable for a subset of patients; consequently, the model’s ability to stratify patients according to survival outcomes could not be evaluated. Finally, radiomic features and the models derived from them have limited biological interpretability, which may hinder their clinical adoption. Future studies integrating multi-omics approaches are needed to provide deeper biological insights into radiomic features and to enhance the interpretability and translational potential of these models.

CONCLUSION

In this study, multiple models were developed to predict NICT response in resectable locally advanced ESCC patients based on pretreatment contrast-enhanced chest CT images and clinical data. Among them, the radiomics-DL model (RF-Vgg11) demonstrated the best performance, achieving AUC values of 0.998, 0.855, and 0.779 for the training, internal validation, and external validation cohorts, respectively. This model showed promising predictive accuracy and generalizability, suggesting its potential as a noninvasive decision-support tool. However, given the limited external validation sample size in the external cohort, these findings should be interpreted with caution. Prospective investigations with larger and more heterogeneous patient populations are essential to substantiate the clinical utility and robustness of this model before it can be integrated into routine practice.

ACKNOWLEDGEMENTS

We would like to thank everyone who has helped with our research.