Ethnic genomic differences in esophageal squamous cell carcinoma: Whole-exome sequencing of Han and Kazakh populations in China

Abstract

BACKGROUND

Esophageal squamous cell carcinoma (ESCC) is a cancer with a poor prognosis, characterized by distinct geographical distribution and family clustering.

AIM

To investigate if ethnic differences (Han vs Kazakh) cause molecular variations in ESCC patients via genomic sequencing 299 samples.

METHODS

Here, we sequenced samples from 299 ESCC patients collected from Henan Key Laboratory for Esophageal Cancer Research and National Key Laboratory of Metabolic Dysregulation and Esophageal Cancer Prevention and Treatment, The First Affiliated Hospital of Zhengzhou University, including Han and Kazakh ethnic groups, and performed a genomic comparative analysis of these two ethnic cohorts.

RESULTS

ESCC patients of Kazakh ethnicity present with a later age of onset compared to Han. Kazakh patients exhibit a slightly higher tumor mutation burden compared to their Han counterparts. Three genes GIGYF1, CACNA1D, and ACOT11 exhibited mutation frequencies threefold higher in Kazakh patients than in Han. This enrichment may be associated with Kazakhs’ adaptation to cold climates and consumption of high-calorie diets. Among Han patients, the apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide (APOBEC)-associated single base substitutions (SBS) 13 mutational signature is more prevalent, whereas SBS6, indicative of DNA mismatch repair deficiency, is more common in Kazakh patients. Additionally, Han Chinese patients with APOBEC-enriched tumors exhibit a significantly higher mutation load than those without. Moreover, patients lacking the APOBEC signature demonstrate superior survival probability compared to the APOBEC-enriched group.

CONCLUSION

Living environment and diet are major factors in the development of ESCC. Genomic difference may provide guidance for the formulation of clinical treatment plans for ESCC from different ethnics regions.

Key Words: Esophageal squamous cell carcinoma; Han; Kazakh; Ethnic difference; Whole-exome sequencing; Diet; Apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide; Single base substitutions 6

Core Tip: This study aimed to explore if ethnic differences (Han vs Kazakh) cause molecular variations in Chinese esophageal squamous cell carcinoma (ESCC) patients, affecting treatment response. It sequenced 299 ESCC samples, found ethnic differences in difference genes (e.g., GIGYF1), tumor mutation burden, mutational signatures (apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide/single base substitutions 6), linked to environment/diet. Results may guide ethnic-specific ESCC treatment.

INTRODUCTION

Over the past few centuries, many well-designed cancer treatment methods have been introduced. These include surgical treatment, radiotherapy, chemotherapy, phototherapy, and molecular targeted cancer therapy after 1970[1]. Cancer treatment strategies are emerging in an endless stream. However, due to the different geographical distribution of patients, epidemiological characteristics, molecular biological backgrounds, clinical stage distributions, and accessibility to medical resources, there are still differences in patient outcomes[2]. The intrinsic genetic mutations and activation of signaling pathways in tumor cells show significant racial specificity, which directly affects the response to targeted therapy[3]. The stromal components such as fibroblasts and vascular endothelial cells in the tumor microenvironment exhibit race-specific activation patterns[4]. The composition of the immune microenvironment (including T cells, myeloid-derived suppressor cells, macrophages, etc.) and the expression of immune checkpoints also show significant racial heterogeneity[5]. Due to cultural differences and economic constraints, minority patients with tumors may have difficulty accessing timely surgeries or standardized treatments[6]. Moreover, differences in socioeconomic status can also lead to variations in the implementation of treatment guidelines[7].

Disease phenotypes not only result from interactions between different genes within the host but also from the interplay between genes and environment[8]. Consequently, research on gene-environment interactions is crucial in understanding how genetic heterogeneity is influenced by various environmental agents[8,9], giving more enlightenment on the biological nature of cancer, improving the capacity to identify the susceptible gene that interacts with other factors[10,11], helping to single out why some individuals exposed to a specific agent develop cancer and others do not, and predicting the population subgroups that are more susceptible to cancer[12,13]. Current research suggests that low intake of fruits and vegetables, smoking and drinking, poor nutritional status, and the use of hot beverages may be the main risk factors contributing to the high incidence of esophageal squamous cell carcinoma (ESCC)[14,15]. Additional ESCC risk factors include comorbid gastrointestinal conditions, prior head/neck malignancies, high-risk human papilloma virus infection, tylosis, socioeconomic disadvantage, and suboptimal oral hygiene. The relative contribution of specific risk factors to disease development may vary significantly across geographic regions and ethnic populations[16,17].

Esophageal cancer is an aggressive and invasive disease[18], ranks the sixth leading cause of cancer-related death worldwide[19] with ESCC accounting for over 90% of total esophageal cancer in China. ESCC has poor prognosis with the 5-year survival rates around 20%-30% probably due to the difficulty in early diagnosis and the lack of effective therapy[19-21]. However, the survival rate of the Kazakh ethnic group is even significantly lower than that of the Han ethnicity (5%-13% vs 20%-30%)[22,23].

Since 2012, there have been dozens of investigations published using the whole-genome sequence or whole-exome sequence (WES) strategy to explore the genetics of ESCC. These studies depicted the general mutational landscape of ESCC, including the significantly mutated genes such as TP53, CDKN2A, EP300, PIK3CA, and NOTCH1, the commonly influenced pathways such as phosphatidylinositol 3-kinase (PI3K)-protein kinase B (AKT) axis, cell cycle, and histone modification, and the commonly identified age-related and apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide (APOBEC) enzymes-related mutational signatures[24,25].

However, what are the reasons for the significant differences in survival between the Kazakh ethnic group and the Han ethnic group? Is the current treatment guideline for ESCC in China applicable to patients of other ethnic groups with ESCC? Analysis of clinical variables-related genomic features is essential for translational research. Here, we compared the genomic features difference between Han and Kazakh ESCC patients in China, to distinguish the molecular basis underlying the striking racial disparities in ESCC, and discuss the reason for the difference.

MATERIALS AND METHODS

Patients

A total of 299 ESCC patients were included in this study, with 290 from Han region and 9 from Kazakh region. All the patients were collected from the database of Henan Key Laboratory for Esophageal Cancer Research and National Key Laboratory of Metabolic Dysregulation and Esophageal Cancer Prevention and Treatment of the First Affiliated Hospital, Zhengzhou University. Patient inclusion criteria: (1) Underwent radical surgical treatment; (2) Complete clinicopathological information available; and (3) Postoperative pathological diagnosis confirmed as primary ESCC with a definite histopathological type. This study was approved by the Ethics Committee of Zhengzhou University (Ethics approval number: No. 2024-KY-2336-002).

Patient basic information was collected, including sex and age; clinicopathological data encompassed tumor differentiation grade, lymph node metastasis status, and pathological stage (Union for International Cancer Control, 2017). Follow-up was conducted for 299 patients with detailed addresses or contact information on record, utilizing telephone interviews and home visits. Follow-up began at diagnosis, with evaluations every 3 months initially, transitioning to annual assessments starting from the second year. The follow-up period concluded in May 2024. Death was defined as the endpoint event.

DNA extraction, library preparation, and WES

WES was performed on 299 paired fresh-frozen ESCC tumor tissues and matched blood samples. Genomic DNA was extracted from both tissue and blood samples using the DNeasy blood and tissue kit (Qiagen, Netherlands; Cat. 69504), according to the manufacturer’s instructions. DNA concentration was precisely quantified using Qubit 2.0. Samples meeting the criteria of ≥ 20 ng/μL concentration and a minimum total mass of 0.6 μg were used for library preparation by Novogene Bioinformatics Technology Co. Ltd (Beijing, China).

Whole exome libraries were constructed as described previously[26]. Briefly, DNA underwent fragmentation, end repair, 3’ end adenylation, and adapter ligation. Following size selection and polymerase chain reaction amplification using the NEB Next^® Ultra DNA library prep kit for Illumina (New England Biolabs, United States; Cat. E7530 L), the initial libraries were amplified with a high-fidelity polymerase to generate a sufficient quantity of the exome library. Upon completion of library construction, quantification was performed using Qubit 3.0 (Invitrogen, United States) and the NGS3K/Caliper system. Libraries were then subjected to cluster generation and sequenced on the Illumina NovaSeq 6000 platform, generating 150 bp paired-end reads.

Statistical analysis

WES was performed on DNA extracted from all tumor tissues and matched blood samples, utilizing a capture system covering 64 Mb of the human genome (GRCh37). Stacked bar plots depicting clinicopathological characteristics were generated using the ggplot2 package in R (version 4.1.2). Analysis of variance (ANOVA) and χ² tests were employed to compare the clinicopathological features of ESCC patients across different geographic regions. Waterfall plots were constructed using the maftools package in R[27]. The mafCompare function was used to identify differentially mutated genes between patient groups, defined as genes with a mutation frequency ≥ 5% and P < 0.05. The somatic interactions function was employed to analyze the mutual exclusivity and co-occurrence of mutations. Mutation sites within genes were visualized using the lollipopPlot2 function. APOBEC enrichment scores were estimated using the trinucleotide matrix function, and samples were classified as either APOBEC-enriched or non-APOBEC-enriched. These two groups were then compared to identify differentially altered genes.

Signature analysis comprised the following steps: (1) Estimate signatures: This function runs non-negative matrix factorization (NMF) over a range of potential signature numbers and assesses the goodness of fit using the cophenetic correlation coefficient; (2) Plot cophenetic: This function generates an elbow plot to aid in determining the optimal number of signatures. The optimal number typically corresponds to the point where the cophenetic correlation exhibits a significant drop; (3) Extract signatures: This function applies NMF to decompose the mutational matrix into n signatures, where n is determined based on steps (1) and (2). If a reliable estimate of n is already available, these preceding steps can be omitted; (4) Compare signatures: The extracted signatures from step 3 are compared against known signatures from the Catalogue of Somatic Mutations in Cancer (COSMIC) database (single base substitutions (SBS) 1, SBS2, SBS3 etc.). Cosine similarity is calculated to identify the best-matching COSMIC signature(s) for each extracted signature; and (5) Plot signatures: Visualizes the extracted mutational signatures.

RESULTS

Patients’ distribution

In this study, we generated WES data using the Illumina platform for both Han (n = 290) and Kazakh (n = 9) ethnic groups with ESCC. To robustly identify ethnic-biased mutated genes, we first assessed the distributions of key clinicopathological factors across these populations. Patient age at diagnosis, geographic region, and tumor differentiation exhibited significant differences among the groups (assessed by ANOVA for age and χ² tests for region and differentiation; Figure 1 and Supplementary Figure 1).

Figure 1 Distribution of five biological factors in two patient groups, justifying a need to balance these confounders. A: Diagnosis; B: Sex; C: Area; D: Smoking; E: Alcohol drinking; F: Differentiation. Analysis of variance test was used to calculate the P value for age at diagnosis, and χ² test was used for the other factors. HIA: High incidence area; LIA: Low incidence area.

Significantly mutated gene of ESCC

In the exonic coding regions of Han Chinese ESCC patients, a total of 242460 somatic mutations were identified, comprising 191962 single-nucleotide polymorphisms (SNPs), 18510 insertions, and 31898 deletions. Among Kazakh ESCC patients, 18000 somatic mutations were detected in exonic coding regions, consisting of 14311 SNPs, 1006 insertions, and 2683 deletions. Figure 2 presents key visualizations generated using the Maftools package. Figure 2A, an oncoplot, displays the top 20 significantly mutated genes identified by the MutSigCV algorithm. Genes are ranked by mutation frequency, while samples are sorted and grouped by ethnicity, revealing differences in the mutational spectrum both between and within ESCC subtypes. Of these 20 genes, 13 were reported as significantly mutated genes in at least one previous ESCC genome study, including TP53 (82%)[25,28], TTN (40%)[28], MUC16 (25%)[28], KMT2D (23%)[25,28], NOTCH1 (19%)[25,28], FLG (14%)[29], CSMD3 (14%)[28,29], PCLO (13%)[30], CDKN2A (13%)[25,31], LRP1B (12%)[28,30], CSMD1 (12%)[25], AHNAK2 (12%)[32], FAT1 (12%)[28,29] have been reported in other ESCC. In our study, seven of these genes have not been reported in other studies on esophageal cancer, including to OBSCN (14%), SYNE1 (14%), ZFHX4 (14%), PER3 (14%), RYR2 (13%), USH2A (13%), CCDC168 (13%). Among them, ZFHX4 and RYR2 are frequently mutated in colon adenocarcinoma, and its high expression in colon adenocarcinoma is significantly associated with shortened progression-free survival[33]. Therefore, it has the potential to serve as a predictive marker for the efficacy of immunotherapy. The other five genes have also been identified in other digestive tract tumors[33], but no corresponding targeted therapies have been proposed yet. And all the 11 genes (GIGYF1, PABPC1, CACNA1D, DNHD1, NF1, ACOT11, CDC27, CGA, NFATC3, PATL1, TRIP12) are significantly enriched in Kazakh patients. TRIP12, NF1, PABPC1 and CACNA1D are core genes that directly or indirectly participate in epigenetic regulation. Targeting TRIP12 can enhance the sensitivity of epigenetic drugs (such as DNA methyltransferase inhibitors), and combined use may reverse chemotherapy resistance[34]. Patients with NF1 mutation and esophageal cancer respond better to histone deacetylase inhibitors (such as vorinostat), which can restore the expression of tumor suppressor genes[34]. GIGYF1 can combine with AKT to form macromolecular complex, the AKT signaling pathway plays an important role in regulating cell apoptosis and autophagy[35], and calcium channel blockers have the potential to serve as adjunctive agents to epigenetic drugs. PABPC1 plays a critical role in Chinese ESCC malignant progression by interacting with eukaryotic translation initiation factor 4G (eIF4G) to regulate IFI27 messenger RNA stability and promote angiogenesis via exosomal miR-21-5p/CXCL10[36]. CACNA1D specifically had low expression in esophageal tumors in other study, however high expression in Kazakh ESCC patients[37]. DNHD1 was not found in other esophageal tumors. Loss-of-function mutation in NF1 was found in ESCC patients[38].

Figure 2 Significantly mutated genes in esophageal squamous cell carcinoma. A: Significantly mutated genes (SMGs) identified by MutSigCV on a combined cohort of Han and Kazakh ethnics groups, Chinese whole-exome sequence samples. Each column denotes an esophageal squamous cell carcinoma patient, and each row is a gene. On top is the number of somatic mutations per sample. On the right are the mutation frequencies of each SMG. The bar plot on the top shows the composition of mutations in the gene. Genes are ordered by their mutation frequencies. The table on the right shows the number of mutation population of total, Han and Kazakh; B and C: Visualizing mutation summary. This maftools plot shows a summary of the multiple alignment format file. Highlighting the most mutated genes, single-nucleotide variant class, and variant classification distributions within a tumor type for Han and Kazakh patients; D and E: Fraction of pathway and sample affected for Han and Kazakh. SNV: Single nucleotide variant; SNP: Single nucleotide polymorphism; INS: Insertion; DEL: Deletion; PI3K: Phosphatidylinositol 3-kinase; TGF: Transforming growth factor; NRF2: Nuclear factor erythroid 2-related factor 2.

As the most important tumor suppressor, the non-silent mutation frequency of TP53 was 82% in 299 tumors, which was consistent with previous reports[39-41]. The types and relative positions of confirmed somatic mutations are shown in the transcript of TP53[40]. Red dots, nonsense mutations, and green dots, missense mutations, are the main variant classification for Han ESCC patients. However, green dots, missense mutation and blue dots, frame shift del are the main variant classification for Kazakh ESCC patients.

Recognizing somatic mutations as the molecular basis of tumor mutational burden (TMB), we characterized genome-wide variations by analyzing somatic mutation data from ESCC patients. An overview of the analytical strategy is presented in Figure 2B and C. Analysis revealed that missense mutations, splice site mutations, and frameshift deletions were the three most frequent variant types in Han Chinese patients (Figure 2B). In contrast, for Kazakh patients, splice site mutations ranked as the second most frequent variant type (Figure 2C). Single-nucleotide variants (SNVs) constituted the vast majority of all variant types in both ethnic groups (Figure 2B and C), with C > T substitutions being the most common SNV class (Figure 2B and C). Moreover, we displayed the number of mutated bases in each of the patients, with a median value of 85 for Han group, 210 for Kazakh (Figure 2B and C). The top 10 mutated genes in Han were TP53, TTN, MUC16, KMT2D, PER3, OBSCN, NOTCH1, SYNE1, CSMD3, and ZFHX4 (Figure 2B). The top 10 mutated genes in Kazakh were TTN, MUC17, MUC16, TP53, MUC12, DNHD1, PABPC1, NF1, GIGYF1 and FAT4 (Figure 2C)[42].

Genes were mapped to known oncogenic pathways, and pathway enrichment analysis was performed using the oncogenic pathways function in Maftools. This analysis quantifies the frequency of altered genes within each pathway and calculates the fraction of samples harboring mutations in genes belonging to that pathway. Pathways are then ranked by the number of mutated genes they contain. For Han Chinese patients, the most frequently altered pathways (with 100% of their constituent genes mutated in at least one sample) were transforming growth factor-beta, TP53, and nuclear factor erythroid 2-related factor 2 (NRF2). However, these pathways were mutated in 12%, 82%, and 9% of samples, respectively. Other enriched pathways included MYC (mutated in 11% of samples), cell-cycle (21%), PI3K (27%), and Hippo (58%) (Figure 2D). For Kazakh Chinese patients, the transforming growth factor-beta, TP53, and NRF2 pathways were also commonly altered. Enrichment levels for other pathways were: MYC (mutated in 11% of samples), cell-cycle (22%), PI3K (44%), and Hippo (78%). Within the pathway plots (Figure 2D and E), tumor suppressor genes are denoted in red font, while oncogenes are shown in blue.

The median variant allele frequency (VAF) of mutations in the top 20 genes for Han Chinese and Kazakh ESCC patients is presented in Supplementary Figure 2A and B, respectively. Among Han patients, genes with a median VAF ≥ 25% included NOTCH1, CDKN2A, and TP53. In contrast, only TP53 met this threshold among Kazakh patients. Conversely, the remaining mutated genes in both groups exhibited a median VAF < 25%, suggesting their presence within a subpopulation of tumor cells.

Next, we investigated the co-occurrence of mutated genes. Significant associations (P < 0.05) were identified in Han Chinese patients (Supplementary Figure 2C). Specifically, mutations in MUC16 co-occurred with mutations in RYR2, CCDC168, OBSCN, FLG, SYNE1, and KMT2D. Mutations in KMT2D co-occurred with mutations in AHNAK2, FAT1, RYR2, and SYNE1. Mutations in ZFHX4 co-occurred with mutations in FAT1, CSMD1, and USH2A. Mutations in OBSCN co-occurred with mutations in KMT2C, FAT1, and CCDC168 (P < 0.05, Supplementary Figure 2C). In Kazakh patients (Supplementary Figure 2D), significant co-occurrence associations (P < 0.05) were observed. Mutations in TP53 and mutations in MYCBP2 and DNAH11. Mutations in DNAH11 and mutations in NF1 and FAT4. Mutations in FAT4 and mutations in NF1 (P < 0.05, Supplementary Figure 2D).

Difference in mutation type

The transition/transversion (Ti/Tv) ratio is commonly used to assess the quality of variant calling. The transition and transversion plots (Figure 3A and B) categorize SNVs into six subtypes. We observed that the Ti/Tv ratio was approximately 2.0 in Han patients, but reached 3.0 in Kazakh patients. This difference may be attributed to a higher proportion of C > T transitions in the Han cohort compared to the Kazakh cohort. The elevated Ti/Tv ratio in the Han population may be associated with cytosine methylation in CpG islands, which promotes C > T transitions[43,44]. TMB was significantly higher in Kazakh ESCC patients compared to Han patients (Figure 3C and D). This elevated TMB implies that tumors in Kazakh patients are likely to produce a greater diversity and quantity of neoantigens. Consequently, these neoantigens have a higher probability of being recognized by the host immune system. Following immune checkpoint inhibitor therapy, which activates the body’s endogenous anti-tumor immune response, this enhanced neoantigen recognition may translate to a greater probability of tumor cell elimination[45]. No clustering of mutations was identified for Han and Kazakh patients (Supplementary Figure 3A for Han patients, and Supplementary Figure 3B for Kazakh patients).

Figure 3 Variant calling and mutation burden. A and B: The transition and transversion plot illustrates the distribution of single-nucleotide variants (above) in esophageal squamous cell carcinoma, categorized into six transition and transversion events. The core purpose is to describe the distribution characteristics (such as central tendency, dispersion degree, outliers, etc.) of a single set of data (for example, a certain transformation/transposition event). The stacked bar chart (bottom) displays the mutational spectrum distribution of each sample in the multiple alignment format file. Present the “proportion distribution of categorical data” in a visual way to assist in understanding the heterogeneity of the mutation spectrum among different samples; C and D: Mutation burden. The core is to describe the distribution characteristics of tumor mutational burden among the samples, including the central tendency (such as the median) and the degree of dispersion (such as the range of variation). TMB: Tumor mutation burden; MB: Mega base.

Difference gene/hub gene /functional domain

Clinical parameters contribute to tumor heterogeneity within a single cancer type. We utilized the mafCompare function to identify differentially mutated genes between the two cohorts, comparing the mutation frequency of each gene using Fisher’s exact tests. This comparison identified 11 differentially mutated genes between Han and Kazakh patients (P < 0.001, Figure 4A; P < 0.005, Supplementary Figure 4A). And all the 11 genes (GIGYF1, PABPC1, CACNA1D, DNHD1, NF1, ACOT11, CDC27, CGA, NFATC3, PATL1, TRIP12) are significantly enriched in Kazakh patients.

Figure 4 Difference genes and protein domain. A: Difference gene between the two ethnics groups (P < 0.001); B: Gene-gene interaction network between difference and significant mutation genes (P < 0.05, and the mutation frequency ≥ 5%); C-H: Schematics of protein changes in ethnics-biased gene products. Scattered distribution of somatic mutations in MUC16, SETD1B, DNHD1, DNAH11, DNAH17, DYNC1H1. OR: Odds ratio; CI: Confidence interval.

We then used the 37 differential genes and significantly mutated genes (P < 0.05 and mutation frequency ≥ 5%) to construct gene-gene association networks using the STRINGdb package (Figure 4B). Based on node centrality and interaction degree within the network, the interaction network was clustered and decomposed, yielding five distinct sub-networks (Supplementary Figure 4B). Among them, sub-network I consists of ZFH3 and SYNE2; sub-network II comprises KCNMA1, RYR2, OBSCN, and MYO18B; sub-network III includes OTOG, MUC12, MUC16, and MUC17, among which MUC16 is the central gene of the network; sub-network IV is composed of SRCAP, SETD1B, CHD4, and KMT2A, among which SETD1B is the central gene of this sub-network; sub-network V is made up of HRNR, DYNC1H1, NF1, DNAH17, DNHD1, and DNAH11, among which DYNC1H1, DNAH17, DNHD1, and DNAH11 are the central genes of this sub-network. Notably, sub-network III contains MUC16, a gene implicated in esophageal cancer onset. This sub-network may represent an important regulatory module potentially involved in esophageal cancer pathogenesis[46].

MUC16 mutation sites were uniformly distributed across the gene without enrichment in functional domains. Frameshift deletions represented the predominant mutation type in both Han and Kazakh patient cohorts (Figure 4C). Similar to MUC16, we also observed differences in different mutation sites in functional domains in SETD1B, DNHD1, DNAH11, DNAH17, and DYNC1H1 (Figure 4D-H).

Mutation signature

From Figure 5A and B, we can see that SBS6 is associated with DNA mismatch repair deficiency in both groups, but the cosine similarity is higher in the Kazakh group (0.858 vs 0.843), indicating that the SBS6 characteristics of the Kazakh group are more similar to the known DNA mismatch repair deficiency characteristics. SBS5 exists in both groups, but the reason is unknown. The SBS5 characteristics of the Han ethnicity have a higher similarity to the known characteristics (0.829 vs 0.686). SBS13 and SBS18 only occur in the Han ethnicity, respectively related to APOBEC cytidine deaminase and possible reactive oxygen species (ROS) damage. SBS30 only occurs in the Kazakh group and is associated with base excision repair deficiency caused by NTHL1 gene inactivation mutation. Figure 5C shows that there is a significant difference in mutation burden between APOBEC-enriched samples and non-APOBEC-enriched samples in Han Chinese esophageal cancer patients (P < 0.001). In APOBEC-enriched samples, the mutation frequencies of genes such as TP53, PIK3CA, USH2A, and MYO1F are higher, while in non-APOBEC-enriched samples, the mutation frequencies of these genes are lower. Non-APOBEC-enriched Han patients has better survival than APOBEC-enriched Han patients (Figure 5D). Supplementary Figure 5A shows best possible value of Han (n = 4) and Kazakh patients (n = 3), and Supplementary Figure 5B shows that the cosine similarity between the mutation characteristics of different samples and the verified characteristics varies, indicating that there are multiple different mutation characteristics in Han Chinese patients. Supplementary Figure 5C also shows that there are multiple different mutation characteristics in Kazakh patients, but the similarity of certain characteristics is different compared to Han Chinese patients.

Figure 5 Signature plot. A: Han; B: Kazakh patients. Each row presents a specific single base substitutions (SBS) signature (Catalogue of Somatic Mutations in Cancer reference signatures with known etiologies) and its best match with known signatures. The distribution of each signature shows the relative frequencies of the six single-base substitution types. Different SBS signatures are associated with different pathogenic mechanisms, such as apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide (APOBEC) cytidine deaminase, unknown mechanisms, possible reactive oxygen species damage, and DNA mismatch repair deficiency, etc.; C: Differences between APOBEC enriched and non-enriched patients for Han patients (none of the samples are enriched for APOBEC in Kazakh patients); D: Survival plot of Han patients. SBS: Single base substitution; APOBEC: Apolipoprotein B messenger RNA-editing enzyme catalytic polypeptide.

DISCUSSION

In general, cancer control strategies that work in high-income areas often don’t work in low-income areas. Because of differences in disease characteristics, health system capacity, sociocultural, treatment completion rates, and drug availability[47]. Our study found that diagnosis age was higher in Kazakh than in Han, and the difference was statistically significant (P = 0.016). The Kazakh ethnic group in China is mainly distributed in Xinjiang Uygur Autonomous Region. Due to the relatively low level of economic development in Xinjiang Uygur Autonomous Region, the income of the Kazakh ethnic group is generally lower than the average level of China. However, the relevant policies of the Chinese government to support the economic development of ethnic minorities, including financial support, tax incentives, poverty alleviation, special industries and employment help, infrastructure construction, legal protection, education and personnel training, rural revitalization and so on[48-52]. These measures have greatly improved the economy and living standards of the Kazakh people in Xinjiang Uygur Autonomous Region.

GIGYF1 gene is located in 7q22, with a total length of about 10kb, and binding to growth factor receptor-bound protein 10, and interacts with both the insulin and IGF1 receptors. Reports have proved that GIGYF1 is significantly expressed in autism[53], type 2 diabetes[54,55], and gastric cancer[35]. We first seen GIGYF1 highly expressed in ESCC patients of Kazakh. GIGYF1 promotes tumor growth by activating the extracellular regulated protein kinases (ERK) and AKT signaling pathways. Its low expression can significantly inhibit the phosphorylation of ERK1/2 and AKT, and block the downstream pro-survival signals[35]. This indicator can be used as a basis for evaluating the feasibility of targeted intervention for gastric cancer[56]. Molecular alterations or potential targets may serve as candidates for targeted therapy in ESCC. PABPC1, poly (A) binding protein cytoplasmic 1, PABPC1 expression was upregulated in ESCC tissues compared with in normal esophageal epithelial tissues. Elevated PABPC1 expression was correlated with tumor cell differentiation and poor prognosis in patients[57,58]. PABPC1 plays a critical role in ESCC malignant progression by interacting with eIF4G to regulate IFI27 messenger RNA stability and promote angiogenesis via exosomal miR-21-5p/CXCL10[36].

Previous studies proved that the dietary calcium intake was inversely associated with the risk of esophageal cancer [odds ratio (OR) = 0.80, 95% confidence interval: 0.71-0.91]. And the protective function of dietary calcium intake was observed in esophageal squamous cell cancer[59,60]. In our study, the occurrence probability of CACNA1D gene in the patients with esophageal cancer of the Kazakh ESCC group is 3 times that of the Han ESCC group (OR = 35.893, P < 0.001). CACNA1D, the gene encoding the pore-forming α1-subunit of Cav1.3 calcium ion (Ca²⁺) channels, highly mutated in Kazakh ESCC patients, which may be associated with low calcium intake in Kazakh because of poorly infrastructure, weak medicines purchasing, and far away distribution of pharmacies. But its actual mechanism of its action needs to be explored through other molecular experiments. L-type calcium channel blockers (such as amlodipine and nifedipine) can block downstream carcinogenic signals by inhibiting intracellular Ca²⁺ influx. CACNA2D1 (also a calcium channel protein) is highly expressed in gastric cancer and promotes tumor growth. Its inhibitor (such as amlodipine) combined with chemotherapy can enhance the therapeutic effect[61,62]. L-type calcium channel blockers (such as amlodipine and nifedipine) can block downstream carcinogenic signals by inhibiting the influx of Ca²⁺[63]. Therefore, it provides feasibility for exploring the combination of calcium channel blockers with standard treatments in clinical trials[64].

Lipid metabolism, which is associated with the tumor microenvironment[65], may affect the responsiveness of the cancer cells to immunotherapy[66]. Previous study explored the mechanisms involved in the regulation of lipid metabolism in ESCC, identified differentially expressed genes and developed a prognosis predictive model for aiding clinical decision-making, which included ACOT11[67]. Its connection to ESCC and lipid metabolism in this context is an area ripe for exploration. ACOT11 is an enzyme involved in the hydrolysis of acyl-CoA compounds, influencing fatty acid metabolism[67]. It is highly expressed in brown adipose tissue (BAT) and is affected by environmental temperature and food consumption[68,69]. In higher organisms, this gene is mostly distributed in the cytoplasm, mitochondria, peroxisome and endoplasmic reticulum. In mammalian tissues, this gene is widely distributed in the brain, liver, kidneys, heart, lungs, steroids, etc.[70,71]. BAT is a kind of adipose tissue used to break down excessive energy and thereby reduce energy storage. It contains a relatively large number of mitochondria. Among them, the energy metabolism of UCP-1[72] and the production of adenosine triphosphate provide a basis for non-shivering thermogenesis[73]. Kirkby demonstrated that in the BAT of mice, the expression of ACOT11 was upregulated at low temperatures[71]. The main function of ACOT11 is to reduce energy consumption and preserve heat in the body, which may be owing to the low-temperature environment and high-energy diets in Kazakh. Therefore, we can treat patients with ACOT11 mutations through two methods[2,74]. The first approach is to inhibit the fatty acid oxidation pathway. This can be achieved by knocking down acetyl coenzyme A carboxylase 1 to suppress fatty acid synthesis, and in combination with ACOT11 inhibition, it can effectively block lipid flow[75,76]. The second approach is immune microenvironment modulation. The APOBEC-ACOT11 axis, where ACOT11 mutations may enhance APOBEC activity, leading to an increase in TMB and an increase in the release of new antigens. Alternatively, by activating the cyclic guanosine monophosphate-adenosine monophosphate synthase-stimulator of interferon genes pathway, it promotes the secretion of type I interferons and the infiltration of T cells. Combined with immune checkpoint inhibitors[77-79]. The response rate of programmed death 1 (PD-1) inhibitors was significantly increased in ESCC patients with high TMB (27% vs 6% in patients with low TMB)[30,80].

The analysis results of mutation signatures in the Han and Kazakh populations with ESCC revealed significant differences in the mechanisms of esophageal cancer occurrence between the two groups, which might be related to genetic background, environmental exposure, and lifestyle factors. In both the Han and Kazakh populations, the SBS6 signature related to DNA mismatch repair deficiency was detected, and the cosine similarity was higher in the Kazakh population (0.858 vs 0.843). This indicates that DNA mismatch repair deficiency may play an important role in the occurrence of esophageal cancer. Defects in the DNA mismatch repair system can lead to genomic instability, increasing the risk of mutation accumulation, thereby promoting the occurrence of cancer[81]. The SBS5 signature exists in both populations, but the reason is unknown. The SBS5 signature in the Han population has a higher similarity to the known signature (0.829 vs 0.686). This may suggest that the SBS5 signature has a more significant role in esophageal cancer in the Han population, but its specific pathogenic mechanism still requires further study. Future research can combine genomics and functional experiments to explore the potential pathogenic mechanism of the SBS5 signature[82]. The Han-specific SBS13 and SBS18 signatures may be related to the specific environmental exposure or genetic factors of this ethnic group. The SBS13 signature is related to the APOBEC cytidine deaminase (C > G), and the APOBEC family enzymes have been found in various cancers, and their activity may be affected by viral infections, inflammatory responses, and other factors[83]. In Han Chinese patients with ESCC, the mutation characteristics mediated by APOBEC3A/APOBEC3B (such as SBS2/SBS13) may be directly related to the chronic inflammation-induced ROS burst. ROS causes DNA single-strand breaks by oxidizing bases (such as 8-oxoguanine), and under replication stress, activates the cytidine deaminase activity of APOBEC3A, thereby catalyzing the transition from TC to TT/TA mutations. This process forms a positive feedback loop of oxidative damage-APOBEC activation-mutation accumulation[84]. Chronic infections of the upper digestive tract (such as Helicobacter pylori) may promote neutrophil infiltration by activating the nuclear factor kappa-B pathway, resulting in excessive production of mitochondrial ROS. This microenvironment is associated with the high expression of APOBEC3B in gastric adenocarcinoma, and a similar mechanism may exist in Han Chinese ESCC: Helicobacter pylori infection-gastric acid reflux-accumulation of ROS in esophageal mucosa-enhancement of APOBEC mutation characteristics[85-88]. For the high-risk Han population, the antioxidant intervention proposed by Liu et al[89] and Ou et al[90] can be adopted: Such as berberine (targeting NOX4) or green tea polyphenols (cleaving ROS), which can block the activation chain of APOBEC. APOBEC inhibitors (such as small molecule antagonists of APOBEC3G) combined with PD-1 antibodies may enhance the efficacy. The SBS18 signature may be related to ROS damage. ROS are produced during cellular metabolism, and excessive ROS can cause DNA damage and promote the occurrence of cancer[91]. The Kazakh-specific SBS30 signature is related to the inactivation mutation of the NTHL1 gene, which leads to a base excision repair defect. The NTHL1 gene encodes a protein involved in base excision repair, and its inactivation mutation can lead to DNA repair defects and increase the risk of mutation accumulation. This suggests that the occurrence of esophageal cancer in the Kazakh population may be related to specific genetic mutations or environmental factors.

Furthermore, through large-scale genomic data, neural network learning[92], generative adversarial networks learning[93], optimization algorithms, and the verification of important markers, applying generative adversarial networks to integrate the methylation and single-cell transcriptome data of Han/Kazakh ethnic groups, to analyze the epigenetic memory of environmental exposure, it can provide important references for the clinical treatment and diagnostic guidelines of esophageal cancer patients of different races.

CONCLUSION

Living environment and diet are major factors in the development of ESCC. Genomic differences may provide guidance for the formulation of clinical treatment plans for ESCC from different ethnics regions.

ACKNOWLEDGEMENTS

We would like to thank Wu HL (University of Chinese Academy of Sciences) and Xia CZ (School of life sciences, Zhengzhou University) for kindly providing analysis support.