Random forest classifier and molecular subtype classifier for molecular typing of esophageal squamous carcinoma and construction method of random forest classifier and molecular subtype classifier

Abstract

The invention belongs to the technical field of bioinformatics and computational biology, and relates to a random forest classifier and a molecular subtype classifier for esophageal squamous cell carcinoma molecular typing and a construction method. The number of molecular subtypes is four, and the molecular subtypes are respectively ECMS 1, ECMS 2, ECMS 3 and ECMS 4. The method is characterized in that: ECMS 1: when a metabolic pathway is abnormal and NFE2L2 is activated, an NFE2L2 inhibitor is selected as a drug; eCMS 2: tumor classical signal pathway is up-regulated, and methylation degree is low; eCMS 3 has the advantages of few copy number variation events, low tumor mutation load and high PD-1 expression, and is beneficial to immunosuppressor treatment; eCMS 4: the epithelial-mesenchymal transition pathway is activated, the overall survival rate and disease-free survival rate are low, and epithelial-mesenchymal transition and VEFG inhibitor treatment are benefited. According to the invention, heterogeneity of esophageal squamous carcinoma is analyzed, and each typing system is unified, so that molecular typing is realized, and popularization is easier.

Description

Technical Field

[0001] The present invention belongs to the technical field of bioinformatics and computational biology, and specifically relates to a random forest classifier and a molecular subtype classifier for molecular typing of esophageal squamous cell carcinoma and a construction method thereof. Background Art

[0002] Esophageal cancer is one of the most common cancers in the world, with a poor prognosis and survival rate, which seriously threatens the health of the people. According to histological classification, it is divided into esophageal adenocarcinoma (EAC) or esophageal squamous cell carcinoma (ESCC). In China, ESCC is the main subtype of esophageal cancer, accounting for 90% of cases. The incidence rate varies across China. Because the early symptoms of ESCC are not obvious, ESCC patients in my country are often in the advanced stage when they are first diagnosed, with a low surgical resection rate and poor treatment effect, with a 5-year survival rate of only 10-15%. Compared with other types of tumors, ESCC lacks effective and specific targeted therapeutic drugs. Although traditional surgery, radiotherapy and chemotherapy have made great progress, their efficacy is limited. Local recurrence and distant metastasis are still important reasons for the failure of ESCC treatment. In addition, there is a high degree of heterogeneity within esophageal squamous cell carcinoma and between patients. The diversity of its molecular biology and clinical manifestations means that traditional classification methods based on histopathological characteristics, tumor size, stage, etc. cannot provide intrinsic molecular biological connotations, making it difficult to determine whether patients can benefit from a specific treatment plan. This is an important cause of treatment tolerance, recurrence and metastasis, and has become a bottleneck restricting esophageal cancer treatment research. The treatment benefits of esophageal squamous cell carcinoma patients are limited, resulting in a low survival rate for patients with esophageal squamous cell carcinoma.

[0003] Improving the diagnosis and treatment of esophageal squamous cell carcinoma and prolonging the survival time of patients are of great significance to improving the health level and quality of life of the people, and have great social and economic value. Molecular typing is the most effective strategy to decipher tumor heterogeneity and is also a key strategy to achieve precision treatment. Tumor patients with different molecular subtypes generally have their own specific clinical characteristics, such as prognosis, outcome or drug sensitivity. So far, molecular typing of tumors has been used to guide clinical treatment. For example, in the diagnosis and treatment guidelines for breast cancer, different treatment strategies are recommended for five different molecular subtypes. It is urgent to find molecular typing methods for esophageal squamous cell carcinoma and explain the regulatory mechanisms specific to esophageal squamous cell carcinoma subtypes, and on this basis, design personalized precision treatment plans in a targeted manner. In recent years, studies on molecular typing related to esophageal squamous cell carcinoma have been conducted to study the heterogeneity of ESCC by exploring gene expression, copy number variation or / and mutation spectrum. For example, two tumor analysis subtypes identified based on RNA expression profiles in Asian populations have been reported. Subtype I is associated with immune response, while subtype II is mainly associated with ectoderm development, glycolysis process and cell proliferation. However, these studies are based on single-omics data. Single-omics data ignore the heterogeneity of other omics data. The Cancer and Tumor Genome Atlas (TCGA) project conducted a multi-omics molecular classification study on esophageal squamous cell carcinoma. Three subtypes (C1-3) were identified through comprehensive multi-omics analysis of 90 ESCC patients (including only one Asian population, Vietnamese). C1 is associated with changes in the NRF2 pathway. C2 is associated with high-frequency somatic mutations in NOTCH1 or ZNF750. C3 is characterized by PI3K pathway activation. However, this study did not include the Chinese population, and C1 tends to have more Vietnamese patients, which indicates that the classification system is biased and the classification subtype is not suitable for the Chinese population. In addition, there are currently two major aspects that hinder the application of subtyping systems in clinical applications. First, different algorithms and patient cohorts lead to different results in different classification studies. There is an urgent need to establish a consensus molecular subtyping system in the existing subtyping system. Summary of the invention

[0004] The present invention provides a random forest classifier for molecular typing of esophageal squamous cell carcinoma. The classifier includes 314 characteristic genes, which are represented by Ensembl database numbers: ENSG00000177272, ENSG00000122224, ENSG00000078589, ENSG00000132185, ENSG00000105369, ENSG00000177455, ENSG00000147138, ENSG0000000177272, ENSG00000122224, ENSG00000078589, ENSG00000132185, ENSG00000105369, 0012124, ENSG00000163564, ENSG00000117091, ENSG00000156738, ENSG00000117215, ENSG00000161405, ENSG00000169442, ENSG00000113088, ENSG 00000137078, ENSG00000254838, ENSG00000175463, ENSG00000162739, ENSG00000164691, ENSG00000185905, ENSG00000182162, ENSG00000048462, ENSG00000100079, ENSG00000163534, ENSG00000166501, ENSG00000128218, ENSG00000132465, ENSG00000005844, ENSG00000160654, ENSG00000079 263, ENSG00000143297, ENSG00000113263, ENSG00000173200, ENSG00000168918, ENSG00000124406, ENSG00000136573, ENSG00000183918, ENSG0000 0115085, ENSG00000122188, ENSG00000122122, ENSG00000081237, ENSG00000117090, ENSG00000110777, ENSG00000137265, ENSG00000026751, ENSG 00000115165, ENSG00000167984, ENSG00000110448, ENSG00000035403, ENSG00000162777, ENSG00000167077, ENSG00000122550, ENSG00000185811,ENSG00000100351,ENSG00000106948,ENSG00000167286,ENSG00000134242,ENSG00000143851,ENSG00000110934,ENSG00000073792,ENSG00000153283,ENSG00000247774,ENSG00000198821,ENSG00000117322,ENSG00000010671,ENSG00000182866,ENSG00000196684,ENSG00000128815,ENSG00000269404,ENSG00000167208,ENSG00000172673,ENSG00000076662,ENSG00000124256,ENSG00000170476,ENSG00000197943,ENSG00000181847,ENSG00000009790,ENSG00000124772,ENSG00000099958,ENSG00000072858,ENSG00000153885,ENSG00000172578,ENSG00000057019,ENSG00000182963,ENSG00000108405,ENSG00000198851,ENSG00000145649,ENSG00000172215,ENSG00000164283,ENSG00000179583,ENSG00000185862,ENSG00000187498,ENSG00000013725,ENSG00000176105,ENSG00000095585,ENSG00000105122,ENSG00000164938,ENSG00000047457,ENSG00000042980,ENSG00000073849,ENSG00000101082,ENSG00000173905,ENSG00000146676,ENSG00000153563,ENSG00000105374,ENSG00000106263,ENSG00000187678,ENSG00000163577,ENSG00000168421,ENSG00000143067,ENSG00000166428,ENSG00000111711,ENSG00000205744,ENSG00000171388,ENSG00000173198,ENSG00000101333,ENSG00000147689,ENSG00000117308,ENSG00000168398,ENSG00000163220,ENSG00000173156,ENSG00000162946,ENSG00000165272,ENSG00000070018,ENSG00000213160,ENSG00000114790,ENSG00000153048,ENSG00000105355,ENSG00000162366,ENSG00000175591,ENSG00000143546,ENSG00000242808,ENSG00000005884,ENSG00000162409,ENSG00000157379,ENSG00000060558,ENSG00000111424,ENSG00000148344,ENSG00000142669,ENSG00000189334,ENSG00000178038,ENSG00000153294,ENSG00000163596,ENSG00000147697,ENSG00000110047,ENSG00000030582,ENSG00000146267,ENSG00000144749,ENSG00000186847,ENSG00000073910,ENSG00000115008,ENSG00000164520,ENSG00000092820,ENSG00000158315,ENSG00000142627,ENSG00000159348,ENSG00000115009,ENSG00000062038,ENSG00000171236,ENSG00000158552,ENSG00000124145,ENSG00000124107,ENSG00000134070,ENSG00000088726,ENSG00000177191,ENSG00000176788,ENSG00000187583,ENSG00000146072,ENSG00000128596,ENSG00000011009,ENSG00000116871,ENSG00000099992,ENSG00000182612,ENSG00000100097,ENSG00000150093,ENSG00000133816,ENSG00000078098,ENSG00000206538,ENSG00000166922,ENSG00000145431,ENSG00000196923,ENSG00000115414,ENSG00000167123,ENSG00000072110,ENSG00000165617,ENSG00000122786,ENSG00000072682,ENSG00000161638,ENSG00000122884,ENSG00000120708,ENSG00000087116,ENSG00000164294,ENSG00000070404,ENSG00000128805,ENSG00000122870,ENSG00000163430,ENSG00000142173,ENSG00000120149,ENSG00000198910,ENSG00000113721,ENSG00000133110,ENSG00000080573,ENSG00000157227,ENSG00000162745,ENSG00000105472,ENSG00000117385,ENSG00000113140,ENSG00000164574,ENSG00000151135,ENSG00000142552,ENSG00000143502,ENSG00000164111,ENSG00000149257,ENSG00000166250,ENSG00000177570,ENSG00000143226,ENSG00000184012,ENSG00000107819,ENSG00000172638,ENSG00000120820,ENSG00000083444,ENSG00000142156,ENSG00000108602,ENSG00000163297,ENSG00000149380,ENSG00000130635,ENSG00000152049,ENSG00000112769,ENSG00000100196,ENSG00000115363,ENSG00000142192,ENSG00000101825,ENSG00000110723,ENSG00000106333,ENSG00000079150,ENSG00000163520,ENSG00000111817,ENSG00000128833,ENSG00000166130,ENSG00000105974,ENSG00000026025,ENSG00000138316,ENSG00000151617,ENSG00000140682,ENSG00000221968,ENSG00000131459,ENSG00000140416,ENSG00000164465,ENSG00000218336,ENSG00000168542,ENSG00000154096,ENSG00000134013,ENSG00000182667,ENSG00000142227,ENSG00000145817,ENSG00000151388,ENSG00000127249,ENSG00000165915,ENSG00000196968,ENSG00000164692,ENSG00000087303,ENSG00000184828,ENSG00000143387,ENSG00000177119,ENSG00000091136,ENSG00000134668,ENSG00000128595,ENSG00000183876,ENSG00000124225,ENSG00000196344,ENSG00000158186,ENSG00000137573,ENSG00000258947,ENSG00000088836,ENSG00000131236,ENSG00000171401,ENSG00000163359,ENSG00000129009,ENSG00000123130,ENSG00000174705,ENSG00000153561,ENSG00000150630,ENSG00000170558,ENSG00000163209,ENSG00000106624,ENSG00000182492,ENSG00000168487,ENSG00000145113,ENSG00000126458,ENSG00000135127,ENSG00000149557,ENSG00000139211,ENSG00000016602,ENSG00000129514,ENSG00000035862,ENSG00000113083,ENSG00000114353,ENSG00000140279,ENSG00000167397,ENSG00000162576,ENSG00000116132,ENSG00000106397, ENSG00000039560, ENSG00000006327, ENSG00000214711, ENSG00000136694, ENSG00000185873, ENSG00000189377, ENSG00000183844. ,

[0005] The present invention also provides a molecular subtype classifier of esophageal squamous cell carcinoma constructed by using the random forest classifier, and the obtained consensus molecular subtypes are 4, namely: ECMS 1, ECMS 2, ECMS 3, and ECMS 4; the characteristics of the four molecular subtypes are: ECMS 1: abnormal metabolic pathway, NFE2L2 activation, and NFE2L2 inhibitors are suitable for drug selection; ECMS 2: upregulation of classical tumor signaling pathways and low methylation level; ECMS 3: few copy number variation events, low tumor mutation load, high expression of PD-1, and benefit from immunosuppressant treatment; ECMS 4: activation of the epithelial-mesenchymal transition pathway, low overall survival rate and disease-free survival rate, and benefit from epithelial-mesenchymal transition and VEFG inhibitor treatment.

[0006] The present invention also provides a method for constructing the esophageal squamous cell carcinoma molecular subtype classifier, comprising the following steps: Obtain the whole genome DNA, transcriptome, and small RNA of esophageal squamous cell carcinoma pathology, perform whole genome sequencing, and preprocess the high-throughput sequencing data to obtain the esophageal squamous cell carcinoma pathology omics data; Based on the obtained omics data, namely copy number, transcriptome, small RNA and methylation data, the omics data were integrated for subtype recognition and identification of molecular subtypes of esophageal squamous cell carcinoma; Collect the esophageal squamous cell carcinoma classification labels from the public dataset, perform consensus molecular classification analysis, and obtain the formula molecular classification subtypes; Hypergeometric test, Fisher test, and t test were applied to identify subtype-specific molecular features; Based on differential gene analysis and gene set enrichment analysis (GSEA), the relationship between each subtype and immunity, metabolism, classical tumor pathways, and epithelial-mesenchymal transition was identified; Immune microenvironment analysis was performed using the R packages deconvolution and ESTIMATE, and survival analysis was performed using Kaplan-Meier; Based on the analysis results, the molecular subtypes and characteristics of ESCC were identified, and an ESCC molecular subtype classifier was constructed.

[0007] Furthermore, the high-throughput sequencing data were preprocessed, specifically: for whole-genome DNA sequencing, the original sequencing data were quality-controlled and filtered using FASTP software to filter low-quality reads, remove sequencing adapter sequences, and remove shorter reads; the sequencing reads were aligned with the human reference genome hg38 using BWA software to obtain SAM files; the SAM files were converted to BAM files using samtools, and the BAM files were sorted according to coordinate positions using Sambamba; the MarkDuplicates mark in the Picard toolkit was used to mark duplicate reads; and the cnv_facets tool was used to detect somatic copy number variations (SCNVs) in tumor samples and paired normal samples.

[0008] For transcriptome sequencing, FASTP software was used to perform quality control filtering on the raw sequencing data, filter low-quality reads, remove sequencing adapter sequences, remove shorter reads, and correct some bases based on the overlapping regions of paired reads; the reads after quality control were aligned to the human reference genome hg38 sequence by the STAR tool and the data were BAM files; the gene expression abundance was estimated by the RSEM tool based on the 98th version of the Ensembl reference annotation and using the unique aligned reads, and the gene expression values ​​were standardized to TPM; The read lengths of whole-genome bisulfite sequencing were first quality-controlled and filtered using FASTP software, and then aligned using the Bismark tool. The resulting BAM file was stripped of duplicate read lengths, and the methylation signals were extracted using the Bismark methylation extractor. For small RNA sequencing data, the adapter sequences were first removed using FASTP software, and the short sequences were compared with the small RNA database miRBase and the human reference genome hg38 sequence using miRDeep2. The expression abundance of each small RNA was obtained using the quantitative module of miRDeep2.

[0009] The copy number data: the log ratio values ​​of the data were extracted from the cnv_facets results, and then the redundant copy number regions were removed based on the CNregions function in the R package iClusterPlus; The methylation data were obtained by merging the methylation levels of all sites of all samples into a unified data matrix using the R package methylKit, and calculating the methylation levels of 2 kb upstream and 0.5 kb downstream of the gene transcription start site. The specific calculation method was to divide the number of methylated sequences detected at all CpG sites in the promoter region of the gene by the total number of sequences. If 20% of the samples were missing values ​​or had low coverage in the promoter region, the data in the promoter region were not considered. The gene expression values ​​were first normalized to the number of transcripts per million reads (TPM), and the gene expression levels of all samples were imported using the R package tximport to obtain a data matrix on gene expression; The expression of small RNA was standardized to the number of small RNAs per 1 million reads CPM, and the small RNA expression values ​​of all samples were read to obtain the small RNA data matrix; The data format of each omics is exemplified as follows: data of a specific omics, columns are for specific samples, a total of n samples, and a total of four omics data matrices; If the matrix is ​​a copy number data matrix, the feature is the copy number variation region, the data is the log ratio value obtained by comparing the tumor sample with the normal sample in the region, and the number of features is the number of regions after redundancy removal; If the matrix is ​​a gene expression matrix, the features are genes, the data are the expression levels of the genes, specifically log2(TPM+1), and the number of features is the number of genes; If the matrix is ​​a methylation data matrix, the feature is the promoter region, the data is the methylation level of the promoter region, ranging from 0 to 1, and the number of features is the number of available promoters; If the matrix is ​​a small RNA expression matrix, the feature is the expression level of the small RNA, specifically log2(CPM+1), and the number of features is the total number of small RNAs; Further filtering was performed based on the median absolute deviation to select the top 1000 genes with large variations, the top 500 small RNAs, the top 1000 promoters, and all copy number regions after redundancy removal; the data of each omics were Z-value standardized, and the Euclidean distance between two samples was calculated at the level of single omics data; the R package SNF was used to construct a similarity network at the level of single omics, and the network was fused under the conditions of hyperparameters alpha of 0.75, K of 30, and T=20; Based on the analysis results, the method for identifying the molecular subtypes and characteristics of esophageal squamous cell carcinoma is as follows: the molecular typing label of each sample is obtained by using the spectral clustering algorithm; to evaluate the contribution of each omics to the subtype, the edges with the top 5% similarity weight are selected; if the difference in variable weights of multiple data types is less than 10%, it is considered that the edge between the two samples in the network is contributed by multiple omics data; if the edge with the highest weight is greater than the edge with the second highest weight, it is considered that the evidence of the edge between the two samples in the network comes from a specific single omics; through the similarity network fusion method based on multi-omics data, four multi-omics molecular subtypes were obtained: MESCC1, MESCC 2, MESCC 3, and MESCC 4.

[0010] The typing labels of seven groups of published typing studies (PubMed ID: 29127303, 32398863, 32580137, 31048097, 28052061, 28052061, 36584672) were further collected to construct the consensus molecular typing input data. 85% of the samples were randomly selected, and the hypergeometric test was performed on the correlation between the typing systems. The p value was less than 0.05, which was considered to be correlated between subtypes. A connected 26-node association network (each node is a molecular subtype of a specific typing study) was constructed, and the Markov clustering algorithm (MCL) was used to identify clusters in the network. The algorithm was repeated 100 times, and the frequency of molecular subtypes in different typing systems being divided into unified clusters was counted to construct a consensus frequency matrix (26*26). MCL was further used to identify the consensus molecular typing ECMS 1, ECMS 2, ECMS 3, and ECMS 4. A hypergeometric test was performed to identify core ECMS samples (P value less than 0.2), and a classifier was trained to predict the consensus molecular subtypes of samples with P value greater than 0.2 to obtain the consensus molecular subtypes of all samples.

[0011] The present invention uses high-throughput sequencing to sequence the whole genome DNA, transcriptome, small RNA, and whole genome methylation of tumor tissues of 152 patients with esophageal squamous cell carcinoma. The high-throughput data is preprocessed on a large scale using quality control software, sequence alignment and processing software, mutation detection software, transcriptome and small RNA quantification software to obtain the omics data of patients with esophageal squamous cell carcinoma. Based on the similarity network fusion method, four omics data are fused for subtype identification, and four multi-omics molecular subtypes of esophageal squamous cell carcinoma (MESCC1-4) are identified. Further, four consensus molecular subtypes (ECMS1-4) are identified through formula molecular typing research. Hypergeometric test, Fisher test, and t-test are used to identify subtype-specific molecular features. Based on differential gene analysis and gene set enrichment analysis (GSEA), each subtype is identified to be related to immunity, metabolism, classical tumor pathways, and epithelial-mesenchymal transition. Subtype ECMS1 is enriched with NFE2L2 amplification events. NFE2L2 in MESCC2 has a higher expression level than the other three subtypes, indicating that NFE2L2 is activated in MESCC2. ECMS3 has fewer copy number events, and the tumor mutation burden (TMB) of ECMS3 is relatively low. There are fewer KMT2D mutations in patients with the ECMS3 subtype. The ECMS3 subtype is characterized by higher immune activity, including myeloid-derived suppressor cells (MDSCs), natural killer cells (NKC), inflammatory helper T cells (Th17 cells), and complement system. ECMS1 is characterized by dysregulation of various metabolic pathways, including cytochrome, glutathione, sucrose, glycoside, and uronic acid metabolic pathways. GSEA revealed that ECMS2 is characterized by upregulation of classical cell cycle-related pathways (WNT, MYC signaling pathways, and KEGG cell cycle pathways). The ECMS4 subtype exhibits mesenchymal characteristics characterized by TGFβ matrix remodeling, VEGF pathway, and epithelial-mesenchymal transition (EMT), and is called the EMT subtype. Compared with the other three subtypes, genes related to epithelial-mesenchymal were significantly overexpressed in the EMT subtype. In addition, the GSEA results based on proteomics data also showed similar GSEA molecular features based on transcriptional sequencing data, indicating that the four consensus molecular subtypes had the same molecular patterns at the protein level.

[0012] Immune microenvironment analysis was performed using the R packages deconvolution and ESTIMATE. ESTIMATE results showed that the immune score of ECMS3 was significantly higher than that of the other three subtypes, while the stromal score of ECM4 was significantly higher than that of ECMS1 and ECMS2. In ECMS3 and ECMS4, the tumor purity was lower, indicating a potential interaction between the tumor and the tumor microenvironment during progression. TIMER immune infiltration analysis also showed that MESCC1 had a higher immune infiltration content, such as B cells, T cell CD4+, and T cell CD8+ components in ECMS3. Tumor-associated macrophages can promote EMT in ESCC, and ECMS4 and ECMS3 have higher macrophage infiltration compared with the classic subtype ECMS2.

[0013] Survival analysis was performed using Kaplan-Meier, and the results of survival analysis of the four subtypes showed that patients with the MESCC4 subtype showed the worst overall survival (OS) and disease-free survival (DFS). The present invention constructs a classifier for molecular typing based on the characteristic genes and gene expression levels of the four subtypes at the transcriptome level. Using the classifier of the present invention, the clinical relevance of the four subtypes was verified. Based on the data analysis of the clinical trial, the ECMS3 subtype was significantly associated with the treatment benefit group in immune checkpoint inhibitor treatment, and the expression level of PDCD1 (PD-1) was significantly higher than that of the other three groups, indicating that patients with ECMS3 may benefit from immune checkpoint inhibitor treatment. In addition, the association analysis between drug sensitivity and gene expression showed that the data came from the cell line data of CTRP (823 cell lines). Negative correlation indicates possible sensitivity to drugs. VER-155008 showed a significant negative correlation with NFE2L2 expression, indicating that Necrosulfonamide may be a potential therapeutic strategy for targeting NFE2L2 in ECMS1 patients. In addition, cerulenin and PF-3758309 could potentially be used to target MYC and VEGFC in ECMS2 and ECMS4.

[0014] The present invention provides a multi-omics molecular subtype classification method and system for esophageal squamous cell carcinoma, and provides the molecular characteristics, clinical characteristics and classifiers of each subtype to analyze the heterogeneity of esophageal squamous cell carcinoma, as well as a classification method and system for the consensus molecular classification of esophageal squamous cell carcinoma, unifying various classification systems, thereby making molecular classification easier to promote. BRIEF DESCRIPTION OF THE DRAWINGS

[0015] Figure 1 To support the source and distribution of the evidence of genomic similarity among the various subtypes of the MESCC classification system; Figure 2 is the consensus frequency matrix, i.e., the frequency of molecular subtypes being divided into the same cluster; Figure 3Esophageal squamous cell carcinoma molecular subtype network based on 8 classification systems; Figure 4 Molecular characteristics of consensus molecular subtypes; In the figure: a is the mutation and copy number characteristic map of subtypes, b is the statistical distribution of copy number variation between subtypes, c is the distribution of tumor mutation burden between subtypes, d is the KMT2D mutation distribution, e is the BST2 methylation degree distribution, f is the NFE2L2 expression distribution, g is the enrichment score distribution of related pathways; Figure 5 Mutation frequency of genes in different datasets; Figure 6 Heat map showing the distribution differences of epithelial-mesenchymal transition (EMT), VEGF, and cell cycle-related genes; In the figure: a is the heat map showing the distribution differences of EMT-related genes; b is the heat map showing the distribution differences of VEGF-related genes; c is the heat map showing the distribution differences of cell cycle-related genes; Figure 7 Distribution of immune microenvironment in each consensus molecular subtype; In the figure: a-i are immune score, stromal score, tumor purity, B cells, CD4+ T cells, CD8+ T cells, neutrophils, macrophages, and myeloid dendritic cells respectively; Figure 8 Overall survival (a) and disease-free survival (b); Fig. 9 Overall survival of consensus molecular subtypes in dataset GSE53625 and disease-free survival period in dataset TCGA; In the figure: a is the overall survival of consensus molecular subtypes in dataset GSE53625, b is the disease-free survival period of consensus molecular subtypes in dataset TCGA; Fig.10 Distribution of chemotherapeutic drug sensitivity in each consensus molecular subtype; In the figure: a-d are Cisplatin, Vinorelbine, Paclitaxel, and Docetaxel respectively; Fig.11 Correlation between cell line drug sensitivity and gene expression; In the figure: a is the correlation distribution of related drugs and small molecule sensitivities with NFE2L2 expression level, b is the distribution of VER-155008 sensitivity with NFE2L2 expression level, c is the correlation distribution of related drugs and small molecule sensitivities with MYC expression level, d is the distribution of cerulenin sensitivity with MYC expression level, e is the correlation distribution of related drugs and small molecule sensitivities with VEGCF expression level, f is the distribution of PF-3758309 sensitivity with VEGFC expression level; Fig.12 Clinical treatment response in each consensus molecular subtype; In the figure: a is the GSE45670 cohort; b is the internal cohort; pCR: pathological complete remission; <pCR: less than pCR; PR: partial response; SD: disease stable; PD: disease progression. Specific implementation methods

[0016] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below. Obviously, the described embodiments are part of the embodiments of the present invention, rather than all the embodiments; based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0017] Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention belongs, and the disclosure and materials cited therein are hereby incorporated by reference.

[0018] Technical equivalents to the specific embodiments described that are apparent to those skilled in the art using no more than routine experimentation are intended to be encompassed by this application.

[0019] The experimental methods in the following examples are all conventional methods unless otherwise specified. The instruments and equipment used in the following examples are all conventional laboratory instruments and equipment unless otherwise specified; the experimental materials used in the following examples are all purchased from conventional stores unless otherwise specified.

[0020] 1. Data preprocessing: For whole genome DNA sequencing, FASTP software was used to perform quality control filtering on the raw sequencing data, filter low-quality reads, remove sequencing adapter sequences, and remove shorter reads. Using BWA software, the sequencing reads were aligned with the human reference genome hg38 to obtain SAM files. SAM files were converted to BAM files using samtools, and BAM files were sorted according to coordinate positions using Sambamba. Repeated reads were marked with MarkDuplicates in the Picard toolkit. Base quality was corrected using the GATK tool and the known human variation database to minimize the systematic errors caused by the sequencer. Mutect2's single sample mode was used to identify variants in normal samples and form PON files, and then Mutect2's normal-tumor tissue pairing detection mutation mode was used to detect somatic variants (point mutations and insertion and deletion variants) in each sample. Variants were further filtered based on sequence context and sample contamination. ANNOVAR software was used to annotate somatic mutations, and relevant information such as the region where the mutation was located, the gene, the frequency in the public database, and the effect on the amino acid were obtained. Mutations located in non-gene regions or with ExAC_EAS frequency > 0.0.5 or synonymous mutations were filtered in subsequent analysis. The cnv_facets tool was used to detect somatic copy number variations (SCNVs) in tumor samples and paired normal samples. The R package “mafTools” was used to analyze mutations and SCNVs. We also used GISTIC to identify significantly amplified or deleted genomic regions.

[0021] For transcriptome sequencing, FASTP software was used to perform quality control filtering on the raw sequencing data, filter low-quality reads, remove sequencing adapter sequences, remove shorter reads, and correct some bases based on the paired read overlap regions. The reads after quality control were aligned to the human reference genome 38 sequences by the STAR tool and the data was BAM file. Gene expression abundance was estimated by the RSEM tool based on the Ensembl (98th version) reference annotation and using the unique aligned reads, and the gene expression values ​​were standardized to TPM.

[0022] The read lengths of whole-genome bisulfite sequencing were first quality-controlled and filtered using FASTP software, and then aligned using the Bismark tool. The resulting BAM file was stripped of duplicate read lengths, and the methylation signals were extracted using the Bismark methylation extractor.

[0023] For small RNA sequencing data, the adapter sequences were first removed using FASTP software, and the short sequences were compared with the small RNA database miRBase and the human reference genome hg38 sequence using miRDeep2. The expression abundance of each small RNA was obtained using the quantitative module of miRDeep2.

[0024] 2. Identification of molecular subtypes based on a method of fusion of multi-omics data and similarity networks: The present invention identifies molecular subtypes based on four types of omics data, namely copy number, transcriptome, small RNA and methylation data.

[0025] The log ratio values ​​of the copy number data were extracted from the cnv_facets results, and then the redundant copy number regions were removed based on the CNregions function in the iClusterPlus R package.

[0026] For methylation data, the R package methylKit was used to merge the methylation levels of all sites in all samples into a unified data matrix, and the methylation levels in the gene promoter region (2kb upstream and 0.5kb downstream of the gene transcription start site) were calculated. The specific calculation method was the number of methylated sequences measured at all CpG sites in the gene promoter region of the sample divided by the total number of sequences. If 20% of the samples were missing values ​​or had low coverage in the promoter region, the data in the promoter region were not considered.

[0027] The gene expression data was first standardized to the number of transcripts per million reads (TPM), and the gene expression levels of all samples were imported using the R package tximport to obtain a data matrix for gene expression. The expression of small RNA was standardized to the number of small RNAs per million reads (CPM), and the small RNA expression values ​​of all samples were read to obtain a small RNA data matrix. Before data processing, TPM and CPM were log2 processed. The data format of each omics behavior is the data of a specific omics, listed as a specific sample, a total of n samples, and a total of four omics data matrices. If it is a copy number data matrix, the feature is the copy number variation region, the data is the log ratio value obtained by comparing the tumor sample with the normal sample in the region, and the number of features is the number of regions after de-redundancy; if it is a gene expression matrix, the feature is the gene, the data is the expression level of the gene, specifically log2(TPM+1), and the number of features is the number of genes; if it is a methylation data matrix, the feature is the promoter region, the data is the methylation level of the promoter region, ranging from 0 to 1, and the number of features is the number of available promoters; if it is a small RNA expression matrix, the feature is the expression level of small RNA, specifically log2(CPM+1), and the number of features is the total number of small RNAs.

[0028] To further screen biologically significant features, further filtering was performed based on the median absolute deviation, and the top 1000 genes with large variations, the top 500 small RNAs, the top 1000 promoters, and all copy number regions after redundancy were selected. The data of each omics were Z-value standardized, and the Euclidean distance between two samples was calculated at the level of single omics data. The R package SNF was used to construct a similarity network at the level of single omics, and the network was fused under the conditions of hyperparameters alpha of 0.75, K of 30, and T=20. The molecular typing label of each sample was obtained using the spectral clustering algorithm. In order to evaluate the contribution of each omics to the subtype, the edges with the top 5% similarity weight were selected. If the difference in variable weights of multiple data types was less than 10%, the edge between the two samples in the network was considered to be contributed by multiple omics data. If the edge with the highest weight was greater than the edge with the second highest weight, the evidence of the edge between the two samples in the network was considered to come from a specific single omics. Four molecular subtypes (MESCC1-4) were obtained by the similarity network fusion method based on multi-omics data.

[0029] The typing labels of seven groups of published typing studies (PubMed ID: 29127303, 32398863, 32580137, 31048097, 28052061, 28052061, 36584672) were further collected to construct the consensus molecular typing input data. 85% of the samples were randomly selected, and the hypergeometric test was performed on the correlation between the typing systems. The p value was less than 0.05, which was considered to be correlated between subtypes. A connected 26-node association network (each node is a molecular subtype of a specific typing study) was constructed, and the Markov clustering algorithm (MCL) was used to identify clusters in the network. The algorithm was repeated 100 times, and the frequency of molecular subtypes in different typing systems being divided into unified clusters was counted to construct a consensus frequency matrix (26*26). MCL was further used to identify the consensus molecular typing ECMS 1, ECMS 2, ECMS 3, and ECMS 4. A hypergeometric test was performed to identify core ECMS samples (P value less than 0.2), and a classifier was trained to predict the consensus molecular subtypes of samples with P value greater than 0.2 to obtain the consensus molecular subtypes of all samples.

[0030] 3. Characteristic analysis of each molecular subtype: The R package LIMMA was used to perform differential gene analysis on the specific subtype and the other three subtypes, and the differentially expressed genes of each molecular subtype compared with the other three subtypes were obtained. Then, the R package HTSanalyzeR2 was used for gene set enrichment analysis, and the bootstrap was set to 1000 times. The Kaplan–Meier analysis of survival time was performed using the R package survival. The TIMER method was used to analyze the tumor microenvironment of each subtype using the R package deconvolution, including B cells, CD4+, CD8+, neutrophils, macrophages, and dendritic cells. In addition, the immune score and matrix score of the subtype were evaluated using the ESTIMATE tool. The subtype-specific molecular characteristics were identified using the hypergeometric test, Fisher test, and t-test.

[0031] 4. Construct subtype classifier: Through differential analysis, select the top 50 genes with fold differential upregulation for each subtype to represent the characteristic genes of the subtype. The characteristic genes should exist in multiple data sets at the same time, so it is necessary to perform intersection analysis with genes in other data sets as the final characteristic genes.

[0032] The expression values ​​of the characteristic genes were extracted and z-score standardized. Then, based on the random forest algorithm and the prediction accuracy of the random forest model, the samples with correct predictions were selected as core samples, and the random forest classifier was constructed again. The number of random forest classifier trees was set to 500. This classifier can be used to predict the molecular subtypes of esophageal squamous cell carcinoma in other transcriptome datasets.

[0033] The classifier included 314 characteristic genes (indicated by Ensembl number): ENSG00000177272, ENSG00000122224, ENSG00000078589, ENSG00000132185, ENSG00000105369, ENSG00000177455, ENSG00000147138, ENSG000000017312, ENSG00000172724, ENSG00000107317, ENSG00000121895, ENSG00000012124, ENSG00000163564, 117091, ENSG00000156738, ENSG00000117215, ENSG00000161405, ENSG00000169442, ENSG00000113088, ENSG00000137078, ENSG00000254838, ENSG0 0000175462 NSG00000166501, ENSG00000128218, ENSG00000132465, ENSG00000005844, ENSG00000160654, ENSG00000079263, ENSG00000143297, ENSG00000113 263, ENSG00000173200, ENSG00000168918, ENSG00000124406, ENSG00000136573, ENSG00000183918, ENSG00000115085, ENSG00000122188, ENSG0000 0122122,ENSG00000081237,ENSG00000117090,ENSG00000110777,ENSG00000137265,ENSG00000026751,ENSG00000115165,ENSG00000167984,ENSG 00000110448, ENSG00000035403, ENSG00000162777, ENSG00000167077, ENSG00000122550, ENSG00000185811, ENSG00000100351, ENSG00000106948,ENSG00000167286,ENSG00000134242,ENSG00000143851,ENSG00000110934,ENSG00000073792,ENSG00000153283,ENSG00000247774,ENSG00000198821,ENSG00000117322,ENSG00000010671,ENSG00000182866,ENSG00000196684,ENSG00000128815,ENSG00000269404,ENSG00000167208,ENSG00000172673,ENSG00000076662,ENSG00000124256,ENSG00000170476,ENSG00000197943,ENSG00000181847,ENSG00000009790,ENSG00000124772,ENSG00000099958,ENSG00000072858,ENSG00000153885,ENSG00000172578,ENSG00000057019,ENSG00000182963,ENSG00000108405,ENSG00000198851,ENSG00000145649,ENSG00000172215,ENSG00000164283,ENSG00000179583,ENSG00000185862,ENSG00000187498,ENSG00000013725,ENSG00000176105,ENSG00000095585,ENSG00000105122,ENSG00000164938,ENSG00000047457,ENSG00000042980,ENSG00000073849,ENSG00000101082,ENSG00000173905,ENSG00000146676,ENSG00000153563,ENSG00000105374,ENSG00000106263,ENSG00000187678,ENSG00000163577,ENSG00000168421,ENSG00000143067,ENSG00000166428,ENSG00000111711,ENSG00000205744,ENSG00000171388,ENSG00000173198,ENSG00000101333,ENSG00000147689,ENSG00000117308,ENSG00000168398,ENSG00000163220,ENSG00000173156,ENSG00000162946,ENSG00000165272,ENSG00000070018,ENSG00000213160,ENSG00000114790,ENSG00000153048,ENSG00000105355,ENSG00000162366,ENSG00000175591,ENSG00000143546,ENSG00000242808,ENSG00000005884,ENSG00000162409,ENSG00000157379,ENSG00000060558,ENSG00000111424,ENSG00000148344,ENSG00000142669,ENSG00000189334,ENSG00000178038,ENSG00000153294,ENSG00000163596,ENSG00000147697,ENSG00000110047,ENSG00000030582,ENSG00000146267,ENSG00000144749,ENSG00000186847,ENSG00000073910,ENSG00000115008,ENSG00000164520,ENSG00000092820,ENSG00000158315,ENSG00000142627,ENSG00000159348,ENSG00000115009,ENSG00000062038,ENSG00000171236,ENSG00000158552,ENSG00000124145,ENSG00000124107,ENSG00000134070,ENSG00000088726,ENSG00000177191,ENSG00000176788,ENSG00000187583,ENSG00000146072,ENSG00000128596,ENSG00000011009,ENSG00000116871,ENSG00000099992,ENSG00000182612,ENSG00000100097,ENSG00000150093,ENSG00000133816,ENSG00000078098,ENSG00000206538,ENSG00000166922,ENSG00000145431,ENSG00000196923,ENSG00000115414,ENSG00000167123,ENSG00000072110,ENSG00000165617,ENSG00000122786,ENSG00000072682,ENSG00000161638,ENSG00000122884,ENSG00000120708,ENSG00000087116,ENSG00000164294,ENSG00000070404,ENSG00000128805,ENSG00000122870,ENSG00000163430,ENSG00000142173,ENSG00000120149,ENSG00000198910,ENSG00000113721,ENSG00000133110,ENSG00000080573,ENSG00000157227,ENSG00000162745,ENSG00000105472,ENSG00000117385,ENSG00000113140,ENSG00000164574,ENSG00000151135,ENSG00000142552,ENSG00000143502,ENSG00000164111,ENSG00000149257,ENSG00000166250,ENSG00000177570,ENSG00000143226,ENSG00000184012,ENSG00000107819,ENSG00000172638,ENSG00000120820,ENSG00000083444,ENSG00000142156,ENSG00000108602,ENSG00000163297,ENSG00000149380,ENSG00000130635,ENSG00000152049,ENSG00000112769,ENSG00000100196,ENSG00000115363,ENSG00000142192,ENSG00000101825,ENSG00000110723,ENSG00000106333,ENSG00000079150,ENSG00000163520,ENSG00000111817,ENSG00000128833,ENSG00000166130,ENSG00000105974,ENSG00000026025,ENSG00000138316,ENSG00000151617,ENSG00000140682,ENSG00000221968,ENSG00000131459,ENSG00000140416,ENSG00000164465,ENSG00000218336,ENSG00000168542,ENSG00000154096,ENSG00000134013,ENSG00000182667,ENSG00000142227,ENSG00000145817,ENSG00000151388,ENSG00000127249,ENSG00000165915,ENSG00000196968,ENSG00000164692,ENSG00000087303,ENSG00000184828,ENSG00000143387,ENSG00000177119,ENSG00000091136,ENSG00000134668,ENSG00000128595,ENSG00000183876,ENSG00000124225,ENSG00000196344,ENSG00000158186,ENSG00000137573,ENSG00000258947,ENSG00000088836,ENSG00000131236,ENSG00000171401,ENSG00000163359,ENSG00000129009,ENSG00000123130,ENSG00000174705,ENSG00000153561,ENSG00000150630,ENSG00000170558,ENSG00000163209,ENSG00000106624,ENSG00000182492,ENSG00000168487,ENSG00000145113,ENSG00000126458,ENSG00000135127,ENSG00000149557,ENSG00000139211,ENSG00000016602,ENSG00000129514,ENSG00000035862,ENSG00000113083,ENSG00000114353,ENSG00000140279,ENSG00000167397,ENSG00000162576,ENSG00000116132,ENSG00000106397,ENSG00000039560,ENSG00000006327, ENSG00000214711, ENSG00000136694, ENSG00000185873, ENSG00000189377, ENSG00000183844. ,

[0034] V. Analysis of the correlation between drug sensitivity and treatment response: Based on the oncoPredict package, the sensitivity of samples to specific chemotherapy drugs is predicted. For subtype-specific genes, the drug sensitivity data and gene expression data in the CTRP2 database are used to perform correlation analysis between gene expression levels and drug sensitivity. Drug molecules with significant negative correlation are potential candidate drugs.

[0035] Based on the random forest classifier constructed above, the subtypes of esophageal squamous cell carcinoma in the public dataset GSE45670 were predicted based on preoperative chemoradiotherapy.

[0036] Based on the copy number data of patients with their own immunosuppressants, the CNregions function in the R package iClusterPlus was used to remove redundant copy number regions and obtain the copy number data frame of the samples. The value was log ratio, and a random forest classifier for immune response prediction was constructed. The number of trees in the classifier was 500, which can be used to predict the molecular subtypes of esophageal squamous cell carcinoma samples with copy number data, and then to study the correlation between treatment response and subtypes.

[0037] 6. Experimental results: 1. Four subtypes were identified by similarity network fusion of multi-omics data: Based on the similarity network fusion method, four subtypes were identified. The similarity evidence of each subtype came from a single-omics or multi-omics level. The sources and distribution of similarity evidence are as follows: Figure 1 .

[0038] 2. Consensus molecular subtype study: Based on the consensus molecular subtype identification method described in the present invention, a 26*26 consensus frequency matrix is ​​obtained, such as Figure 2 As shown, the rows and columns represent the molecular subtypes of the classification system, and the colors represent the frequency of two molecular subtypes being divided into the same cluster in 100 85% sampling repeated analyses. Based on the consensus frequency matrix (network), Markov clustering identified four subtypes (ECMS1-4), and each consensus molecular classification contains the following subtype types: Figure 3 As shown, each node represents a molecular subtype defined by the corresponding typing system, its size is proportional to the sample size, and the width of the edge represents the consensus frequency between two subtypes.

[0039] In order to obtain the consensus molecular subtype of each sample, the hypergeometric test was used to analyze the significance of each sample belonging to each molecular subtype. Samples with P values ​​less than 0.2 were taken as core samples, and the consensus molecular subtype labels of samples with P values ​​greater than or equal to 0.2 were predicted, thereby obtaining the consensus molecular subtype of all samples.

[0040] 3. Molecular and clinical characteristics of subtypes: Figure 4 As shown in a, similar to previous studies, the gene with the highest mutation and copy number variation frequency is TP53 (91%), followed by CDKN2A (84%), NFE2L2 (48%), KMT2C (39%), KMT2D (30%), NOTCH1 (18%), CSMD3 (16%), PIK3CA (14%) and LRP1B (12%). Among them, CDKN2A and FAT1 contain more copy number loss events, while NFE2L2 contains more copy number amplification events. ECMS1 is enriched in NFE2L2 copy number verification events (P<0.01, 22 / 32 vs 44 / 120, chi-square test), and NFE2L2 also has a higher expression level in ECMS1 ( Figure 4 f). ECMS3 tumors are characterized by fewer copy number variation events ( Figure 4 b) and relatively low tumor mutation burden (TMB) ( Figure 4 c), which is consistent with the CNV2 subtypes defined in the study by Du et al. Nonsynonymous mutations in KMT2D are associated with poor prognosis, and the mutation frequency of this gene is not the same in the four consensus molecular subtypes, with tumors of the ECMS3 subtype containing fewer KMT2D mutations ( Figure 4 d). Comparison of mutation frequencies in the TCGA dataset, the dataset involved in the present invention (SXM-I), and another Chinese esophageal squamous cell carcinoma dataset (SXM-II) ( Figure 5 ). TP53 showed a higher frequency in the Chinese population, but no significant difference was observed between the datasets (P = 0.22). CDKN2A (Chi-square test, P = 0.0019) and FAT1 (P = 0.046, Chi-square test) had higher mutation frequencies in the Chinese population compared with the TCGA-ESCC cohort. In contrast, NEF2L2 had a lower mutation frequency (P = 0.045, Chi-square test).

[0041] Tumors of the ECMS1 subtype are characterized by dysregulation of multiple metabolic pathways, including cytochrome, glutathione, sucrose, glycoside, and uronic acid metabolism pathways ( Figure 4 g). NFE2L2 has been previously shown to be involved in metabolic reprogramming, which may partially explain the metabolic dysregulation phenotype of ECMS1. ECMS2 is characterized by dysregulation of classical tumor pathways such as WNT, MYC, and cell cycle ( Figure 4 g). ECMS3 is associated with immune activation, including high expression of myeloid-derived suppressor cells (MDSC), mature killer cells (NKC), inflammatory helper T cells (Th17 cells), and PD-1 ( Figure 4 g). The ECMS4 subtype showed mesenchymal properties characterized by TGFβ activation, VEGF pathway, and epithelial-mesenchymal transition (EMT) ( Figure 4 g). We consistently found that core signature genes of the EMT and VEGF pathways showed significantly higher expression in ECMS4 compared with the other three subtypes ( Figure 6 a and Figure 6 b), while core signature genes of the cell cycle were overexpressed in ECMS2 ( Figure 6 c).

[0042] The immune and stromal scores were calculated based on the deconvolution of tumor transcriptome data using ESTIMATE. Consistent with the GSEA results, tumors of the ECMS3 subtype showed significantly higher immune scores ( Figure 7 a), while tumors of the ECMS4 subtype had significantly higher stromal content ( Figure 7 b). The tumor purity of ECMS3 and ECMS4 is low ( Figure 7 c), suggesting that these tumors may be more susceptible to the tumor microenvironment. Further analysis using TIMER also demonstrated higher immune infiltration of B cells, CD4+, and CD8+ T cells in the ECMS3 subtype ( Figure 7 di). The same molecular features were also observed based on gene set enrichment analysis based on 90 esophageal squamous cell carcinoma proteomics data.

[0043] The survival analysis results of the four subtypes showed that patients with ECMS4 subtype showed the worst overall survival and disease-free survival. In order to verify this correlation, the present invention constructed a subtype classifier based on random forests, predicted the subtype labels of samples in the independent data set GSE53625, and performed clinical correlation analysis, and found that it was significantly correlated with overall survival ( Figure 8 a) and was found to be significantly associated with disease-free survival in the TCGA dataset ( Figure 8 b).

[0044] 4. Drug sensitivity analysis: Based on the public cell line drug database CTRP data and oncoPredict package, the sensitivity of tumor samples to specific chemotherapy drugs is predicted. The analysis results show that ECMS3 tumors are resistant to cisplatin and docetaxel, while ECMS1 tumors are resistant to paclitaxel. Fig.10). The present invention found that NFE2L2, MYC, PDCD1 (PD-1), and VEGFC were significantly upregulated in the ECMS1-4 subtypes. Similarly, based on the data of esophageal cancer cell lines in the CTRP database, the correlation analysis between drug sensitivity and gene expression showed that a negative correlation indicated possible sensitivity to the drug. VER-155008 showed a significant negative correlation with NFE2L2 expression, indicating that Necrosulfonamide might be a potential therapeutic strategy targeting NFE2L2 in ECMS1 patients ( Fig.11 a). In addition, cerulenin and PF-3758309 might be used to target MYC and VEGFC in ECMS2 and ECMS4 ( Fig.11 b and Fig.11 c).

[0045] Data analysis based on clinical trials showed that the ECMS3 subtype was significantly correlated with no complete pathologic response (<pCR) to neoadjuvant chemoradiotherapy ( Fig.12 a), which was consistent with the results obtained from computational analysis. However, the expression level of PDCD1 (PD-1) in ECMS3 was significantly higher than that in the other three groups ( Figure 4 g), and the ECMS3 subtype was significantly correlated with the treatment benefit group in immunotherapy with immune checkpoint inhibitors ( Fig.12 b), indicating that patients with ECMS3 might benefit from immunotherapy with immune checkpoint inhibitors.

[0046] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those of ordinary skill in the art should understand that they can still modify the technical solutions described in the foregoing embodiments, or perform equivalent replacements for some or all of the technical features; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A random forest classifier for molecular typing of esophageal squamous cell carcinoma, characterized by: A total of 314 characteristic genes were included in the classifier, which are represented by the Ensembl database numbers: ENSG00000177272, ENSG00000122224, ENSG00000078589, ENSG00000132185, ENSG00000105369, ENSG00000177455, ENSG00000147138, ENSG000000017312, ENSG00000172724, ENSG00000107317, ENSG00000121895, ENSG00000012124, ENSG00000163564, and ENSG000000177776. E NSG00000175463, ENSG00000162739, ENSG00000164691, ENSG00000185905, ENSG00000182162, ENSG00000048462, ENSG00000100079, ENSG0000016353 4. ENSG00000166501, ENSG00000128218, ENSG00000132465, ENSG00000005844, ENSG00000160654, ENSG00000079263, ENSG00000143297, ENSG000001 13263, ENSG00000173200, ENSG00000168918, ENSG00000124406, ENSG00000136573, ENSG00000183918, ENSG00000115085, ENSG00000122188, ENSG000 00122122,ENSG00000081237,ENSG00000117090,ENSG00000110777,ENSG00000137265,ENSG00000026751,ENSG00000115165,ENSG00000167984,ENS G00000110448, ENSG00000035403, ENSG00000162777, ENSG00000167077, ENSG00000122550, ENSG00000185811, ENSG00000100351, ENSG00000106948,ENSG00000167286,ENSG00000134242,ENSG00000143851,ENSG00000110934,ENSG00000073792,ENSG00000153283,ENSG00000247774,ENSG00000198821,ENSG00000117322,ENSG00000010671,ENSG00000182866,ENSG00000196684,ENSG00000128815,ENSG00000269404,ENSG00000167208,ENSG00000172673,ENSG00000076662,ENSG00000124256,ENSG00000170476,ENSG00000197943,ENSG00000181847,ENSG00000009790,ENSG00000124772,ENSG00000099958,ENSG00000072858,ENSG00000153885,ENSG00000172578,ENSG00000057019,ENSG00000182963,ENSG00000108405,ENSG00000198851,ENSG00000145649,ENSG00000172215,ENSG00000164283,ENSG00000179583,ENSG00000185862,ENSG00000187498,ENSG00000013725,ENSG00000176105,ENSG00000095585,ENSG00000105122,ENSG00000164938,ENSG00000047457,ENSG00000042980,ENSG00000073849,ENSG00000101082,ENSG00000173905,ENSG00000146676,ENSG00000153563,ENSG00000105374,ENSG00000106263,ENSG00000187678,ENSG00000163577,ENSG00000168421,ENSG00000143067,ENSG00000166428,ENSG00000111711,ENSG00000205744,ENSG00000171388,ENSG00000173198,ENSG00000101333,ENSG00000147689,ENSG00000117308,ENSG00000168398,ENSG00000163220,ENSG00000173156,ENSG00000162946,ENSG00000165272,ENSG00000070018,ENSG00000213160,ENSG00000114790,ENSG00000153048,ENSG00000105355,ENSG00000162366,ENSG00000175591,ENSG00000143546,ENSG00000242808,ENSG00000005884,ENSG00000162409,ENSG00000157379,ENSG00000060558,ENSG00000111424,ENSG00000148344,ENSG00000142669,ENSG00000189334,ENSG00000178038,ENSG00000153294,ENSG00000163596,ENSG00000147697,ENSG00000110047,ENSG00000030582,ENSG00000146267,ENSG00000144749,ENSG00000186847,ENSG00000073910,ENSG00000115008,ENSG00000164520,ENSG00000092820,ENSG00000158315,ENSG00000142627,ENSG00000159348,ENSG00000115009,ENSG00000062038,ENSG00000171236,ENSG00000158552,ENSG00000124145,ENSG00000124107,ENSG00000134070,ENSG00000088726,ENSG00000177191,ENSG00000176788,ENSG00000187583,ENSG00000146072,ENSG00000128596,ENSG00000011009,ENSG00000116871,ENSG00000099992,ENSG00000182612,ENSG00000100097,ENSG00000150093,ENSG00000133816,ENSG00000078098,ENSG00000206538,ENSG00000166922,ENSG00000145431,ENSG00000196923,ENSG00000115414,ENSG00000167123,ENSG00000072110,ENSG00000165617,ENSG00000122786,ENSG00000072682,ENSG00000161638,ENSG00000122884,ENSG00000120708,ENSG00000087116,ENSG00000164294,ENSG00000070404,ENSG00000128805,ENSG00000122870,ENSG00000163430,ENSG00000142173,ENSG00000120149,ENSG00000198910,ENSG00000113721,ENSG00000133110,ENSG00000080573,ENSG00000157227,ENSG00000162745,ENSG00000105472,ENSG00000117385,ENSG00000113140,ENSG00000164574,ENSG00000151135,ENSG00000142552,ENSG00000143502,ENSG00000164111,ENSG00000149257,ENSG00000166250,ENSG00000177570,ENSG00000143226,ENSG00000184012,ENSG00000107819,ENSG00000172638,ENSG00000120820,ENSG00000083444,ENSG00000142156,ENSG00000108602,ENSG00000163297,ENSG00000149380,ENSG00000130635,ENSG00000152049,ENSG00000112769,ENSG00000100196,ENSG00000115363,ENSG00000142192,ENSG00000101825,ENSG00000110723,ENSG00000106333,ENSG00000079150,ENSG00000163520,ENSG00000111817,ENSG00000128833,ENSG00000166130,ENSG00000105974,ENSG00000026025,ENSG00000138316,ENSG00000151617,ENSG00000140682,ENSG00000221968,ENSG00000131459,ENSG00000140416,ENSG00000164465,ENSG00000218336,ENSG00000168542,ENSG00000154096,ENSG00000134013,ENSG00000182667,ENSG00000142227,ENSG00000145817,ENSG00000151388,ENSG00000127249,ENSG00000165915,ENSG00000196968,ENSG00000164692,ENSG00000087303,ENSG00000184828,ENSG00000143387,ENSG00000177119,ENSG00000091136,ENSG00000134668,ENSG00000128595,ENSG00000183876,ENSG00000124225,ENSG00000196344,ENSG00000158186,ENSG00000137573,ENSG00000258947,ENSG00000088836,ENSG00000131236,ENSG00000171401,ENSG00000163359,ENSG00000129009,ENSG00000123130,ENSG00000174705,ENSG00000153561,ENSG00000150630,ENSG00000170558,ENSG00000163209,ENSG00000106624,ENSG00000182492,ENSG00000168487,ENSG00000145113,ENSG00000126458,ENSG00000135127,ENSG00000149557,ENSG00000139211,ENSG00000016602,ENSG00000129514,ENSG00000035862,ENSG00000113083,ENSG00000114353,ENSG00000140279,ENSG00000167397,ENSG00000162576,ENSG00000116132,ENSG00000106397,ENSG00000039560,ENSG00000006327,ENSG00000214711,ENSG00000136694,ENSG00000185873,ENSG00000189377,ENSG00000183844。, 2. The molecular subtype classifier of esophageal squamous cell carcinoma constructed using the random forest classifier according to claim 1, characterized in that: There are 4 molecular subtypes obtained, namely ECMS 1, ECMS 2, ECMS 3, and ECMS 4. The characteristics are: ECMS 1: abnormal metabolic pathways, NFE2L2 activation, and NFE2L2 inhibitors are suitable for drug selection; ECMS 2: upregulation of tumor classical signaling pathways and low methylation levels; ECMS 3: Few copy number variation events, low tumor mutation burden, high PD-1 expression, and benefit from immunosuppressant treatment; ECMS 4: EMT pathway activated, low overall survival and disease-free survival, and benefit from EMT and VEFG inhibitor treatment.

3. A method for constructing the esophageal squamous cell carcinoma molecular subtype classifier according to claim 2, characterized in that: The steps include: Obtain the whole genome DNA, transcriptome, and small RNA of esophageal squamous cell carcinoma pathology, perform whole genome sequencing, and preprocess the high-throughput sequencing data to obtain the esophageal squamous cell carcinoma pathology omics data; Based on the obtained omics data, namely copy number data, transcriptome, small RNA and methylation data, the omics data were integrated for subtype recognition and identification of molecular subtypes of esophageal squamous cell carcinoma; Collect the esophageal squamous cell carcinoma classification labels from the public dataset, perform consensus molecular classification analysis, and obtain the formula molecular classification subtypes; Hypergeometric test, Fisher test, and t test were applied to identify subtype-specific molecular features; Based on differential gene analysis and gene set enrichment analysis (GSEA), the relationship between each subtype and immunity, metabolism, classical tumor pathways, and epithelial-mesenchymal transition was identified; Immune microenvironment analysis was performed using the R packages deconvolution and ESTIMATE, and survival analysis was performed using Kaplan-Meier; Based on the analysis results, the molecular subtypes and characteristics of ESCC were identified, and an ESCC molecular subtype classifier was constructed.

4. The construction method according to claim 3, characterized in that: The high-throughput sequencing data were preprocessed, specifically: for whole-genome DNA sequencing, the original sequencing data were quality-controlled and filtered using FASTP software to filter low-quality reads, remove sequencing adapter sequences, and remove shorter reads; the sequencing reads were aligned with the human reference genome hg38 using BWA software to obtain SAM files; the SAM files were converted to BAM files using samtools, and the BAM files were sorted according to coordinate positions using Sambamba; the MarkDuplicates mark in the Picard toolkit was used to mark duplicate reads; the somatic copy number variation SCNV was detected in tumor samples and paired normal samples using the cnv_facets tool; For transcriptome sequencing, FASTP software was used to perform quality control filtering on the raw sequencing data, filter low-quality reads, remove sequencing adapter sequences, remove shorter reads, and correct some bases based on the overlapping regions of paired reads; the reads after quality control were aligned to the human reference genome hg38 sequence by the STAR tool, and the data was a BAM file; the gene expression abundance was estimated by the RSEM tool based on the 98th version of the Ensembl reference annotation and using the unique aligned reads, and the gene expression value was standardized to TPM; The read lengths of whole genome sequencing were first filtered for quality control using FASTP software, and then aligned using the Bismark tool. The resulting BAM file was stripped of duplicate read lengths, and the methylation signals were extracted using the Bismark methylation extractor. For small RNA sequencing data, the adapter sequences were first removed using FASTP software, and the short sequences were compared with the small RNA database miRBase and the human reference genome hg38 sequence using miRDeep2. The expression abundance of each small RNA was obtained using the quantitative module of miRDeep2.

5. The construction method according to claim 3, characterized in that: The copy number data: the log ratio values ​​of the data were extracted from the cnv_facets results, and then the redundant copy number regions were removed based on the CNregions function in the R package iClusterPlus; For methylation data, the R package methylKit was used to merge the methylation levels of all sites in all samples into a unified data matrix, and the methylation levels of 2 kb upstream and 0.5 kb downstream of the gene transcription start site were calculated. The specific calculation method was the number of methylated sequences measured at all CpG sites in the promoter region of the gene in the sample divided by the total number of sequences; if 20% of the samples were missing values ​​or had low coverage in the promoter region, the data in the promoter region were not considered.

6. The construction method according to claim 4, characterized in that: The gene expression values ​​were first normalized to the number of transcripts per million reads (TPM), and the gene expression levels of all samples were imported using the R package tximport to obtain a data matrix on gene expression; The expression of small RNAs was normalized to the number of small RNAs per 1 million reads CPM, and the small RNA expression values ​​of all samples were read to obtain the data matrix of small RNAs; If the matrix is ​​a copy number data matrix, the feature is the copy number variation region, the data is the log ratio value obtained by comparing the tumor sample with the normal sample in the region, and the number of features is the number of regions after redundancy removal; If the matrix is ​​a gene expression matrix, the features are genes, the data are the expression levels of the genes, specifically log2(TPM+1), and the number of features is the number of genes; If the matrix is ​​a methylation data matrix, the feature is the promoter region, the data is the methylation level of the promoter region, ranging from 0 to 1, and the number of features is the number of available promoters; If the matrix is ​​a small RNA expression matrix, the feature is the expression level of the small RNA, specifically log2(CPM+1), and the number of features is the total number of small RNAs; Further filtering was performed based on the median absolute deviation to select the top 1000 genes with larger variations, the top 500 small RNAs, the top 1000 promoters and all copy number regions after redundancy removal. The data of each omics were Z-value standardized, and the Euclidean distance between two samples was calculated at the level of single omics data. The R package SNF was used to construct a similarity network at the level of single omics, and network fusion was performed under the conditions of hyperparameters alpha of 0.75, K of 30 and T=20.

7. The construction method according to claim 2, characterized in that: Based on the analysis results, the method for identifying the molecular subtypes and characteristics of esophageal squamous cell carcinoma is as follows: the molecular classification label of each sample is obtained by using the spectral clustering algorithm; to evaluate the contribution of each omics to the subtype, the edges with the top 5% similarity weight are selected. If the difference in variable weights of multiple data types is less than 10%, it is considered that the edge between the two samples in the network is contributed by multiple omics data; if the edge with the highest weight is greater than the edge with the second highest weight, it is considered that the evidence of the edge between the two samples in the network comes from a specific single omics; through the similarity network fusion method based on multi-omics data, four multi-omics molecular classification subtypes were obtained: MESCC 1, MESCC 2, MESCC 3, and MESCC4; Seven groups of typing labels with PubMed IDs: 29127303, 32398863, 32580137, 31048097, 28052061, 28052061, and 36584672 of published typing studies were further collected to construct consensus molecular typing input data. 85% of the samples were randomly selected, and hypergeometric tests were performed on the correlations between the typing systems. The Markov clustering algorithm (MCL) was used to identify clusters in the network. This was repeated 100 times, and the frequencies of molecular subtypes in different typing systems being classified into the same cluster were counted to construct a 26*26 consensus frequency matrix. MCL was further used to identify the consensus molecular typing ECMS 1, ECMS 2, ECMS 3, and ECMS 4.

Images