Abstract:Dear Editor, Colorectal cancer (CRC) is considered as a genetic disease, which arises from the stepwise accumulation of genetic and epigenetic alterations.1, 2 We found a novel mutational signature (MS) that could assist clinicians to select patients who are more suitable for immunotherapy; the risk score (RS) combined with pathological TNM stage could provide comprehensive and precise prognostic information for CRC patients. To explore the genomic basis of tumor variability in the tumor microenvironment of CRC, we integrated single nucleotide variation (SNV) and transcriptome data and collected information from 1133 and 588 CRC patients of the Memorial Sloan Kettering Cancer Center (MSKCC) and The Cancer Genome Atlas (TCGA) databases. In the training (MSKCC) cohort, we identified an MS consisting of 27 genomic variant genes and generated a prognostic model (Figures 1A and 1B). The date showed that the high-risk group has poorer overall survival (OS), which was verified in both MSKCC and TCGA cohorts (Figures 1C and 1D). The Kaplan-Meier survival curve and ROC curve were applied to evaluate the predicting power of the model by using R packages “survival” and “survival ROC.” The ROC curve indicated that the classifier had a good predictive ability (Figure 1E).The univariate and multivariate analyses revealed MS is an independent, unfavorable prognostic factor for CRCs (Figures 1F and 1G; Table S1). A flowchart is shown in Figure S1. The coefficients of the 27 mutated genes are shown in Table S4. In the training cohort, a nomogram was generated to predict the OS of CRC patients (Figure S2A). The predictors included tumor location, M stage, TNM stage, and RS, among which the RS had the highest C-index (Figures S2B and S2C). The clinical figures of CRC patients are listed in Table S2 and S3. To explore the differences of genomic alterations in these two groups, we analyzed the data containing somatic mutations from the TCGA database. First, it revealed a significant enrichment of different mutations between low- and high-risk groups (Figure 2A). The data showed that more than 90% of CRC in the low-risk group had more mutations in the MSKCC cohort (Figure 2B), which means most of the genes with more mutations in the low-risk group. Besides, the low- and high-risk groups had different distribution of the top 10 mutated genes (Figure 2C). Significant enrichment of oncogenic alterations in such genes as BRAF, ZFHX3, and MTOR was found in right-sided tumors and MSI (Microsatellite instability) patients, while oncogenic alteration of APC (Adenomatous Polyposis Coli) was primary found in the left-sided tumors and MSS (Microsatellite stability) patients (Figure 2D). And all the results were consistent in the validation cohort (Figures 2E–2H; Figure S5). The MSI status is critical when considering immunotherapy and chemotherapeutic drugs as options for CRC patients.3, 4 The RS was observed to be significantly associated with the status of MSI/dMMR and other clinical features (Figures 3A and 3B; Figure S6). In line with previous observation, the status of MSI was more common in low-risk group (Figure 3C). Furthermore, we observed that low-risk group exhibited a higher mutations number (Figures 3D–3F). Due to the hypermutation or high mutational load, these patients might have increased neoantigens, leading to increased immune infiltration, and thus might be more sensitive to immunotherapy. We speculate that there is a potential connection between the MS model and the immune environment. To further clarify the relationship between MS and immune-phenotyping, we analyzed SNV and transcriptome data in the TCGA database. The immune activity was determined by analyzing 29 immune-related genesets. These genesets were analyzed using the ssGSEA.5 A heatmap of the infiltration levels and scores of each sample of immune cells in the three subtypes is shown in Figure 4A. A higher expression level of the PD-L1 gene was found in the immunity-H cluster, and the immunity-H cluster was correlated with better survival outcome (Figures 4B and S5A-S5D). The immunity-H cluster was significantly enriched in the low-risk group (Figure 4E). To reconfirm the findings above, we also performed consensus molecular subgroups (CMS) classification,6 which give a more profound biological insight into immunity typing, and has a strong prognostic effect. PD-L1 was highly expressed in CMS1, which is defined by upregulation of immune genes and associated with MSI-h. The data suggested that CMS1 cluster was significantly enriched in the low-risk group (Figures 4D and 4F), and the CMS1 cluster is more likely benefit from PD-L1 inhibitor treatment. To further investigate the potential predictive value of the MS for the immune status, we examined the possible associations between the immune status and RS. TMB score was negatively correlated with the RS, but TMB was positively correlated with immune score (Figures S6A and S6B). Therefore, we postulated that the RS was negatively associated with immune score. We found a low immune score was related to a worse OS outcome (Figure S3C), and the low-risk group may also be associated with a better survival. To figure out the infiltrated immune cell composition in the defined risk groups, we analyzed the expression signature matrix of 22 infiltrated immune cell types in tumor samples using the CIBERSORT test (Figures S3D-S3G). Regarding tumor-infiltrating immune cells in CRC microenvironment, the number of CD4 memory T cells decreased, and the macrophage M0 increased in the high-risk group (Figures S6D-S6G). Immunotherapy has been raised as a novel effective treatment against CRC; however, the current guidelines only based on the TNM stage cannot reflect the information of host immune system response.7, 8 In clinic, MSI-H is an especially good indicator for checkpoint blockade immunotherapy in CRC, but only about 45% of MSI-H CRC patients could benefit from immunotherapy.9 In our prediction model, we have identified a novel MS, which can generate a prognostic tool to effectively classify CRC patients with different OS risks. Moreover, the MS classifier can be used to predict patients who are more suitable for immunotherapy, and a nomogram comprising the MS could help medical staff in directing personalized therapeutic treatment selection for CRC patients. The authors declare that there is no conflict of interest that could be perceived as prejudicing the impartiality of the research reported. Study design: Ping Lan, Xiaosheng He, and Liang Xu. Literature research: Yanyun Lin and Xijie Chen. Data acquisition: Liang Xu, Yanyun Lin, and Xijie Chen. Data analysis/interpretation: Yanyun Lin and Guanman Li. Statistical analysis: Xijie Chen and Zengjie Chi. Manuscript preparation: Liang Xu and Bin Zheng. Manuscript definition of intellectual content: Lisheng Zheng and Bin Zheng. Manuscript editing: Yufeng Cheng, Jiancong Hu, Shuang Guo, and Danling Liu. All data generated or analyzed during the present study are available via the corresponding author on reasonable request. FIGURE S1 Flowchart detailing the procedure of analyzing mutational signature and their correlation with clinical characteristics and immunity, as well as prognostic models of mutated genes. FIGURE S2 Univariate and multivariate COX regression analyses of clinical factors and independence associated with prognosis/nomogram A to predict the risk of overall survival in CRC. (A) Nomogram to predict distant metastasis-free survival. (B) Calibration curves of the nomogram to predict overall survival at 3 and 5 years in the MSKCC cohort. (C) Calibration curves of the nomogram to predict overall survival at 3 and 5 years in the TCGA cohort. FIGURE S3 The mutation rate of low- and high-risk group. (A) The proportion of right tumor was significantly increased in the low-risk group both in the training and validation cohorts. (B) The proportion of stage III-IV was significantly increased in the high-risk group both in the training and validation cohorts. (C) The high-risk group has a significantly higher risk score than the low-risk group both in the training and validation cohorts. FIGURE S4 The mutation rate of low- and high-risk group. (A) The proportion of right tumor was significantly increased in the low-risk group both in the training and validation cohorts. (B) The proportion of stage III-IV was significantly increased in the high-risk group both in the training and validation cohorts. (C) The high-risk group has a significantly higher risk score than the low-risk group both in the training and validation cohorts. FIGURE S5 Composition of immune cell. (A) PCA analysis of the three clusters. (B) Comparison of the expression levels of HLA genes between CRC subtypes (ANOVA test). (C) Kaplan-Meier analysis of three immunity clusters. (D) Comparison of the stromal score, immune score, ESTIMATE score, and tumor purity between CRC subtypes (Mann-Whitney U test). FIGURE S6 Composition of immune cells at low- and high-risk tissues in the TCGA cohort. (A) Correlation analysis between risk score and TMB in the TCGA dataset. (B) Correlation analysis between immune score and TMB in the TCGA dataset. (C) Kaplan-Meier analysis between low and high immune score. (D) Fractions of immune cells in 63 high-risk and 63 low-risk groups in the TCGA dataset. (E) Correlation of immune cells in the TCGA dataset. (F) Comparison of immune cells between high- and low-risk groups in the TCGA dataset. (G) GSEA analysis of high- and low-risk group. FIGURE S7 DCA curve for MSKCC and TCGA cohorts. (A) DCA curve at 3 and 5 years in the MSKCC cohort. (B) Enriched GO terms in the “biological process” category. Different colors indicate a different function. (C) Enriched KEGG biological pathways. Table S1 Univariate and multivariate COX regression analyses of clinical factors and independence associated with prognosis. Table S2 Relationship between risk score and clinicopathological factors in the training, external validation cohorts. Table S3 Clinical characteristics. Table S4 Twenty-seven mutated genes and its coefficient. Table S5 Relationship between risk score and gene mutation in the training, external validation cohorts. Please note: The publisher is not responsible for the content or functionality of any supporting information supplied by the authors. Any queries (other than missing content) should be directed to the corresponding author for the article.

Abstract:Clinical and Translational MedicineVolume 10, Issue 2 e104 LETTER TO EDITOROpen Access Analysis of tumor microenvironment-related key mRNAs and construction of a prognosis signature in colon cancer Zhiyuan Zhang, Zhiyuan Zhang orcid.org/0000-0003-0336-7783 Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorGuodong He, Guodong He Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorQingyang Feng, Qingyang Feng Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorPeng Zheng, Peng Zheng Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYang Lv, Yang Lv Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYihao Mao, Yihao Mao Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYuqiu Xu, Yuqiu Xu Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorJianmin Xu, Corresponding Author Jianmin Xu xujmin@aliyun.com Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, China Correspondence Jianmin Xu, Department of General Surgery, Zhongshan Hospital, Fudan University, 180 Fenglin Road, Shanghai 200032, China. Email: xujmin@aliyun.comSearch for more papers by this author Zhiyuan Zhang, Zhiyuan Zhang orcid.org/0000-0003-0336-7783 Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorGuodong He, Guodong He Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorQingyang Feng, Qingyang Feng Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorPeng Zheng, Peng Zheng Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYang Lv, Yang Lv Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYihao Mao, Yihao Mao Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorYuqiu Xu, Yuqiu Xu Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, ChinaSearch for more papers by this authorJianmin Xu, Corresponding Author Jianmin Xu xujmin@aliyun.com Department of General Surgery, Zhongshan Hospital, Fudan University, Shanghai, China Correspondence Jianmin Xu, Department of General Surgery, Zhongshan Hospital, Fudan University, 180 Fenglin Road, Shanghai 200032, China. Email: xujmin@aliyun.comSearch for more papers by this author First published: 20 June 2020 https://doi.org/10.1002/ctm2.104Citations: 3 Zhiyuan Zhang, Guodong He, and Qingyang Feng contributed equally to this work. AboutSectionsPDF ToolsRequest permissionExport citationAdd to favoritesTrack citation ShareShare Give accessShare full text accessShare full-text accessPlease review our Terms and Conditions of Use and check box below to share full-text version of article.I have read and accept the Wiley Online Library Terms and Conditions of UseShareable LinkUse the link below to share a full-text version of this article with your friends and colleagues. Learn more.Copy URL Share a linkShare onFacebookTwitterLinkedInRedditWechat Dear Editor, Tumor microenvironment (TME) serves as a crucial factor during the process of colon cancer, which can provide valuable information for prognosis of colon cancer patients.1 We used bioinformatics analysis flow (Figure S1) to integrated analyze the prognostic TME-related mRNAs and construct a prognostic model to predict the overall survival in colon cancer patients. In total, 2831 mRNAs were identified as prognostic mRNAs and were combined with the results of ESTIMATE analysis for further analysis via WGCNA (weight gene co-expression network analysis).2, 3 The process of WGCNA was shown in Figures S2A, S2B, and 1A(a), eight gene modules were identified (Figure 1A[b]), cells that were correlated with any one of the three scores were screened out (the criterion was correlation value >0.4 and P-value < .05), and the mRNAs involved in these cells were identified as TME-related key mRNAs (the criterion was corresponding |Gene Significance (GS)| > 0.2; |Module Membership (MM)| > 0.7; P-value of GS and MM <.05). Finally, 317 mRNAs were identified as TME-related key mRNAs. GO and KEGG analysis indicated that TME-related key mRNAs mainly enriched in the cell cycle (Figures S3A and S3B). It revealed that cell cycle also plays a vital role in regulating the TME of colon cancer. The PPI network was used to identify the 20 hub mRNAs and to assess the interactions (Figures S2C and 2CD): the cluster that owned the highest score contained 833 interactions but only 43 mRNAs, and the number of interactions was much larger compared to the number of mRNAs involved in the PPI network indicating that the mechanism of TME in colon cancer was complex. FIGURE 1Open in figure viewerPowerPoint WGCNA and construction and validation of the prognostic signature. A, WGCNA based on overall survival-related mRNAs: (a) hierarchical cluster analysis dendrogram to identify co-expression models along with relevant color assignments; (b) correlation analysis between the gene modules and TME-related scores including immune score, stromal score, and ESTIMATE score. Each cell contained a relevant correlation and P-value. The cells were colored by correlation, which decreased in degree from red to blue. Gene Significance (GS) and Module Membership (MM) were computed to assess the relationship between gene expression and ESTIMATE scores (immune score, stromal score, and estimate score) and individual module. B, Construction and validation of prognostic signature. Efficacy of the prognostic signature in the internal training set (a), the internal validation set (b), the external validation set (c), and the entire set (d). The upper panel presented the risk score analysis of the signature, the middle panel presented the Kaplan-Meier survival analysis of the signature, and the cutoff of high-risk and low-risk groups was the median risk score in the internal training set. The lower panel presented the ROC curves of the signature. C, Univariate and multivariate Cox analysis in the entire set. (a) Forest plot depicting the results of univariate Cox analysis. (b) Forest plot depicting the results of multivariate Cox analysis (The risk factors that were identified as statistically significant in univariate Cox regression were further analyzed by multivariate Cox regression). P-value < .05 was identified as statistically significant In the internal training set, least absolute shrinkage and selection operator (LASSO) was performed based on the TME-related key mRNAs.4 We performed cross-validation, and the log(λ) that owned the lowest deviance was selected to screen the mRNAs with the coefficients that were not 0 (Figure S4). In total, 21 mRNAs were identified and further included to construct the prognostic signature to predict overall survival of colon cancer patients via multivariate Cox regression analysis. The signature was as follows: risk score = (exp level of APOL3 * –0.473) + (exp level of FIGNL1 * –0.188) + (exp level of HLX * 0.496) + (expression level of MCM6 * –0.251) + (exp level of RNF114 * –0.645) + (exp level of TIMM13 * –0.692) + (exp level of TNFRSF9 * –0.491). The efficacy of the signature was assessed and validated to be great in the internal training set, the internal validation set, and the external validation (Figure 1B). In the entire set, the area under the curve (AUC) at 1, 3, and 5 years was 0.745, 0.673, and 0.651, respectively. The AUC at 1 year ranged from 0.729 to 0.766 in different sets indicating the great efficiency at 1 year, though decreased to about 0.65 at 3 and 5 year, and the signature still remained good stability and efficiency. Univariate and multivariate analyses proved the risk score derived from the signature to be an independent prognostic factor (Figure 1C). The entire set was classified into different subgroups based on several clinical features; the AUC of different subgroups at 1, 3, and 5 years did not alter much indicating that our signature was stable in different situations (Figure 2). We further checked if our signature and mRNA involved in the signature could have a great correlation with the TME, which was estimated via TME scores and infiltrating immune cells abundance (assessed by ESTIMATE and CIBERSORT algorithm). The results were shown in Figure S5 indicating that our signature and mRNA involved in the signature had a great correlation with TME. FIGURE 2Open in figure viewerPowerPoint Evaluation of the signature in different subgroups classified by different clinical features. The entire set was divided into diverse subgroups according to age (A), gender (B), KRAS status (C), and BRAF status (D); the efficacy was evaluated in diverse groups. The upper panel presented the ROC curves and the lower panel presented the Kaplan-Meier survival analysis of the signature in corresponding groups; the cutoff of high risk and low risk groups was the median risk score in the internal training set In conclusion, we comprehensively analyzed mRNA that might involve in mediating of TME in colon cancer, we identified TME-related key mRNAs, and mining the potential interactions and the functional enrichment of these mRNAs, moreover, we constructed a robust and stable signature in order to predict the overall survival of colon cancer patients. Our analysis could help to strengthen the knowledge of specific mRNAs in the regulation of TME in colon cancer and the signature could further assist to better predict the prognostic outcomes of colon cancer patients. ACKNOWLEDGMENTS We want to thank the GEO database for providing the raw data. The National Natural Science Foundation of China (No. 81702302 and 81472228) granted funding support for our research. CONFLICT OF INTEREST The authors declare no conflict of interest. AVAILABILITY OF DATA AND MATERIALS Data and materials are all available. CONSENT FOR PUBLICATION All authors consent to publish. ETHICS APPROVAL AND CONSENT TO PARTICIPATE This is a study based on data from public datasets, so we want to waive the ethics consent. AUTHOR CONTRIBUTION ZYZ performed and designed the study, and wrote this letter. JMX, GDH, and QYF assisted to improve the design. YL, YHM, YQX, and PZ helped to analyze the data and write the article. Supporting Information Filename Description ctm2104-sup-0001-figureS1.tif217.6 KB Figure S1 Flow chart of the whole research. ctm2104-sup-0002-figureS2.tif560.8 KB Figure S2 WGCNA analysis based on overall survival-related mRNAs and TME scores. A, Identification of network topology of various soft-thresholding powers and appropriate soft-thresholding power to construct a scale-free network. B, Dendrogram of clustering eigengenes with consensus modules. 0.25 was set as the cutoff of the height to merge similar modules which corresponded to a correlation of 0.75. ctm2104-sup-0003-figureS3.tif2.9 MB Figure S3 Functional enrichment analysis and PPI network construction based on TME-related key mRNAs. A, GO analysis of TME-related key mRNAs. B, KEGG analysis of TME-related key mRNAs. C, Top 20 hub mRNAs in the whole PPI network. D, Cluster with the highest score of the whole PPI network. ctm2104-sup-0004-figureS4.tif754.8 KB Figure S4 Results of LASSO Cox regression analysis. A, The correlation between log(λ) and the deviance. B, The correlation between log(λ) and the coefficients of relevant mRNAs. ctm2104-sup-0005-figureS5.tif829.9 KB Figure S5 Correlation analysis between signature, mRNAs involved in the signature and TME-related scores, immune cell infiltrating. Each cell represented the correlation between expression of mRNA or risk score with corresponding TME-related scores or immune cell infiltrating. Cells were colored by correlation, which decreased in degree from red to blue. * and ** represent the p-value of correlation value was < 0.05 and 0.01 respectively. Please note: The publisher is not responsible for the content or functionality of any supporting information supplied by the authors. Any queries (other than missing content) should be directed to the corresponding author for the article. REFERENCES 1Feng Q, Chang W, Mao Y, et al. Tumor-associated macrophages as prognostic and predictive biomarkers for postoperative adjuvant chemotherapy in patients with stage II colon cancer. Clin Cancer Res. 2019; 25(13): 3896- 3907. 2Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008; 9: 559. 3Yoshihara K, Shahmoradgoli M, Martínez E, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013; 4: 2612. 4Tibshirani R. The lasso method for variable selection in the Cox model. Stat Med. 1997; 16(4): 385- 395. Citing Literature Volume10, Issue2June 2020e104 FiguresReferencesRelatedInformation

A Mutational Signature for Colorectal Cancer Prognosis Prediction: Associated with Immune Cell Infiltration

A Mutational Signature That Predicts Prognosis and Benefit of Immune Checkpoint Blockade in Colorectal Cancer

Prognostic Risk Model of Immune-Related Genes in Colorectal Cancer

Identification of a Novel Immune-Related Prognostic Model for Patients with Colorectal Cancer Based on 3 Subtypes.

A New Prognostic and Therapeutic Immune-related Risk Signature of Colorectal Cancer

Predictive Value Of Immunoscore In Patients With Colorectal Cancer

Low Mutation Burden and Differential Tumor-Infiltrating Immune Cells Correlate with Lymph Node Metastasis in Colorectal Cancer

Analysis of Tumor Microenvironment‐related Key Mrnas and Construction of a Prognosis Signature in Colon Cancer

Development of 14‐gene signature prognostic model based on metastasis for colorectal cancer

Mutation patterns in colorectal cancer and their relationship with prognosis

mRNAsi-related metabolic risk score model identifies poor prognosis, immunoevasive contexture, and low chemotherapy response in colorectal cancer patients through machine learning

Molecular subtype identification and prognosis stratification by a metabolism-related gene expression signature in colorectal cancer

A novel risk-scoring model developed based on tertiary lymphoid structure related genes predict the prognosis and immunotherapeutic efficacy in cutaneous melanoma.

An Immune Landscape Based Prognostic Signature Predicts the Immune Status and Immunotherapeutic Responses of Patients with Colorectal Cancer.

Comprehensive Analysis of Genomic and Expression Data Identified Potential Markers for Predicting Prognosis and Immune Response in CRC

A Prognostic Model Based on Nine DNA Methylation-Driven Genes Predicts Overall Survival for Colorectal Cancer

Clinicopathological Features Combined With Immune Infiltration Could Well Distinguish Outcomes in Stage II and Stage III Colorectal Cancer: A Retrospective Study

The prognostic index of m7G-related genes in CRC correlates with immune infiltration

A Novel Prognostic Signatures Based on Metastasis- and Immune-Related Gene Pairs for Colorectal Cancer.

Comprehensive Analysis of Microsatellite-Related Transcriptomic Signature and Identify Its Clinical Value in Colon Cancer

The exploration of mitochondrial-related features helps to reveal the prognosis and immunotherapy methods of colorectal cancer