Clinical Significance and Inflammatory Landscape of aNovel Recurrence-Associated Immune Signature in Stage II/III Colorectal Cancer

Background: A considerable number of patients with stage II/III colorectal cancer (CRC) will relapse within 5 years after surgery, which is a leading cause of death in early-stage CRC. The current TNM stage system is limited due to the heterogeneous clinical outcomes displayed in patients of same stage. Therefore, searching for a novel tool to identify patients at high recurrence-risk for improving post-operative individual management is an urgent need.

Methods: Using four independent public cohorts and qRT-PCR data from 66 tissues, we developed and validated a recurrence-associated immune signature (RAIS) based on global immune genes. The clinical and molecular features, tumor immune microenvironment landscape, and immune checkpoints profiles of RAIS were also investigated.

Results: In five independent cohorts, this novel scoring system was proven to be an independent recurrent factor and displayed excellent discrimination and calibration in predicting the recurrence-risk at 1~5 years. Further analysis revealed that the high-risk group displayed high mutation rate of TP53, while the low-risk group had more abundance of activated CD4+/CD8+ T cells and high expression of PD-1/PD-L1.

Conclusions: The RAIS model is highly predictive of recurrence in patients with stage II/III CRC, which might serve as a powerful tool to further optimize decision-making in adjuvant chemotherapy and immunotherapy, as well as tailor surveillance protocol for individual patients.

Introduction

By 2020, colorectal cancer (CRC) has become the third most prevalent cancer and second leading cause of cancer-related mortality globally (1). The American Joint Committee on Cancer (AJCC) and Union for International Cancer Control (UICC) TNM stage system is usually applied to manage the prognosis and adjuvant chemotherapy of early-stage CRC (2). 5-Fluorouracil (5-FU)-based adjuvant chemotherapy aims to eliminate residual cancer cells after surgical resection, thus reducing the recurrence rate or extending the time to recurrence. It is a conventional therapy for stage III and a subset of stage II CRC patients (e.g., T4, high grade) (3, 4). However, administration of adjuvant chemotherapy is unsatisfactory in clinical practice, with many patients not benefiting from this because they are either cured by surgery or relapse after adjuvant chemotherapy (5). Moreover, current guidelines demonstrate that the present definition of “high-risk” stage II CRC remains insufficient (6). Previous studies have shown that approximately half of stage III CRC patients will relapse within 5 years after surgical resection, while the 5-years recurrence rate of stage II CRC patients is about 12~38% (7–9). In clinical practice, the AJCC stage system alone is limited in stratifying these patients (10), thus, searching for new ways to identify patients at high risk for recurrence who have urgent need for additional therapeutic schemes is imperative.

Indeed, substantial efforts have been put into identifying and developing markers for assessing the recurrence-risk in early-stage CRC. The consensus molecular subtype (CMS) system has been reported to correlate with clinical outcomes in stage II/III CRC, CMS4 tumors have dismal recurrence and overall survival (11). In parallel, the transcriptomic-based CRC intrinsic subtype (CRIS) that reveals cancer-cell intrinsic characteristics has been proven to be an independent factor of prognosis, and CRIS-C tumors are superior at risk of recurrence (12). Patients with high level of microsatellite instability (MSI-H) are more prone to show a dramatically reduced risk of death and recurrence (13). Recently, circulating tumor DNA analyses have been demonstrated to serve as biomarkers of recurrence and benefit of adjuvant therapy in stage III CRC (14). In addition, alterations in various genes, such as BRAF, KRAS, and PIK3CA mutations, SMAD4 deletion, and HER2 amplification, are also significantly associated the recurrence after curative surgery of early-stage CRC (15–19). We have previously reported that TTN/OBSCN “Double-hit” suggests favorable prognosis in CRC (20). Nonetheless, these markers only have a moderate prediction accuracy and limited clinical usefulness (21, 22).

In recent years, tumor immune microenvironment (TIME) markers have displayed encouraging progress in evaluating prognostic and immunotherapy of patients (22–24). In-situ adaptive immune cells infiltration (e.g., CD8+ cytotoxic T lymphocyte) is linked with the time to recurrence (25). Actually, an immunohistochemistry-based scoring pipeline has been established and validated (termed Immunoscore^®), which quantifies the densities of two adaptive immune cells, CD3+ and CD8+ T cells, in the core and invasive margin of tumor (22). Although Immunoscore^® displays a stable predictive power of prognosis in early-stage CRC, and has been introduced as a diagnostic standard in the 5^th edition of the WHO Digestive System Tumors (26), its performance remains at a moderate accuracy of Harrell’s C-statistics ranging from 0.56 to 0.68 in international researches (22). This may be due to the fact that only two adaptive immune cells are considered, but other immune components have been also reported to correlate with the time to death and recurrence, such as macrophages, nature killer (NK) cells, and γδ T lymphocytes (27–29). Therefore, estimating the recurrence-risk of patients with stage II/III CRC based on the global immune milieu might improve the accuracy of model. Traditional techniques such as immunohistochemistry or quantitative real-time polymerase chain reaction (qRT-PCR) are difficult to quantify a broad spectrum of immune genes, but advances in bioinformatics and machine learning have made it possible. With the help of machine learning, such as the least absolute shrinkage and selection operator (LASSO) algorithm, it is easy to identify the most important elements based on the expression profiles of global immune genes and fit a model with strong generalization performance (30).

In this study, using four independent public cohorts with gene expression and recurrence-free survival (RFS) data, we developed and validated an individualized recurrence-associated immune signature (RAIS) for stage II/III CRC after surgical resection. The clinical and molecular features, TIME landscape, and immune checkpoints profiles of RAIS were also investigated. Furthermore, we used 66 frozen tissue samples with qRT-PCR data for experimental verification to prove the stability and reliability of the RAIS model. Initial construction of the RAIS for stratifying the recurrence-risk will enhance the understanding of underlying mechanisms between immune molecules and the recurrence of CRC, and might help optimize decision-making in adjuvant chemotherapy and immunotherapy for patients with stage II/III CRC.

Materials and Methods

Publicly Available mRNA Data and Immune Gene Sets

We retrospectively collected four CRC cohorts from The Cancer Genome Atlas (TCGA, https://portal.gdc.cancer.gov) and Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo) databases, including TCGA-CRC, GSE143985, GSE29621, and GSE92921. The data processing was described in detail in the Supplementary Method. In four cohorts, we only retained CRC patients that met the following criteria: (1) Have mRNA expression data; (2) Have both recurrent status and RFS information; (3) In the AJCC stage II/III; (4) No preoperative chemotherapy or radiotherapy received. A total of 171 patients from TCGA-CRC were used as the training set, and GSE143985 (n =91), GSE29621 (n =40), and GSE92921 (n =59) from the GEO database were used as the validation sets. The corresponding clinical information of four cohorts was also downloaded, and the baseline data were summarized in Table S1. The list of immune-related genes (IRGs) was retrieved from ImmPort (https://www.immport.org), IRIS (http://www.immunegene.org), and InnateDB (https://www.innatedb.com) databases. There were a total of 2504 IRGs that could be detected in all four cohorts.

Signature Generation

First, based on univariate Cox regression, we identified stable recurrence-associated IRGs using the following criteria: (1) P <0.05 in at least half of the enrolled cohorts; (2) all hazard ratio (HR) greater or less than 1 in the cohorts with statistical significance. Second, using the expression of these stable recurrence-associated IRGs in TCGA-CRC, we fitted a LASSO Cox regression model for assessing the RFS of patients with stage II/III CRC. Using the 10-fold cross validations, the optimal lambda was obtained when the partial likelihood deviance reached the minimum value. Third, based on the optimal lambda, the IRGs with nonzero coefficients were selected to construct the prediction signature. The risk score for each patient was calculated with the LASSO model weighting coefficient as follows:

$$Risk score={\displaystyle \sum _{i=1}^{n}Ex{p}_i}\times Coe{f}_i$$

where n is the number of key IRGs, Exp_i is the expression of IRG i, and Coef_i is the LASSO coefficient of IRG i.

TIME Characterization Analysis

Gene expression profiles were utilized to decipher the TIME characterization of CRC samples with multiple bioinformatics tools. The immunophenoscore (IPS) was applied to assess the immune state of each sample. IPS is a scoring scheme that quantifies the immunogenicity of tumor samples using a variety of markers of immune response or immune tolerance. The higher the IPS z-score, the stronger the immunogenicity of the sample (31). To measure the infiltration abundance of immune cell populations in tumor tissues, two different tools were applied. First, we leveraged the R package ESTIMATE to infer the fraction of stromal and immune cells in CRC samples, and generated two scores including the immune and stromal scores (32). Next, to describe a more detailed landscape of immune cell types infiltration, we employed the single sample gene set enrichment analysis (ssGSEA) implemented in R package GSVA. The gene sets for marking each cell were obtained from the research of Charoentong, which stored 28 human immune cell subtypes (31).

In addition, to predict their putative response to immune checkpoint blockade (ICB), CRC samples were scored using T-cell inflammatory signature (TIS) and Tumor Immune Dysfunction and Exclusion (TIDE) approaches. TIS proposed by Ayers et al. was used to predict clinical response to PD-1 blockade. The signature was composed of 18 inflammatory genes associated with antigen presentation, chemokine expression, cytotoxic activity, and adaptive immune resistance (33). The TIDE algorithm (http://tide.dfci.harvard.edu/) integrates the expression signature of two primary mechanisms of immune evasions: T cell dysfunction and T cell exclusion, to model tumor immune evasion. Patients with higher TIDE score suggest the greater potential of tumor immune evasion; thus, these patients would derive worse immunotherapy response (34).

Patients and Tissue Specimens

From January 2015 to December 2016, we collected a total of 66 frozen surgically resected CRC tissues with AJCC stage II/III at The First Affiliated Hospital of Zhengzhou University. Follow up was concluded five years after surgery, and all patients gave written informed consent. Detailed baseline data of CRC patients were displayed in Table S1. Total RNA was isolated from CRC tissues using RNAiso Plus reagent (Takara, Dalian, China) according to the manufacturer’s instructions. RNA quality was evaluated using a NanoDrop One C (Waltham, MA, USA), and RNA integrity was assessed using agarose gel electrophoresis. An aliquot of 1 μg of total RNA was reverse transcribed into complementary DNA (cDNA) according to the manufacturer’s protocol using a High-Capacity cDNA Reverse Transcription kit (TaKaRa BIO, Japan). All cDNA samples were prepared for qRT-PCR. This project was approved by the Ethics Committee Board of The First Affiliated Hospital of Zhengzhou University.

qRT-PCR Analysis

In the qRT-PCR analysis, the enrolled 12 genes in the RAIS signature and feature genes (including PD-L1, PD-1, CD4 and CD8A) were detected. qRT-PCR was performed using SYBR Assay I Low ROX (Eurogentec, USA) and SYBR^® Green PCR Master Mix (Yeason, Shanghai, China). The expression value of the target genes was normalized to GAPDH, and then log2 transformed for subsequent analysis. The primer sequences of the included 12 genes and GAPDH were shown in Table S2.

Statistical Analysis

All data processing, statistical analysis, and plotting were conducted in R 4.0.2 software. The time from surgery to cancer recurrence was defined as RFS. Correlations between two continuous variables were assessed via Spearman’s correlation coefficients. Fisher’s exact test or Pearson’s chi-squared test was applied to compare categorical variables. Continuous variables were compared between two groups through the Wilcoxon rank-sum test or T test. The patients were divided into high and low-risk groups based on the median risk score. The Kaplan-Meier method and the log-rank test were utilized to estimate the different RFS between two groups. The receiver operating characteristic curve (ROC) used to predict binary categorical variables was implemented using the R package pROC. The time-dependent area under the ROC (AUC) for survival variable were conducted by the R package timeROC. The R package rms was applied to plot calibration curves. GO and KEGG enrichment analysis were performed by the R package clusterProfiler. All statistical tests were two-sided. P <0.05 was regarded as statistically significant.

Results

Patients’ Demographic Characteristics

As illustrated in Table S1, we collected a total of 427 patients with stage II/III CRC from five independent multicenter cohorts. In TCGA-CRC, there were 101 stage II and 70 stage III patients with a median RFS of 2.2521 years; the 1-, 2-, and 3-year recurrence rate were 5.3%, 8.2%, and 10.5%, respectively. In GSE143985, there were 55 stage II and 36 stage III patients with a median RFS of 5.8904 years; the 1-, 2-, and 3-year recurrence rate were 9.9%, 12.1%, and 15.4%, respectively. In GSE29621, there were 22 stage II and 18 stage III patients with a median RFS of 3.9379 years; the 1-, 2-, and 3-year recurrence rate were 7.5%, 10.0%, and 15.0%, respectively. In GSE92921, there were 43 stage II and 16 stage III patients with a median RFS of 5.7260 years; the 1-, 2-, and 3-year recurrence rate were 6.8%, 8.5%, and 10.2%, respectively. In qRT-PCR validation cohort, there were 40 stage II and 26 stage III patients with a median RFS of 3.9671 years; the 1-, 2-, and 3-year recurrence rate were 15.2%, 16.7%, and 21.2%, respectively.

Construction and Validation of the RAIS With Stage II/III CRC in Public Datasets

In order to comprehensively evaluate the relationships between IRGs and recurrence in patients with stage II/III CRC, an analysis pipeline was developed to display the process (Figure 1). Univariate Cox results of four cohorts revealed a total of 88 stable IRGs were significantly associated with RFS (Table S3). GO and KEGG enrichment analysis indicated these genes were mainly involved in immune response and inflammatory signaling pathways (Figure S1). Based on the expression of these genes in TCGA-CRC, we fitted a LASSO Cox regression model and identified 12 IRGs that were strongly predictive of RFS, including MFI2, LTA, VEGFA, NPY, SHC3, RAG1, CASP1, NTF3, COCH, NMB, ERN1, and NLRP7 (Figure 2A). Next, a risk score for RAIS was calculated using a formula that including the 12 IRGs weighted by their regression coefficients in a penalized Cox model as follows: Risk score = 0.4906 × Exp(VEGFA) + 0.2516 × Exp(MFI2) + 0.2276 × Exp(RAG1) + 0.1638 × Exp(COCH) + 0.0970 × Exp(NPY) + 0.0369 × Exp(ERN1) - 0.0715 × Exp(CASP1) - 0.0967 × Exp(LTA) - 0.1198 × Exp(NLRP7) - 0.1547 × Exp(SHC3) - 0.3043 × Exp(NTF3) - 0.3210 × Exp(NMB) (Figure 2B). This formula was utilized to calculate the risk score for each patient in four cohorts. The expression heatmap of the 12 selected IRGs and the distribution of risk scores were illustrated in Figure 2C. All patients were assigned into high- and low-risk groups according to the median risk score (Figure 2C). Relative to the low-risk group, patients in the high-risk group had significantly dismal RFS in TCGA-CRC (HR =4.051, 95% confidence interval (95%CI) =2.497~6.571, log-rank P =0.00043), GSE143985 (HR =5.591, 95%CI =2.910~10.743, log-rank P =0.00012), GSE29621 (HR =2.317, 95%CI =1.536~4.854, log-rank P =0.0063), and GSE92921 (HR =6.823, 95%CI =2.568~18.133, log-rank P =0.0078) (Figures 3A, B). After controlling the available clinical characteristics in four cohorts, multivariate Cox regression analysis revealed RAIS remained an independent risk factor for evaluating RFS of stage II/III CRC patients (all P <0.05) (Figure 3B).

FIGURE 1

Figure 1 The flowchart of this study.

FIGURE 2

Figure 2 The development of the RAIS model based on the LASSO algorithm. (A) Ten-fold cross-validations to tune the parameter selection in the LASSO model. The two dotted vertical lines are drawn at the optimal values by minimum criteria (left) and 1−SE (standard error) criteria (right). (B) LASSO coefficient profiles of the candidate genes for RAIS construction. (C) The distribution of risk score, recurrence status, and gene expression panel in four cohort.

FIGURE 3

Figure 3 Survival significance of RAIS in four cohorts. (A) Kaplan-Meier curves of RFS according to the RAIS. (B) Univariate and multivariate Cox regression analysis of the risk score. The bold values mean P <0.05.

Evaluation of the RAIS Model

In this study, we evaluated this model from two perspectives: discrimination and calibration. The discrimination was assessed by ROC and Harrell’s C-index, and the calibration was assessed by calibration plots. The results showed that the AUCs for predicting RFS at 1~5 years was 0.783, 0.841, 0.858, 0.859, and 0.951 in TCGA-CRC, 0.941, 0.922, 0.880, 0.878, and 0.877 in GSE143985, 0.892, 0.879, 0.840, 0.735, and 0.785 in GSE29621, and 0.962, 0.973, 0.961, 0.960, and 0.959 in GSE92921, respectively (Figure 4A). The C-index were 0.814 (95%CI: 0.734~0.895), 0.860 (95%CI: 0.790~0.930), 0.825 (95%CI: 0.703~0.947), and 0.938 (95%CI: 0.883~0.993) in four cohorts, respectively (Figure 4B). The above indicated the high predictive accuracy of this model. Furthermore, the RAIS displayed excellent calibration, with the predicted probabilities of RFS at 1~5 years accurately, describing the true risk observed in all four cohorts (Figure 4C). The RAIS also can accurately separate the recurrence and recurrence-free CRC with tumor stage II/III after surgical resection. As illustrated in Figure 4D, patients in the high-risk group displayed a significantly higher fraction of recurrence (high-risk vs. low-risk: 22% vs. 6% in TCGA-CRC, 31% vs. 2% in GSE143985, 35% vs. 0% in GSE29621, and 21% vs. 0% in GSE92921; all P <0.05). The ROC analysis further suggested the RAIS possessed high accuracy for identifying CRC patients with recurrence in all four cohorts (Figure 4E). Taken together, the RAIS signature presented stable and excellent performance in evaluating RFS in patients with stage II/III CRC after surgical resection.

FIGURE 4

Figure 4 Evaluation of the RAIS model in four cohorts. (A) Time-dependent ROC analysis for predicting RFS at 1~5 years. (B) The Harrell’s C-index of RAIS. (C) Calibration plots for comparing the actual probabilities and the predicted probabilities of RFS at 1~5 years. (D) Comparison of recurrence rate between the high-risk and low-risk groups. (E) ROC analysis of the RAIS model for predicting the recurrence event of patients.

Clinical and Molecular Characteristics of RAIS

In order to characterize the clinical significance of RAIS, we explored the relationship between clinical features and RAIS. As shown in Figures S2–S5, most clinical features, including age, sex, tumor stage, microsatellite instability, adjuvant chemotherapy, and KRAS mutation, were not significantly different between the high and low-risk groups in all four cohorts. In CRC, the TP53 mutation is a key step driving the transition from adenoma to adenocarcinoma, and it is associated with adverse clinical outcomes (35). Obviously, in four cohorts, patients in the high-risk group had consistently higher proportion of TP53 mutation compared with the low-risk group (all P <0.05) (Figure S2, S3, and S5), and the risk score was predominantly higher in patients with TP53 mutation (all P <0.01) (Figure S6).

TIME Landscape and Immune Checkpoints Profiles of RAIS

Since the establishment of RAIS was based on immune-related genes, we hypothesized that there were differences in the immune characterization between two recurrence-risk groups. First, the IPS score was utilized as a general index of immune activation in tumor tissues. Patients in the low-risk group displayed a higher IPS z-score compared with the high-risk group (all P <0.05) (Figure S7A). Next, ESTIMATE software was employed to infer the fraction of stromal and immune cells. In line with the IPS results, the low-risk group scored better in the immune category (all P <0.05) (Figure S7B). To gain more detailed insights into this issue, we applied the ssGSEA approach to quantify the infiltration abundance of different immune cell populations. Overall, the immune cells infiltration was more abundant in the low-risk group (Figure S8A), suggesting that their immune activity and immune response would be more active. Also, we could observe the risk score had broad negative relationships with different immune cell types (Figure 5A). Of note, there were significant and stable correlations between risk score and activated CD4+ T cell and CD8+ T cell in four cohorts (all P <0.01) (Figure 5A). Patients in the low-risk group possessed higher infiltration abundance of activated CD4+ T cell and CD8+ T cell (all P <0.05) (Figure 5B).

FIGURE 5

Figure 5 TIME landscape and immune checkpoints profiles of RAIS in four cohorts. (A) The correlation analysis between RAIS and 28 immune cells infiltration abundance. (B) The distribution difference of activated CD4+/CD8+ T cells infiltration between the high-risk and low-risk groups. (C) Four heatmaps of 27 immune checkpoints profiles in high-risk and low-risk groups. *P < 0.05, **P < 0.01, ***P < 0.001.

Furthermore, we extended our analysis to encompass 27 immune checkpoint members, including B7-CD28 family (PD-L1, PD-L2, PD-1, CTLA4, CD276, HHLA2, ICOS, ICOSLG, TMIGD2, and VTCN1) (36), the TNF superfamily (BTLA, CD27, CD40, CD40LG, CD70, TNFRSF18, TNFRSF4, TNFRSF9, and TNFSF14) (37), and several other molecules (ENTPD1, FGL1, HAVCR2, IDO1, LAG3, NCR3, NT5E, and SIGLEC15) (38, 39), and the results for four cohorts were presented in Figure 5C. In total, we observed that two molecules, including PD-1 and PD-L1, were significantly upregulated in the low-risk group in four cohorts (all P <0.05) (Figure S8B). Consistent with this, the risk score was negatively correlated with the expression of PD-1 and PD-L1 (all P <0.05) (Figure S8C). Previous study has demonstrated patients with high expression of PD-1 and PD-L1 would benefit more from pembrolizumab and nivolumab (40, 41). Hence, to predict their putative response to ICB, two bioinformatics tools, TIS and TIDE, were used. As shown in Figure S9, patients in the low-risk group displayed the higher TIS score and the lower TIDE score (all P <0.05), suggesting that they were more likely to yield considerable clinical benefit from ICB therapy.

Validation of RAIS in an Independent Cohort From Frozen Tissue Samples

To further verify the performance of our 12-gene RAIS model into a clinically translatable tool, we next evaluated the expression of these genes in a clinical cohort of 66 CRC patients by conducting qRT-PCR assay. The expression heatmap of the 12 selected IRGs and the distribution of risk scores were illustrated in Figure S10A. Consistently, the Kaplan-Meier survival analysis demonstrated that patients with high score displayed the dramatically poor RFS (HR =2.299, 95%CI =1.478~3.577, log-rank P <0.0001) (Figures 6A, B). After controlling for confounding variables (including age, sex, stage, and postoperative chemotherapy), the RAIS remained the statistical significance (HR = 2.190, 95%CI =1.043~4.598, P <0.05) (Figure 6B). ROC analysis showed the pinpoint accuracy of RAIS: the AUCs for predicting RFS at 1~5 years was 0.959, 0.950, 0.924, 0.891, and 0.900, respectively (Figure 6C). Similarly, the C-index reached 0.939 (95%CI =0.884~0.994). The calibration plot further displayed the predicted probabilities of RFS at 1~5 years accurately describing the true risk observed (Figure S10B). In addition, we also found that patients in the low-risk group might have more abundance of T cells (CD4+ and CD8+), and higher expression of PD-1 and PD-L1 (Figure S10C). Taken together, the results from a clinical in-house cohort supported that our discovery and in-silico validation cohort findings, which validated and confirmed that our RAIS model was quite robust, and can serve as an independent predictor of recurrence in stage II/III CRC.

FIGURE 6

Figure 6 Validation of our discovery in a clinical in-house cohort. (A) Kaplan-Meier curves of RFS according to the RAIS. (B) Univariate and multivariate Cox regression analysis of the risk score. (C) Time-dependent ROC analysis for predicting RFS at 1~5 years.

Discussion

CRC is a highly heterogeneous tumor that possesses complex biological processes, for which post-operative surveillance and therapeutic regimens are necessary to be tailored to generate an optimal outcome for each patient. Nevertheless, a considerable proportion of stage II/III CRC patients not only derive benefit from 5-FU-based adjuvant chemotherapy but also display drug reactions (5). A limitation of the current AJCC stage system is that patients in the same stage have distinct clinical outcomes, which leads to latent under- or over-treatment. Herein, developing a novel classifier that can be routinely implemented into clinical practice is critical for identifying those early-stage patients who are at high recurrence-risk and who might thus benefit from adjuvant chemotherapy. Considering various immune responses have close connections with the occurrence and progression of CRC (23, 25), we believe that introducing immune-based parameters in the clinical management of early-stage CRC is a high priority. Despite an international study recently proposed an index termed Immunoscore^® that presented stable and independent predictive of prognosis, it only considered two adaptive immune cells and showed a moderate accuracy. We hypothesized that a signature with high performance could be developed according to the global immune milieu. With the development of artificial intelligence and bioinformatics, an advanced machine learning algorithm can identify several key indicators that are most meaningful to predict clinical outcomes from a large number of genes (30), which is actually in line with the biological scale-free network which was dominated by a few hub nodes (42). Therefore, for the first time, we developed a novel signature (termed RAIS) to evaluate the recurrence-risk of patients with stage II/III CRC in multicenter cohorts based on the RNA expression of global immune genes. The reproducibility and powerful performance of RAIS in multiple independent cohorts and external qRT-PCR data not only prove that it is a robust and highly accurate model, but also is promising to be routinely implemented into clinical practice due to the advantages of high sensitivity and specificity, simplicity, and low cost of qRT-PCR.

In this study, we fitted a recurrence model consisting of 12 IRGs, including MFI2, LTA, VEGFA, NPY, SHC3, RAG1, CASP1, NTF3, COCH, NMB, ERN1, and NLRP7. All genes have been reported to be involved in the progression and TIME cross-talking of tumor. For example, VEGFA induces the expression of transcription factor TOX to drive T cell exhaustion (43), and the expression of CASP1 is able to be repressed by G9A and further promotes tumor immune escape (44). Based on the 12 enrolled genes, we developed the RAIS model, which performed stably in predicting recurrent-risk of patients with stage II/III CRC. The prognostic meta-analysis showed that RAIS was an extremely vicious indicator of recurrence and was proven to be an independent factor after adjusting multiple clinical clinicopathologic features. More importantly, in four cohorts, RAIS demonstrated a high discrimination and calibration in predicting the recurrence-risk at 1~5 years. To prevent false positive results from sequencing data, we conducted another validation according to qRT-PCR results from 66 frozen CRC tissues with tumor stage II/III, confirming our prior findings and evaluating their practicality in different centers. As reported previously, patients with a high-risk score suggested dismal RFS, and thus might need to adjust therapy strategies or add additional adjuvant chemotherapy. For example, current guidelines recommend that a subset of stage II patients without “high-risk” traits do not require adjuvant chemotherapy (6), but when these patients show a high-risk score, using additional adjuvant chemotherapy might be essential.

Afterwards, we conducted a comprehensive analysis in the relationships between RAIS and clinical and molecular traits in four cohorts. RAIS was found to have no significant association with most of available traits. Of note, the risk score was dramatically higher in the TP53 mutant CRC compared with wild type tumors. Based on the previous literature, patients with TP53 mutation tended to display more aggressiveness and higher recurrence rate, which was in line with our finding (35). Strikingly, the interesting finding in this study was the relationship between RIAS and the prevalent immunotherapeutic biomarkers in stage II/III CRC. As is well-known, cancer immunotherapy represented by ICB has revolutionized the treatment of solid tumors, including a subset of CRC. Two monoclonal antibodies targeting PD-1, nivolumab and pembrolizumab, have demonstrated great efficacy in CRC with MSI-H mismatch repair deficiency (dMMR), and have been approved by FDA (40). In this study, patients in the low-risk group suggested an “immune-hot” subtype, displaying the rich infiltration of activated CD4+/ CD8+ T cells and the higher expression PD-1/PD-L1. High abundance of tumor-infiltrating lymphocytes (especially T cells) is not only a strong prognostic indicator, but also provides the backup resource for immunotherapy (45). Similarly, previous studies have reported that the high expression of PD-1/PD-L1 correlates with favorable prognosis and better immunotherapeutic response in early-stage CRC (40, 41, 46). Hence, patients with the low-risk score suggested a potential benefit from immunotherapy. Two bioinformatic algorithms, TIS and TIDE, also showed the low-risk group might be more sensitive to immunotherapy relative to the high-risk group in four cohorts. Overall, patients with high-risk score might be not suitable for immunotherapy, due to potential ineffectiveness and immune-related adverse events (irAEs).

To the best of our knowledge, this is the first and most comprehensive study to date validating the prognostic accuracy of an immune signature in patients with stage II/III CRC undergoing surgical resection, based on the global immune genes. Prior to this study, a few reports established molecular signatures for predicting prognostic risk of CRC (47–50). In comparison with these studies, our work has several advantages and novelties: (1) The RAIS model was developed based on the recurrence rather than overall survival in patients with stage II/III CRC, which allowed it to accurately identify high-risk patients with early-stage CRC; (2) We performed comprehensive statistical approaches to evaluate the discrimination and calibration of the RAIS model, and our model remained stable and highly accurate performance at 1~5 years; (3) qRT-PCR was used to validate the performance of RAIS to ensure its robustness and clinical feasible; (4) We also demonstrated the TIME profiles and immune checkpoint landscape of RAIS, revealing its potential predictive value of immunotherapy. Despite the RAIS model is promising, some limitations should be acknowledged. First, all the samples from five centers were retrospective, and future validation of the RAIS model should be conducted in prospective fresh samples. Second, some clinical characteristics on public datasets were very inadequate, which thus had concealed the potential associations between RAIS and some clinical traits. Third, patients treated with immunotherapy were not examined in this study, so the performance of RAIS for predicting immunotherapeutic response was investigated indirectly. Further prospective study is still necessary.

In conclusion, we established a reproducible and powerful model for evaluating the recurrence-risk of patients with stage II/III CRC. Our study provides novel implications regarding immune profiles and stage II/III CRC recurrence. More importantly, the RAIS model may be a promising tool to optimize decision-making in adjuvant chemotherapy and immunotherapy, as well as tailor surveillance protocol for individual patients with stage II/III CRC.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

Ethics Statement

The studies involving human participants were reviewed and approved by Ethics Committee of The First Affiliated Hospital of Zhengzhou University [December 19, 2019, TRN: 2019-KW-423]. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

ZL and XH designed this work. ZL, QD, ZS, JL, KX, and DJ integrated and analyzed the data. ZL, TL, LW, CG, and LL wrote this manuscript. ZL, JL, KX, ZS, and XH edited and revised the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This study was supported by the National Natural Science Foundation of China (81972663), Henan Province Young and Middle‐Aged Health Science and Technology Innovation Talent Project (YXKC2020037), and Henan Provincial Health Commission Joint Youth Project (SB201902014).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2021.702594/full#supplementary-material

References

1. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin (2021) 71:209–49. doi: 10.3322/caac.21660

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Weiser MR. AJCC 8th Edition: Colorectal Cancer. Ann Surg Oncol (2018) 25(6):1454–5. doi: 10.1245/s10434-018-6462-1

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Kannarkatt J, Joseph J, Kurniali PC, Al-Janadi A, Hrinczenko B. Adjuvant Chemotherapy for Stage II Colon Cancer: A Clinical Dilemma. J Oncol Pract (2017) 13(4):233–41. doi: 10.1200/JOP.2016.017210

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Taieb J, Gallois C. Adjuvant Chemotherapy for Stage III Colon Cancer. Cancers (Basel) (2020) 12(9):2679–95. doi: 10.3390/cancers12092679

CrossRef Full Text | Google Scholar

5. Moertel CG, Fleming TR, Macdonald JS, Haller DG, Laurie JA, Goodman PJ, et al. Levamisole and Fluorouracil for Adjuvant Therapy of Resected Colon Carcinoma. N Engl J Med (1990) 323(3):197–8. doi: 10.1056/NEJM199007193230312

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Benson AB, Venook AP, Al-Hawary MM, Cederquist L, Chen YJ, Ciombor KK, et al. NCCN Guidelines Insights: Colon Cancer, Version 2.2018. J Natl Compr Canc Netw (2018) 16(4):359–69. doi: 10.6004/jnccn.2018.0021

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Carlsson U, Lasson A, Ekelund G. Recurrence Rates After Curative Surgery for Rectal Carcinoma, With Special Reference to Their Accuracy. Dis Colon Rectum (1987) 30(6):431–4. doi: 10.1007/BF02556491

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Dekker E, Tanis PJ, Vleugels JLA, Kasi PM, Wallace MB. Colorectal Cancer. Lancet (2019) 394(10207):1467–80. doi: 10.1016/S0140-6736(19)32319-0

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Osterman E, Hammarstrom K, Imam I, Osterlund E, Sjoblom T, Glimelius B. Recurrence Risk After Radical Colorectal Cancer Surgery-Less Than Before, But How High Is It? Cancers (Basel) (2020) 12(11):3308–38. doi: 10.3390/cancers12113308

CrossRef Full Text | Google Scholar

10. Benson AB 3rd, Hamilton SR. Path Toward Prognostication and Prediction: An Evolving Matrix. J Clin Oncol (2011) 29(35):4599–601. doi: 10.1200/JCO.2011.37.8646

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Guinney J, Dienstmann R, Wang X, de Reynies A, Schlicker A, Soneson C, et al. The Consensus Molecular Subtypes of Colorectal Cancer. Nat Med (2015) 21(11):1350–6. doi: 10.1038/nm.3967

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Allen WL, Dunne PD, McDade S, Scanlon E, Loughrey M, Coleman H, et al. Transcriptional Subtyping and CD8 Immunohistochemistry Identifies Poor Prognosis Stage II/III Colorectal Cancer Patients Who Benefit From Adjuvant Chemotherapy. JCO Precis Oncol (2018) 2018:1–15. doi: 10.1200/PO.17.00241

CrossRef Full Text | Google Scholar

13. Petrelli F, Ghidini M, Cabiddu M, Pezzica E, Corti D, Turati L, et al. Microsatellite Instability and Survival in Stage II Colorectal Cancer: A Systematic Review and Meta-Analysis. Anticancer Res (2019) 39(12):6431–41. doi: 10.21873/anticanres.13857

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Tie J, Cohen JD, Wang Y, Christie M, Simons K, Lee M, et al. Circulating Tumor DNA Analyses as Markers of Recurrence Risk and Benefit of Adjuvant Therapy for Stage III Colon Cancer. JAMA Oncol (2019) 5:1710–7. doi: 10.1001/jamaoncol.2019.3616

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Margonis GA, Buettner S, Andreatos N, Kim Y, Wagner D, Sasaki K, et al. Association of BRAF Mutations With Survival and Recurrence in Surgically Treated Patients With Metastatic Colorectal Liver Cancer. JAMA Surg (2018) 153(7):e180996. doi: 10.1001/jamasurg.2018.0996

PubMed Abstract | CrossRef Full Text | Google Scholar

16. McCuaig S, Barras D, Mann EH, Friedrich M, Bullers SJ, Janney A, et al. The Interleukin 22 Pathway Interacts With Mutant KRAS to Promote Poor Prognosis in Colon Cancer. Clin Cancer Res (2020) 26(16):4313–25. doi: 10.1158/1078-0432.CCR-19-1086

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Kato S, Iida S, Higuchi T, Ishikawa T, Takagi Y, Yasuno M, et al. PIK3CA Mutation is Predictive of Poor Survival in Patients With Colorectal Cancer. Int J Cancer (2007) 121(8):1771–8. doi: 10.1002/ijc.22890

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Wasserman I, Lee LH, Ogino S, Marco MR, Wu C, Chen X, et al. SMAD4 Loss in Colorectal Cancer Patients Correlates With Recurrence, Loss of Immune Infiltrate, and Chemoresistance. Clin Cancer Res (2019) 25(6):1948–56. doi: 10.1158/1078-0432.CCR-18-1726

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Richman SD, Southward K, Chambers P, Cross D, Barrett J, Hemmings G, et al. HER2 Overexpression and Amplification as a Potential Therapeutic Target in Colorectal Cancer: Analysis of 3256 Patients Enrolled in the QUASAR, FOCUS and PICCOLO Colorectal Cancer Trials. J Pathol (2016) 238(4):562–70. doi: 10.1002/path.4679

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Liu Z, Wang L, Guo C, Liu L, Jiao D, Sun Z, et al. TTN/OBSCN ‘Double-Hit’ Predicts Favourable Prognosis, ‘Immune-Hot’ Subtype and Potentially Better Immunotherapeutic Efficacy in Colorectal Cancer. J Cell Mol Med (2021) 7:3239–51. doi: 10.1111/jcmm.16393

CrossRef Full Text | Google Scholar

21. Dienstmann R, Villacampa G, Sveen A, Mason MJ, Niedzwiecki D, Nesbakken A, et al. Relative Contribution of Clinicopathological Variables, Genomic Markers, Transcriptomic Subtyping and Microenvironment Features for Outcome Prediction in Stage II/III Colorectal Cancer. Ann Oncol (2019) 30(10):1622–9. doi: 10.1093/annonc/mdz287

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Pages F, Mlecnik B, Marliot F, Bindea G, Ou FS, Bifulco C, et al. International Validation of the Consensus Immunoscore for the Classification of Colon Cancer: A Prognostic and Accuracy Study. Lancet (2018) 391(10135):2128–39. doi: 10.1016/S0140-6736(18)30789-X

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Mlecnik B, Bindea G, Kirilovsky A, Angell HK, Obenauf AC, Tosolini M, et al. The Tumor Microenvironment and Immunoscore are Critical Determinants of Dissemination to Distant Metastasis. Sci Transl Med (2016) 8(327):327ra26. doi: 10.1126/scitranslmed.aad6352

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Liu Z, Zhang Y, Shi C, Zhou X, Xu K, Jiao D, et al. A Novel Immune Classification Reveals Distinct Immune Escape Mechanism and Genomic Alterations: Implications for Immunotherapy in Hepatocellular Carcinoma. J Transl Med (2021) 19(1):5. doi: 10.1186/s12967-020-02697-y

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Mlecnik B, Tosolini M, Kirilovsky A, Berger A, Bindea G, Meatchi T, et al. Histopathologic-Based Prognostic Factors of Colorectal Cancers Are Associated With the State of the Local Immune Reaction. J Clin Oncol (2011) 29(6):610–8. doi: 10.1200/JCO.2010.30.5425

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Quezada-Marin JI, Lam AK, Ochiai A, Odze RD, Washington KM, Fukayama M, et al. Gastrointestinal Tissue-Based Molecular Biomarkers: A Practical Categorisation Based on the 2019 World Health Organization Classification of Epithelial Digestive Tumours. Histopathology (2020) 77(3):340–50. doi: 10.1111/his.14120

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Tachibana T, Onodera H, Tsuruyama T, Mori A, Nagayama S, Hiai H, et al. Increased Intratumor Valpha24-Positive Natural Killer T Cells: A Prognostic Factor for Primary Colorectal Carcinomas. Clin Cancer Res (2005) 11(20):7322–7. doi: 10.1158/1078-0432.CCR-05-0877

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Corvaisier M, Moreau-Aubry A, Diez E, Bennouna J, Mosnier JF, Scotet E, et al. V Gamma 9V Delta 2 T Cell Response to Colon Carcinoma Cells. J Immunol (2005) 175(8):5481–8. doi: 10.4049/jimmunol.175.8.5481

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Edin S, Wikberg ML, Rutegard J, Oldenborg PA, Palmqvist R. Phenotypic Skewing of Macrophages In Vitro by Secreted Factors From Colorectal Cancer Cells. PloS One (2013) 8(9):e74982. doi: 10.1371/journal.pone.0074982

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Tibshirani R. The Lasso Method for Variable Selection in the Cox Model. Stat Med (1997) 16(4):385–95. doi: 10.1002/(sici)1097-0258(19970228)16:4<385::aid-sim380>3.0.co;2-3

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Charoentong P, Finotello F, Angelova M, Mayer C, Efremova M, Rieder D, et al. Pan-Cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep (2017) 18(1):248–62. doi: 10.1016/j.celrep.2016.12.019

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Yoshihara K, Shahmoradgoli M, Martinez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring Tumour Purity and Stromal and Immune Cell Admixture From Expression Data. Nat Commun (2013) 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Ayers M, Lunceford J, Nebozhyn M, Murphy E, Loboda A, Kaufman DR, et al. IFN-Gamma-Related mRNA Profile Predicts Clinical Response to PD-1 Blockade. J Clin Invest (2017) 127(8):2930–40. doi: 10.1172/JCI91190

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Jiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, et al. Signatures of T Cell Dysfunction and Exclusion Predict Cancer Immunotherapy Response. Nat Med (2018) 24(10):1550–8. doi: 10.1038/s41591-018-0136-1

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Raskov H, Soby JH, Troelsen J, Bojesen RD, Gogenur I. Driver Gene Mutations and Epigenetics in Colorectal Cancer. Ann Surg (2020) 271(1):75–85. doi: 10.1097/SLA.0000000000003393

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Janakiram M, Chinai JM, Zhao A, Sparano JA, Zang X. HHLA2 and TMIGD2: New Immunotherapeutic Targets of the B7 and CD28 Families. Oncoimmunology (2015) 4(8):e1026534. doi: 10.1080/2162402X.2015.1026534

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Ward-Kavanagh LK, Lin WW, Sedy JR, Ware CF. The TNF Receptor Superfamily in Co-Stimulating and Co-Inhibitory Responses. Immunity (2016) 44(5):1005–19. doi: 10.1016/j.immuni.2016.04.019

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Chretien S, Zerdes I, Bergh J, Matikas A, Foukakis T. Beyond PD-1/PD-L1 Inhibition: What the Future Holds for Breast Cancer Immunotherapy. Cancers (Basel) (2019) 11(5):628. doi: 10.3390/cancers11050628

CrossRef Full Text | Google Scholar

39. Wang J, Sanmamed MF, Datar I, Su TT, Ji L, Sun J, et al. Fibrinogen-Like Protein 1 Is a Major Immune Inhibitory Ligand of LAG-3. Cell (2019) 176(1-2):334–47 e12. doi: 10.1016/j.cell.2018.11.010

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Ganesh K, Stadler ZK, Cercek A, Mendelsohn RB, Shia J, Segal NH, et al. Immunotherapy in Colorectal Cancer: Rationale, Challenges and Potential. Nat Rev Gastroenterol Hepatol (2019) 16(6):361–75. doi: 10.1038/s41575-019-0126-x

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Payandeh Z, Khalili S, Somi MH, Mard-Soltani M, Baghbanzadeh A, Hajiasgharzadeh K, et al. PD-1/PD-L1-Dependent Immune Response in Colorectal Cancer. J Cell Physiol (2020) 235(7-8):5461–75. doi: 10.1002/jcp.29494

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Zhang B, Horvath S. A General Framework for Weighted Gene Co-Expression Network Analysis. Stat Appl Genet Mol Biol (2005) 4:17. doi: 10.2202/1544-6115.1128 Article17.

CrossRef Full Text | Google Scholar

43. Kim CG, Jang M, Kim Y, Leem G, Kim KH, Lee H, et al. VEGF-A Drives TOX-Dependent T Cell Exhaustion in Anti-PD-1-Resistant Microsatellite Stable Colorectal Cancers. Sci Immunol (2019) 4(41):eaay0555. doi: 10.1126/sciimmunol.aay0555

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Huang T, Zhang P, Li W, Zhao T, Zhang Z, Chen S, et al. G9A Promotes Tumor Cell Growth and Invasion by Silencing CASP1 in non-Small-Cell Lung Cancer Cells. Cell Death Dis (2017) 8(4):e2726. doi: 10.1038/cddis.2017.65

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Havel JJ, Chowell D, Chan TA. The Evolving Landscape of Biomarkers for Checkpoint Inhibitor Immunotherapy. Nat Rev Cancer (2019) 19(3):133–50. doi: 10.1038/s41568-019-0116-x

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Tome M, Pappalardo A, Soulet F, Lopez JJ, Olaizola J, Leger Y, et al. Inactivation of Proprotein Convertases in T Cells Inhibits PD-1 Expression and Creates a Favorable Immune Microenvironment in Colorectal Cancer. Cancer Res (2019) 79(19):5008–21. doi: 10.1158/0008-5472.CAN-19-0086

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Salazar R, Roepman P, Capella G, Moreno V, Simon I, Dreezen C, et al. Gene Expression Signature to Improve Prognosis Prediction of Stage II and III Colorectal Cancer. J Clin Oncol (2011) 29(1):17–24. doi: 10.1200/JCO.2010.30.1077

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Jorissen RN, Gibbs P, Christie M, Prakash S, Lipton L, Desai J, et al. Metastasis-Associated Gene Expression Changes Predict Poor Outcomes in Patients With Dukes Stage B and C Colorectal Cancer. Clin Cancer Res (2009) 15(24):7642–51. doi: 10.1158/1078-0432.CCR-09-1431

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Li J, Zhang J, Hu H, Cai Y, Ling J, Wu Z, et al. Gene Expression Signature to Predict Prognosis and Adjuvant Chemosensitivity of Colorectal Cancer Patients. Cancer Manag Res (2020) 12:3301–10. doi: 10.2147/CMAR.S243490

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Kim SK, Kim SY, Kim JH, Roh SA, Cho DH, Kim YS, et al. A Nineteen Gene-Based Risk Score Classifier Predicts Prognosis of Colorectal Cancer Patients. Mol Oncol (2014) 8(8):1653–66. doi: 10.1016/j.molonc.2014.06.016

PubMed Abstract | CrossRef Full Text | Google Scholar