Characteristics of genomic mutations and signaling pathway alterations in thymic epithelial tumors.

BACKGROUND: To elucidate the mechanisms of thymic epithelial tumor (TET) canceration by characterizing genomic mutations and signaling pathway alterations.
METHODS: Primary tumor and blood samples were collected from 21 patients diagnosed with TETs (thymoma and thymic cancer), 15 of whom were screened by nucleic acid extraction and whole exon sequencing. Bioinformatics was used to comprehensively analyze the sequencing data for these samples, including gene mutation information and the difference of tumor mutation burden (TMB) between thymoma and thymic carcinoma groups. We performed signaling pathway and functional enrichment analysis using the WebGestalt 2017 toolkit.
RESULTS: ZNF429 (36%) was the gene with the highest mutation frequency in thymic carcinoma. Mutations in BAP1 (14%), ABI1 (7%), BCL9L (7%), and CHEK2 (7%) were exclusively detected in thymic carcinoma, whereas ZNF721 mutations (14%) and PABPC1 (14%) were found exclusively in thymoma. The mean TMB values for thymic carcinoma and thymoma were 0.722 and 0.663 mutations per megabase (Mb), respectively, and these differences were not statistically significant. The ErbB signaling pathway was enriched in the thymoma and intersection groups, and pathways of central carbon metabolism in cancer, longevity regulating and MAPK signaling were only found in the thymoma group, while pathways in cancer (hsa05200) was found in the thymoma and thymic carcinoma groups.
CONCLUSIONS: Multiple differences in somatic genes and pathways have been identified. Our findings provide insights into differences between thymoma and thymic carcinoma that could aid in designing personalized clinical therapeutic strategies. 2021 Annals of Translational Medicine. All rights reserved.

Introduction

As the most common anterior mediastinal tumors in adults, thymic epithelial tumors (TETs) account for approximately 0.2% to 1.5% of all malignancies (1,2). According to the World Health Organization (WHO) histologic classification of TETs, these tumors can be divided into types A, AB, B1, B2, B3, Thymic carcinoma, Squamous cell carcinoma, Undifferentiated carcinoma, Large-cell neuroendocrine carcinoma, NOS and Micronodular Thymoma (MNT) (3,4). Compared with thymoma, thymoma has the characteristics of inconspicuous lobular structure, larger nucleus, unbalanced nucleus to pulp ratio, and obvious heterogeneity (5). The most common tissue types of thymic carcinoma are squamous cell carcinoma and undifferentiated carcinoma (6).

Generally, thymomas are prone to local recurrence rather than metastasis. At present, the main treatment for thymoma is surgical resection with postoperative radiotherapy and chemotherapy (7,8). However, even with surgery, chemotherapy, and radiotherapy, thymic carcinoma is associated with a high risk of recurrence and death (9). As a result, the 5-year median survival rate of thymoma is better than that of thymic carcinoma, standing at 69% and 36%, respectively (10). Despite the mutation information and signaling pathways of thymic epithelial tumors have been studied (11), the differences in molecular characteristics and signaling pathway changes between these two tumor types are currently unknown.

Dozens of biomarkers related to checkpoint inhibitors have been identified to date. Among them, tumor mutation burden (TMB), programmed cell death-1 (PD-1) ligand 1 (PD-L1), and microsatellite instability/deficient mismatch repair (MSI/dMMR) have been verified by phase III clinical trials and are widely used in clinical practice (12-14). Although it has been reported that TMB-high (TMB-H) alone is not suitable for predicting the immunotherapeutic effect against solid tumors (15), TMB can serve as a biomarker for predicting PD-1/PD-L1 immune response (16,17). In adults, the TMB of thymic tumors is generally 0.48 mutations per megabase (18). PD-L1 is commonly expressed in TETs, including in 23–92% of thymomas and 36–80% of thymic carcinomas (19-25).

Multiple studies have investigated the relationship between PD-L1 expression and survival/prognosis in patients with TETs; however, their results have been conflicting. In their study, Padda et al. reported patients with high PD-L1 expression to have poor overall survival (19). However, Weissferdt et al. failed to find a correlation between PD-L1 expression and prognosis (25), while Arbour et al. and Yokoyama et al. reported a good prognosis in patients with high PD-L1 expression (22,26). Studies have also shown complex signaling pathway crossover in the tumor canceration process (27), although few studies have investigated signaling pathways in the process of thymoma canceration.

To better understand the potential mechanisms and biological characteristics of TETs, next-generation sequencing (NGS) and signal pathway analyses of TET samples were conducted to investigate somatic gene mutations, TMB, the transition/transversion (TI/TV) ratio, and mutated gene signaling pathways in patients with thymoma and thymic carcinoma. Our findings may provide effective guidance for clinical intervention, thereby improving the long-term prognosis of patients with TETs. We present the following article in accordance with the MDAR checklist (available at https://dx.doi.org/10.21037/atm-21-5182).

Methods

Patient specimen acquisition

Twenty-one patients with TETs including 9 with thymomas and 12 with thymic carcinomas were enrolled into this study. PD-L1 expression was detected in all 21 patients, and whole-exome sequencing (WES) was performed for 15 patients.

Slide-mounted tumor specimens were classified by hematoxylin and eosin (H&E) staining according to the 2010 WHO classification of digestive system tumors. The included tumor tissue samples were histologically confirmed as adenocarcinoma by two molecular pathologists and shown to consist of more than 70% tumor cells.

This study was conducted in accordance with the Helsinki Declaration (as revised in 2013) (28), and the protocol was approved by the Ethics Committee of the First Affiliated Hospital of Sun Yat-sen University (No. [2017]174). Informed written consent was obtained from each patient before their inclusion in the study.

DNA extraction and qualification

Tumor tissue sections were subjected to DNA extraction using a DNeasy Blood and Tissue Kit (69504; Qiagen, Venlo, Netherlands) in line with the manufacturer’s instructions. The extracted DNA was characterized using an Agilent Bioanalyzer (Santa Clara, CA, USA). Germline mutations were determined by reference to matched leukocyte DNA.

Target gene sequencing and WES

A NadPrep DNA Library Preparation Kit (#1002101; Nanodigmbio, Nanjing, China) was used to prepare libraries for the Illumina sequencing platform, which is equipped with a unique dual index (UDI) adapter that can mitigate label-skipping and label-mismatch issues. A NadPrep Hybrid Capture Reagents kit (#1005101; Nanodigmbio) was used for hybridization capture–based target enrichment.

Variant annotation analysis and visual mapping

False-positive results were excluded through realigning split-reads to hg19/GRch37 using Blast. The computational efficiency of WES data capture was assessed using the Picard tool (http://broadinstitute.github.io/picard/command-line-overview.html#CollectHs-Metrics), which filters low-quality bases, repetitive bases, bases that deviate from the target, and ends that overlap with reads at both ends due to short inserts. Using this tool, the strictest depth distribution was calculated. For the purpose of quality control, sequencing data were considered to be qualified if the mean bait coverage exceeded 100. To improve the alignment of individual reads, further local rearrangements were performed using SpeedSeq (29). Somatic mutations were identified using MuTect (30). ANNOVAR and Oncotator were used to annotate genetic variant data, which were converted to a MAF (Mutation Annotation Format) file using the maf tool (31). Cancer driver genes were analyzed using the Integrative OncoGenomics (IntOgen) platform, running Oncodrive FM and Oncodrive CLUST (32). The landscape of the top driver mutation spectrum for tumors predicted by IntOgen, including mutation rate and mutation subclass/subtypes (Oncodrive FM P≤0.1), was visualized with R Script.

Pathway enrichment analysis

Differences in the frequency of somatic gene mutations between the thymoma and thymic carcinoma groups were illustrated using Wayne diagrams, as some gene mutations occurred only in one group, while others occurred in both groups. WebGestalt (WEB-based Gene SeT AnaLysis Toolkit, http://www.webgestalt.org) was used to analyzed the canonical pathways associated with genes containing single-nucleotide variations (SNVs). P values were calculated based on a hypergeometric distribution, and the Benjamini-Hochberg method was adopted for false discovery rate correction (33).

Statistical analysis

Correlations of clinical and biological variables were analyzed using the Statistical Package for the Social Sciences (SPSS) 22.0 (IBM Corp., Armonk, NY, USA). Where necessary, non-parametric tests, such as Welch’s t-test and Wilcoxon signed-rank test, in the R package ggpubr were used (34).

Results

Patient characteristics

To better understand the molecular mechanisms underlying the genesis of TETs, fresh frozen tumor tissues and corresponding blood samples were selected from 15 of the 21 patients with TETs in the study and were subjected to bioinformatics analyses according to strict standards. Clinicopathological assessment indicated that there were 9 cases (S01, S04, S05, S06, S10, S12, S13, S14, and S15) of thymic carcinoma, of which one case (S06) also had a lung metastasis. The other 12 cases (S02, S03, S07, S08, S09, S11, S16, S17, S18, S19, S20, and S21) were thymomas. Cases S16, S17, S18, S19, S20, and S21 had a tumor cell content of less than 50%; WES was not performed on these samples. Detailed clinical characteristics of the patients are shown in . The PD-L1 positive rate was 33.3% (3/9) and 58.3% (7/12) in the thymic carcinoma and thymoma groups, respectively ( and ). The detailed clinical characteristics of the patients are shown in Table S1.

Table 1

Patient characteristics

Characteristic	No. of cases
Total number	n=21
Age, years, mean [range]	51 [21–69]
Sex, n (%)
Male	15 (71.4)
Female	6 (28.6)
Pathological type, n (%)
Thymic carcinoma	9 (42.9)
Thymoma	12 (57.1)
PD-L1 positivity, n (%)
Thymic carcinoma	3 (33.3)
Thymoma	7 (58.3)

Figure 1

Pictures of hematoxylin and eosin staining and programmed death-ligand 1 (PD-L1) detection of thymic epithelial tumor samples (all images were ×400 magnification fields).

Identifying somatic mutations in 15 patients with TETs

NGS sequencing of DNA from 15 primary tumor tissues and matched blood samples was performed, and several somatic mutations were annotated using MuTect and Somatic Indel Detector. The average depth of WES of the tumor samples was 346.06×, and the average depth of sequencing of the control samples was 127.96×. A total of 1,779 nonsynonymous SNVs were detected, including 912 in thymic carcinomas and 867 in thymomas (https://cdn.amegroups.cn/static/public/ATM-21-5182-1.xlsx).

As shown in , case S08 had the most SNVs, followed by case S03. The top 30 driver genes were ranked based on the frequency of their somatic mutation, which showed that ZNF429 (36%) had the highest mutation frequency. The most common mutation type was missense mutations, followed by nonsense, splice site, and other mutations. Among the top 30 driver genes, BAP1 (14%), ABI1 (7%), BCL9L (7%), and CHEK2 (7%) mutations were detected exclusively in thymic carcinoma, while ZNF721 (14%) and PABPC1 (14%) mutations were found exclusively in thymoma. ZNF429, ZNF208, BAP1, ERBB4, GNAQ, ABI1, BCL9L were common somatic mutations in thymic carcinoma subgroup compared with healthy people.

Figure 2

Somatic mutations in thymic epithelial tumor (TET) samples. Top: variations in the numbers of mutations in each sample; middle: mutated genes and mutation types; right: frequency of mutations in somatic genes in TET samples.

The TMB was calculated using only somatic nonsynonymous mutations. The mean TMB values were 0.722 and 0.663 mutations per Mb in the thymic carcinoma and thymoma groups, respectively, and the difference between the groups was not statistically significant (P=0.817, Welch’s t-test; ). Notably, the TMB value in this study (0.698) was significantly higher than that reported in The Cancer Genome Atlas thymoma (TCGA-THYM) (0.488) (P=0.0038, Wilcoxon signed-rank test; ).

Figure 3

Comparative analysis of differences in tumor mutational burden (TMB). (A) TMB in thymic carcinoma and thymoma samples. (B) Comparison of TMB between the data in this study and The Cancer Genome Atlas thymoma (TCGA-THYM) database. (C) Comparative analysis of the relative transition/transversion (TI/TV) ratio between this study and the TCGA-THYM study. (D) Transition and transversion proportions of six nucleotide changes. The stacked proportion bar chart is sorted according to the increase in TI/TV fraction.

In principle, all mutation types (e.g., substitutions, indels, and rearrangements) and any accessory mutation characteristic (e.g., the sequence context of the mutation or the transcriptional chain where the mutation occurred) can be incorporated into the set of features by which a mutational signature is defined. Mutational signatures were extracted using base substitutions, and six classes of substitutions were assigned according to the pyrimidine of the mutated Watson-Crick base pair: C > A, C > G, C > T, T > C, T > A, and T > G. The six mutation types were then compared with the TCGA-THYM database, and the proportions of these six mutation replacement types were found to be roughly the same. Among all the substitutions, C > T had the highest mutation percentage, which differed significantly from that for the TCGA-THYM database (P=0.0014). The T > G, C > A, and T > A substitutions in this study also differed significantly from those in the TCGA-THYM study ().

Comparison of signaling pathway alterations between thymoma and thymic carcinoma

For further characterization the functions of mutational genes and pathways involved in TET, the PANTHER classification system, an ontology-based pathway database, was coupled with data analysis tools (35). Using pathological diagnosis results, a Venn diagram was created, which divided the driver mutant genes detected by WES into three clusters: thymic carcinoma (38 genes), thymoma (38 genes), and the intersection of the two (15 genes) ().

Figure 4

Analysis of pathways involving nonsynonymous genes in thymic epithelial tumors. (A) Venn diagram showing gene distribution in three clusters. (B) Volcano plot of uniquely mutated driver gene enrichment pathways in thymic carcinoma. (C) Volcano plot of uniquely mutated driver gene enrichment pathways in thymoma. (D) Volcano plot of shared driver gene enrichment pathways in the intersection group.

In thymic carcinoma, mutant somatic genes were mainly enriched in chronic myeloid leukemia (hsa05220), the ErbB signaling pathway (hsa04012), and the T-cell receptor signaling pathway (hsa04660) (). Meanwhile, altered somatic genes in thymoma were mainly enriched in central carbon metabolism in cancer (hsa05230), acute myeloid leukemia (hsa05221), MAPK signaling pathways (hsa04010), and longevity-regulating pathways (hsa04213) (). The shared driver gene enrichment pathways at the intersection of the thymic carcinoma and thymoma clusters mainly included lysine degradation (hsa00310), endocrine and other factor-regulated calcium reabsorption (hsa04961), long-term depression (hsa04730), long-term potentiation (hsa04720), gastric acid secretion (hsa04971), and the ErbB signaling pathway (hsa04012) ().

Discussion

TET, a rare malignant tumor of the chest, is the most common anterior mediastinal neoplasm in adults, with an incidence rate of 1.5–3.2 cases per million people per year (36). Among TETs, thymomas are often uniquely associated with autoimmune diseases, whereas thymic carcinomas exhibit more invasive characteristics in clinic. A study by Radovich et al. reported a strong association of broad histological subtypes (A, AB, B, and thymic carcinoma) with multiple classes of aberrations occurring at different levels (18). The etiology of TETs is unknown, largely due to the limited knowledge of the genomic underpinnings of thymoma and thymic carcinoma.

In this study, WES of 15 tumor-normal paired samples revealed some high-frequency mutation driver genes that have been reported in previous studies on TETs, including BAP1, ERBB4, HRAS, and KIT (37,38). BAP1 and BCL9L, as oncogenes, promote the growth, invasion and lung metastasis of breast cancer by interacting with β-catenin and Pygopus (39-41). ABI1 inhibits epithelial-mesenchymal transformation by inhibiting atypical WNT signal activation, providing a new mechanism for the inhibition of prostate tumors (42). Analysis further showed that the ZNF429 (36%) and ZNF208 (29%) genes had the highest frequency of mutations (). The ZNF429 gene was reported to be significantly associated with acute postoperative pain in patients with breast cancer, while ZNF208 polymorphisms were found to be associated with ischemic stroke in a southern Chinese Han population (43,44). Research has found that PABPC1 [poly(A) binding protein cytoplasmic 1] is directly involved in the tumorigenesis of gastric cancer through its promotion of cell proliferation; consistent with this, downregulated PABPC1 expression is associated with tumor progression and a poor prognosis in colorectal cancer (45,46). Our analysis of the TCGA-THYM dataset further revealed that the genes with the highest mutation frequency were GTF2I (50%), HRAS (8%), and MUC16 (7%) (). Radovich et al. showed that mutated GTF2I is a thymoma-specific oncogene, and that cells expressing mutant GTF2I also highly express genes involved in the Wnt and Sonic hedgehog (SHH) signaling pathways (18). Perplexingly, the current study did not detect mutations in the GTF2I gene in thymoma or thymic carcinoma samples, which is possibly due to insufficient sample size. MUC16, a mucin expressed at high levels on the surface of epithelial ovarian tumor cells, exerts similar effects to MUC1, which binds or aggregates neutrophils and other cell types to provide antitumor immune protection (47,48).

Figure 5

Mutational landscape of somatic alterations in The Cancer Genome Atlas thymoma (TCGA-THYM) database.

PD-L1 expression has been used as a tumor cell marker in a number of clinical trials and has been approved by FDA (the food and drug administration) (49). Positive expression of PD-L1 in tumors is considered to indicate patients who are more likely to respond to immune checkpoint inhibitors treatment (50). However, some patients who test positive for PD-L1 may fail to respond to treatment, and—more importantly—some patients who test negative may respond, making this an imperfect biomarker (51). Multiple studies have reported conflicting data on the relationship between PD-L1 expression and survival/prognosis in patients with TETs (19,22,25,26). In the current study, the rates of PD-L1 positivity in thymoma and thymic carcinoma (33.3% and 58.3%, respectively) were consistent with the positivity-rate range reported in the literature (19-25).

High TMB is a characteristic associated with responsiveness to immunotherapy (52). However, a high TMB value (≥20 mutations per Mb) is not suitable for predicting the immunotherapeutic effect in all solid tumors (14). Establishing the high TMB threshold for different cancer types may require additional clinical studies and information statistics of a larger number of patients. The TMB value in this study was significantly higher than that in the TCGA-THYM dataset, which suggests that the patients in our study would be more likely to benefit from immunotherapy ().

Because of the differences between our findings regarding TMB and those reported by the TCGA-THYM dataset, we also analyzed the types of mutation, as defined by the TI/TV ratio. Significant differences existed between our study and the TCGA data with respect to C > A, C > T, T > A, and T > G mutations; however, no significant difference was found for C > G and T > C mutations (). Mutational signatures can be understood as different mutation processes that often generate different combinations of mutation types (53). Thousands of somatic mutations can be identified in a single cancer sample, thus making it possible to decipher the mutant signature, even if several mutations are operative (54). The C > A mutational signature has been shown to be associated with tobacco smoking, whereas C > T mutations are induced by ultraviolet light (55).

Single or multiple genes alternation and expression level change usually response to pathway change and complex cross-talk between different signaling pathways simultaneously, such as the angiogenesis and notch signaling pathways (56). Cancers are not driven by single gene mutations or expression changes but by coordinated changes affecting multiple signaling pathways (57,58). ErbB (erythroblastic oncogene B), which is also called the epidermal growth factor receptor (EGFR), is generally expressed on the surface of epithelial cells, but is often overexpressed in certain tumor cells (59). Overexpression of EGFR is associated with tumor cell metastasis and invasion, and a poor prognosis (60). In the current study, we found enrichment of the ErbB and T-cell signaling pathways in a high proportion of samples in the thymic carcinoma group and that were also found at the intersection of gene clusters (). In the thymoma group, prominent pathways included longevity-regulating pathways (hsa04213) and central carbon metabolism in cancer (hsa05230) (). Cell immortalization is an important stage in tumorigenesis (61). In mammals, calorie restriction (CR) is considered to be the most reliable and repeatable longevity-promoting strategy. Among the pathways involved in the regulation of the CR effect are the insulin-like growth factor (IGF-1)/insulin signaling pathway, the sirtuin pathway, the adenosine monophosphate (AMP)-activated protein kinase (AMPK) pathway. and the target of rapamycin (TOR) pathway (62). These pathways are believed to respond to CR via a common mechanism by which autophagy, stress defense mechanisms, and survival pathways are activated, and proinflammatory mediators and cell growth become weaker, thereby promoting cell health and, ultimately, longevity (63). Anti-PD-l/PD-Ll and other immunotherapeutic drugs have good efficacy and safety in the treatment of solid tumors such as melanoma, non-small cell lung cancer, and head and neck cancer (64). At present, some domestic studies suggested that the expression of PD-L1 was related to the pathological classification of thymus, especially in type B2/B3 thymoma and thymic carcinoma, which indicated that the expression of PD-L1 was related to the prognosis of thymic cancer (65). However, whether patients with thymic cancer can be determined to benefit from immunotherapy requires more clinical research. In short, the carcinogenesis of thymic epithelial cells may involve benign or malignant mutations of multiple genes or cross feedback of signaling pathways. In the treatment of thymic epithelial tumors, we must not only use surgical resection or adjuvant radiotherapy and chemotherapy based on the pathological staging of the patient, but also pay attention to the comprehensive treatment based on the molecular classification of the patient.

Conclusions

The present study has analyzed somatic gene mutations, TMB, and TI/TV ratios in thymoma and thymic carcinoma using bioinformatics, and has also expounded the differences between our findings and those of the TCGA-THYM dataset. Drive genes of ZNF429, BAP1, ABI1, BCL9L, CHEK2, ZNF721 and PABPC1 are expressed to varying degrees in thymoma or thymic carcinoma, respectively. Enrichment of the longevity-regulating, MAPK, and ErbB signaling pathways in thymoma still remains to be further investigated. An in-depth understanding of the differences in somatic genes and pathways between thymoma and thymic cancer will help to design personalized clinical therapy for TET patients.

The article’s supplementary files as