This study was a single-center retrospective study registered in the Chinese Clinical Trial Registry (approval no. ChiCTR2000041147). Our study was approved by the Ethics Committee of Shanghai Jiao Tong University Affiliated Sixth People’s Hospital (approval no. 2020-253). The cohort from Shanghai Jiao Tong University Affiliated Sixth People’s Hospital was used to develop and validate the model. All patients were pathologically diagnosed with stage II–III CRC between January 2009 and December 2016. Written informed consent was obtained for this study. Patients fulfilling the criteria patients were excluded: 1) age <18 years; 2) had emergency surgery; 3) had multiple primary carcinoma; 4) with incomplete clinical data; 5) died within 30 days; 6) lost to follow-up; and 7) underwent preoperative adjuvant therapy.

In our study, clinical features such as gender, age, pTNM stage, tumor location, tumor cross-sectional area (CSA), tumor long axis, tumor differentiation, lymphatic infiltration, vascular infiltration, nerve infiltration, neutrophil-to-lymphocyte ratio (NLR), platelet-to-lymphocyte ratio (PLR), preoperative tumor markers [carcinoembryonic antigen (CEA), carbohydrate antigen 19-9 (CA19-9), and CA125], and the MMR status were collected. Patient clinical data were mainly provided by examination of their medical history and by the Electronic Medical Record Department. Pathological staging was based on the 8th AJCC criterion for CRC. NLR and PLR were calculated as (neutrophil count)/(lymphocyte count) and (platelet count)/(lymphocyte count), respectively. Preoperative tumor markers were examined within 1 week before surgery. All patients were followed up according to the current National Comprehensive Cancer Network (NCCN) guidelines, including analysis of serum tumor markers, colonoscopy, chest X-ray, and CT (or MRI). Patient follow-up data were updated by telephone, email, and medical history. Disease-free survival (DFS) was defined as the time from surgery to cancer metastasis or recurrence. Overall survival (OS) was defined as the time from surgery to death.

Immunostaining of CD3⁺ and CD8⁺ T cells was performed on formalin-fixed paraffin-embedded sections. Antigen retrieval was conducted with an EDTA buffer (pH 9.0) for 90 s, followed by quenching of endogenous peroxidase activity by 3% H₂O₂ for 30 min at room temperature. Sections were incubated at 4°C with primary antibodies: rabbit anti-human monoclonal antibody against CD3 (EP41; ZSGB-BIO, Beijing, China) and rabbit anti-human monoclonal antibody against CD8 (SP16; ZSGB-BIO, Beijing, China). Revelation with the Ultra DAB IHC Detection Kit (Maxim, Fuzhou, China) and counterstaining with Harris hematoxylin were performed. Counterstained slides were scanned at ×40 magnification (NanoZoomer S360, Hamamatsu, Japan) to generate a whole slide imaging file in NDPI format. CT was the core of the tumor, and the invasive margin (IM) was defined as a region of 500-μm width surrounding the CT. The CT and IM regions were manually marked on the whole slide using QuPath software (20), in which hematoxylin/eosin-stained sections were used to help CT/IM labeling. Two independent pathologists, who were blinded to the patients’ clinical information, participated in the analysis to avoid the interference of necrotic areas and to verify the location of the CT/IM. Positive CD3 and CD8 cells within the CT and IM areas were obtained via QuPath software (20), and the densities of CD3 and CD8 were quantified by the number of cells per square millimeter in both CT and IM. The concordance in the semi-quantitative evaluation between CD3 and CD8 was determined by two independent pathologists.

For every patient, the densities of CD3⁺ and CD8⁺ cells in the CT and IM regions (CD3_CT, CD3_IM, CD8_CT, and CD8_IM) were converted into percentiles (0%–100%) based on our cohort, as described by Galon et al. (9). The mean of the four percentiles (CD3_CT, CD3_IM, CD8_CT, and CD8_IM) was then calculated and converted into a percentile Immunoscore. In a three-category Immunoscore analysis, a 0%–25% density was scored as low, 25%–70% density was scored as intermediate, and 70%–100% density was scored as high (9). In our study, we found that a low (0%–25%) and an intermediate (25%–70%) Immunoscore in the three-category Immunoscore had similar clinical outcomes (DFS). Consequently, we combined the low-Immunoscore (0%–25%) and intermediate-Immunoscore (25%–70%) groups as the low-Immunoscore group (0%–70%) in a two-category Immunoscore, in which a 0%–70% density was scored as low and 70%–100% was scored as high. Samples were excluded from the analysis if counts were missing from a tumor region, if there was improper histology (e.g., broken tissue, atypical CT/IM, excessive necrotic cavity, excessive mucous area, etc.), or if the staining intensity was regarded as low.

The tumor DNA MMR status was determined by immunohistochemical analysis of MMR proteins (MLH1, MSH2, MSH6, and PMS2) on formalin-fixed paraffin-embedded sections. The conditions of having deficient MMR (dMMR; loss of at least one MMR protein) and proficient MMR (pMMR) were denoted as microsatellite instability (MSI) and microsatellite stability (MSS), respectively. Two independent pathologists, who were blinded to the patients’ clinical information, participated in the analysis to verify the MMR status.

The association between the clinicopathological characteristics and the Immunoscore was analyzed via a chi-squared test. All numeric variables were tested for normality using the Shapiro–Wilk test.

To develop the prognostic model, we performed univariate analysis of all variables using Cox proportional hazards regression. Subsequently, the significant variables (p < 0.05) were analyzed with the multivariate Cox proportional hazards regression. After removing the non-significant covariates in the multivariate analysis, a final multivariable Cox regression model was constructed. A nomogram was constructed to predict the 3- and 5-year DFS/OS probabilities with the total points of all variables. The risk score was the linear predictor of the Cox model built on our cohort with selected variables. To evaluate the predictive accuracy of the different variables or models, we used the integrated area under the receiver operating characteristic (ROC) curve (iAUC) with 1,000× bootstrap resampling. The performances of the models were compared using likelihood ratio tests, when the models were nested. The relative importance of each variable to the risks of recurrence and death was estimated using the χ ² from Harrell’s rms R package (version 6.0-1).

Model performance was evaluated with the concordance index (C-index) and corrected 1,000 times by bootstrapping. The calibration curves of the nomogram were drawn for 3- and 5-year DFS/OS to evaluate the accuracy of the model by comparing the DFS/OS probabilities between observations and predictions. A time-dependent ROC was used to compare the discrimination between our nomogram and pTNM. Patients were classified into three risk groups (high, intermediate, and low) according to the risk scores obtained from the nomogram: the 30% with the highest scores were designated the “high” risk group, the 30% with the lowest scores the “low” group, and the remaining 40% as the “intermediate” group. The Kaplan–Meier (K-M) method was applied to estimate the survival probabilities. Hazard ratios (HRs) and 95% confidence intervals (CIs) were estimated, and a log-rank test was used to determine the statistical differences between different groups. Decision curve analysis was conducted using the ggDCA package in R (version 1.2) to determine the clinical usefulness of the nomogram via quantifying the net benefits at different threshold probabilities (21).

Data processing, data analysis, and figures were performed and produced in R language (version 4.0.3). All analyses were two sided, and p < 0.05 was considered statistically significant.

A total of 1,048 patients with stage II–III CRC were collected between January 2009 and December 2016 from Shanghai Jiao Tong University Affiliated Sixth People’s Hospital ( Supplementary Figure S1 ). After clinical quality control, there were 735 eligible patients with complete clinical data. Formalin-fixed paraffin-embedded sections of tumor samples from 350 out of the 735 patients were collected and the Immunoscore data were retrieved between 2020 and 2021. Within the 350 samples, 96 were excluded due to mismatch in the quality control, among which 67 patients were excluded after histology quality control, 20 patients were excluded after staining quality control, and 9 patients were excluded due to missing staining data. Subsequently, only eligible patients with qualified immunohistochemical data samples (n = 254) were finally included in the development and validation of the prognostic models.

The characteristics of our study population are shown in Table 1 . In total, 60.0% of patients were males, and the median age of all patients was 66.0 years (IQR = 56–76 years). One hundred fifty-one (59.0%) patients had stage II and 103 (41.0%) had stage III CRC. Colon tumors located on the left and right sides were 73 (29.0%) and 84 (33.0%), respectively, and 97 (38.0%) patients had rectum tumors. The degree of differentiation in more than half of the tumors was identified as moderate or well (147, 58.0%), and 107 (42.0%) patients had a poor level. Of the patients, 235 (93.0%) showed MSS and only 19 (7.0%) showed MSI. Seventy-seven (30.0%) patients had a relapse, and 74 (29.0%) patients died. The median follow-up time for all patients was 53.0 months.

The median densities of CD3⁺ T cells in CT and IM were 227/mm² (14–2,061/mm²) and 629/mm² (38–2,359/mm²), respectively, and those of CD8⁺ T cells were 119/mm² (8–1,832/mm²) and 340/mm² (13–1,365/mm²), respectively. More than half of CRC patients had a low Immunoscore (189, 74.0%), while 26.0% (65) of patients had a high Immunoscore ( Table 1 ).

To study the association between the Immunoscore and other characteristics in the tumor microenvironment, we performed chi-squared test analysis ( Supplementary Table S1 ). We did not find any relationship between the Immunoscore and other characteristics, except for the microsatellite status. The results showed that a high Immunoscore was found more frequently than a low Immunoscore in tumors with dMMR (14.0% vs. 5.0%, p = 0.0468). Although tumor location was not significantly associated with the Immunoscore (p = 0.0689), it showed some trends, and a high Immunoscore was found to be less frequent in tumors in right-sided colon cancer ( Supplementary Table S1 ).

Representative images of CD3⁺ and CD8⁺ T-lymphocyte immunostaining on formalin-fixed paraffin-embedded sections are provided in Figure 1 . The CT and IM areas of the tumor were manually marked on the whole slide using QuPath software (20), in which hematoxylin/eosin-stained sections were used to help in CT/IM labeling. We validated the two-level categorical Immunoscore ( Supplementary Figure S2 ) whose prognostic impact was previously shown in an international validation study in TNM stage I–III colon cancers (9). When tumors were categorized into predetermined low (0%–70%) and high (70%–100%) groups, a low Immunoscore was associated with a statistically significant and poorer DFS (p = 0.0390) and OS (p = 0.0070). The 3-year DFS for low vs. high Immunoscore was 68.7% vs. 82.5%, and the 3-year OS was 75.9% vs. 87.4%. The 5-year DFS for low vs. high Immunoscore was 62.6% vs. 77.3%, and the 5-year OS was 63.1% vs. 82.1% ( Supplementary Figure S2 ). For consistency with prior work, the associations between DFS and OS were also shown for percentile Immunoscore and three-level categorical Immunoscore (9) ( Supplementary Table S2 ).

To screen the prognostic factors, we performed univariate analysis. Our results suggested that the pT stage, pN stage, lymphatic infiltration, CEA, and the Immunoscore were significantly associated with DFS and OS in CRC patients (p < 0.05) ( Supplementary Table S2 ). Vascular infiltration affected DFS, but not OS, whereas age, NLR, and CA125 affected OS, but not DFS (p < 0.05) ( Supplementary Table S2 ).

Compared to patients with T1, T2, and T3 tumors, patients with T4 tumors had significantly worse DFS (HR = 2.05, 95% CI = 1.11–3.80, p = 0.0226) and OS (HR = 2.58, 95% CI = 1.32–5.04, p = 0.0054) ( Supplementary Table S2 ). Patients with N2 had statistically significant and worse DFS (HR = 3.98, 95% CI = 2.34–6.78, p < 0.0001) and OS (HR = 4.27, 95% CI = 2.47–7.37, p < 0.0001) than patients with N0. Patients with N1 also had worse OS (HR = 2.05, 95% CI = 1.11–3.80, p = 0.0226), but not DFS, than patients with N0 ( Supplementary Table S2 ). Patients with lymphatic infiltration displayed a statistically significant association with DFS (p < 0.05) and OS (p < 0.05) ( Supplementary Table S2 ) and were therefore included in the multivariable Cox models. An increased level of CEA was associated with shorter DFS (p < 0.0001) and OS (p < 0.01), which achieved statistical significance ( Supplementary Table S2 ). In the univariate analysis, the Immunoscore was analyzed as a two-level categorical variable, and a high Immunoscore was associated with a statistically significant and better DFS (HR = 0.54, 95% CI = 0.30–0.98, p = 0.0421) and OS (HR = 0.41, 95% CI = 0.21–0.80, p = 0.0092) ( Supplementary Table S2 ). For consistency with prior studies, the associations with DFS and OS were also shown for the percentile Immunoscore and three-level categorical Immunoscore (9), and both were protective prognostic factors for DFS (p < 0.05) and OS (p < 0.05) ( Supplementary Table S2 ). The MMR status was not prognostic either in DFS (p = 0.5753) or OS (p = 0.3085) ( Supplementary Table S2 ). Subsequently, all variables statistically significant in the univariate analysis were entered into the multivariate Cox proportional hazards regression analysis.

Multivariate Cox proportional hazards regression analysis for DFS and OS revealed that the pN stage, vascular infiltration, CEA, and the Immunoscore were independently associated with DFS (p < 0.05) ( Table 2 ), and age, pN stage, CEA, CA125, and the Immunoscore were independently associated with OS (p < 0.05) ( Table 2 ) in our cohort. The predictive accuracy of the Immunoscore was evaluated by determining the time-dependent AUC ( Supplementary Figure S3 ). For DFS, the predictive accuracy of the Immunoscore was found to be similar to that of the PLR (p > 0.05) and was superior to that of gender, age, tumor location, tumor CSA, tumor long axis, tumor differentiation, nerve infiltration, NLR, CA19-9, CA125, or MMR (p < 0.05), however was lower than that of pT stage, pN stage, lymphatic infiltration, vascular infiltration, and CEA (p < 0.05) ( Supplementary Figure S3A ). For OS, the predictive accuracy of the Immunoscore was found to be similar to that of pT stage (p > 0.05) and was superior to that of gender, age, tumor location, tumor CSA, tumor long axis, tumor differentiation, lymphatic infiltration, vascular infiltration, nerve infiltration, NLR, PLR, CA19-9, or MMR (p < 0.05), however was lower than that of pN stage, CEA, and CA125 (p < 0.05) ( Supplementary Figure S3B ). Furthermore, adding preoperative serum tumor markers (CEA, CA19-9, and CA125) or the Immunoscore to a model that combined all clinical variables (gender, age, pT stage, pN stage, tumor location, tumor CSA, tumor long axis, tumor differentiation, lymphatic infiltration, vascular infiltration, and nerve infiltration) significantly improved both DFS (likelihood ratio: p = 0.0052 and p = 0.0276, respectively) and OS (likelihood ratio: p = 0.0004 and p = 0.0117, respectively) prediction ( Supplementary Figures S3A, B ). Therefore, clinical variables, preoperative serum tumor markers, and the Immunoscore were all required to optimize the determination of patient prognosis.

Variables that were statistically significant in the multivariate Cox analysis were used to develop the final prognostic models, which included independent variables that were associated with DFS and OS. Finally, four indicators were selected for the prognostic model of DFS in CRC, including pN stage, vascular infiltration, CEA, and the Immunoscore ( Table 2 ), and five indicators were selected for OS prediction in CRC, including age, pN stage, CEA, CA125, and the Immunoscore ( Table 2 ).

The final model was then used to generate a nomogram to predict the DFS ( Figure 2A ) and OS ( Figure 2B ) rates for individual patients in clinical practice. The nomogram assigns points to each variable and allows predicting the DFS/OS probabilities at 3 and 5 years using the total points of all variables. The risk scores were generated as a linear predictor of the Cox model and were calculated as follows: 1) risk score for DFS in CRC patients: (0*N0 + 0.1488*N1 + 1.2511*N2) + 0.7320*(vascular infiltration) + 0.3625*ln(CEA value) + −0.7759*(Immunoscore); 2) risk score for OS in CRC patients: 0.0336*age + (0*N0 + 0.5231*N1 + 1.7046*N2) + 0.2192*ln(CEA value) + 0.4896*ln(CA125 value) + −0.8166*(Immunoscore).

We analyzed the relative importance of all variables in the final multivariable model, and the results revealed that the pN stage (41.64%) had the largest impact on DFS, followed by CEA (34.26%), vascular infiltration (12.18%), and then the Immunoscore (11.92%) ( Figure 3A ). For OS, the pN stage (50.08%) had the largest impact, followed by age (16.94%), CA125 (16.12%), CEA (8.28%), and then the Immunoscore (8.59%) ( Figure 3B ).

The C-index of the nomogram was 0.6941, corrected with 1,000 permutations, for DFS in CRC in our cohort ( Table 2 ). The C-index of our OS model was 0.7138, corrected with 1,000 permutations ( Table 2 ). Notably, the C-index of pTNM based on the 8th edition of AJCC was 0.6456 for DFS and was 0.6647 for OS. The calibration curves for CRC based on the nomograms showed very good agreement between the predicted and observed probabilities of DFS and OS at 3 and 5 years ( Figures 4A–D ). Consistently, our nomogram also showed a slightly higher prognostic accuracy than the pTNM stage from 30 to 70 months for both DFS (3-year AUC: nomogram = 0.74, pTNM = 0.69; 5-year AUC: nomogram = 0.75, pTNM = 0.71) ( Figure 4E ) and OS (3-year AUC: nomogram = 0.75, pTNM = 0.68; 5-year AUC: nomogram = 0.78, pTNM = 0.74) ( Figure 4F ) in the time-dependent ROC analysis.

According to the risk scores obtained from the nomogram, the patients were categorized into three risk groups: 30% of patients with the highest scores classified as “high”, 30% with the lowest scores as “low”, and the rest (40%) as “intermediate”. Consequently, the cutoff values were 0.269/0.887 for DFS and 3.435/4.344 for OS in CRC. K-M curves were applied to compare the survival differences. The K-M curve analysis showed statistically significant differences among the different risk groups (p < 0.0001) ( Figures 5A, B ). In the subgroup analysis, it was found that the correlation between nomogram-based risk stratification and DFS/OS was significant both in the subsets of stage II (DFS: p = 0.0017; OS: p = 0.0014) ( Figures 5C, D ) and stage III (DFS: p = 0.0028; OS: p = 0.0031) ( Figures 5E, F ) patients. Therefore, the risk scores generated based on our prognostic models efficiently distinguished the prognosis of patients with stage II or III CRC.

Based on our prognostic models, some patients previously believed to have a high risk of relapse or death rate were found to be at low risk. The classification of patients with stage III CRC into the low-risk (T_1–3N₁) and high-risk (T₄ or N₂) groups is routinely used to guide the treatment of adjuvant FOLFOX (folinic acid–fluorouracil–oxaliplatin) or CAPOX (capecitabine–oxaliplatin) in clinical practice (22). Based on this classification, our cohort included 17 low-risk (T_1–3N₁) and 86 high-risk (T₄ or N₂) patients with stage III CRC. As the number of patients in the low-risk group was too low, we only analyzed the high-risk group using our prognostic models. The analysis showed that our nomogram could significantly identify a group of patients with good OS, but not DFS, within clinically high-risk stage III CRC ( Supplementary Figures S4A, B ). Furthermore, we performed a similar analysis in patients with stage II CRC. High-risk stage II CRC had the characteristics of positive biomarkers for vascular infiltration, lymphatic infiltration, or nerve infiltration (VILINI⁺) or T₄ stage II, whereas low-risk CRC was negative for VILINI markers (VILINI^–) and T_1–3 stage II (22). Our cohort included two low-risk (VILINI^– and T_1-3) and 149 high-risk (VILINI⁺ or T₄) patients with stage II CRC. As the number of patients in the low-risk (VILINI^– and T_1–3) group was too low, we only analyzed the high-risk (VILINI⁺ or T₄) group using our prognostic models. The analysis showed that our nomogram could significantly identify a group of patients with very good DFS and OS within clinically high-risk stage II CRC ( Supplementary Figures S4C, D ). These results suggest that our new multidimensional models could improve patient prognosis prediction.

To determine the clinical usefulness of our nomograms, a decision curve analysis for the nomogram based on our model and pTNM stage was performed, as shown in Figure 6 . By applying our prognostic models, a higher net benefit than that for the strategy of accepting or rejecting interventions for every patient could be achieved when the risk thresholds for DFS range from 12% to 100% at 3 years ( Figure 6A ) and from 16% to 100% at 5 years ( Figure 6B ) and for OS from 6% to 80% at 3 years ( Figure 6C ) and from 15% to 100% at 5 years ( Figure 6D ). Especially, at 3 years, if the threshold probability ranged from 0.12 to 0.30 and from 0.43 to 1.00, our prognostic models for DFS showed a better net benefit than that of pTNM stage ( Figure 6A ); at 5 years, if the threshold probability ranged from 0.16 to 0.36 and from 0.53 to 1.00, our prognostic models for DFS showed a better net benefit than that of pTNM stage ( Figure 6B ). As for the prognostic models for OS at 3 years, if the threshold probability ranged from 0.06 to 0.22 and from 0.27 to 0.80, our nomogram showed a better net benefit than that of pTNM stage ( Figure 6C ); at 5 years, our nomogram showed a better net benefit than that of pTNM stage if the threshold probability ranged from 0.15 to 0.37 and from 0.47 to 1.00 ( Figure 6D ).

Based on the two prognostic models of the 254 patients from our cohort, we developed a website for predicting the risks of recurrence (DFS) and death (OS) of stage II–III CRC patients (http://www.biostatistics.online/liuyuan2). The scoring system based on our model was built into the website, and the prediction results, including risk stratification (high/intermediate/low) and recurrence/survival probabilities at different times, could be obtained after inputting the model variable values. It is easy to operate and friendly to clinicians, which is very helpful for the generalization and application of our models.