Clonal Evolution Dynamics in Primary and Metastatic Lesions of Pancreatic Neuroendocrine Neoplasms

Background: Data on inter-tumoral heterogeneity and clonal evolution of pancreatic neuroendocrine neoplasms (panNENs) with liver metastasis are limited. The aim of this study was to explore different patterns of clonal evolution of pancreatic neuroendocrine neoplasms with liver metastasis and the possible distinctive signaling pathways involved between G2 neuroendocrine tumors (NETs) and neuroendocrine carcinomas (NECs).

Methods: Tumor tissues of five patients (10 samples) with pancreatic neuroendocrine neoplasms with synchronous liver metastasis were analyzed using next-generation sequencing. PyClone, Gene Ontology, and Reactome pathway enrichment analysis were also applied.

Results: Mutated genes varied in individuals, reflecting the inter-tumoral heterogeneity of panNENs. The distribution of subclones varied during tumor metastasis, and different clonal evolution patterns were revealed between NETs and NECs. Gene Ontology and Reactome analyses revealed that in both NETs and NECs, signaling pathways and biological processes shared similarities and differences in the primary and metastatic lesions. In addition, the signaling pathway features were different between NETs and NECs. In the primary lesions, epigenetic changes and post-transcriptional modifications participated in NETs, while FGFR signaling, EGFR signaling, and NTRK2 signaling were largely involved in NECs. Although DNA repair and TP53 regulation were both involved in the metastatic lesions, most of the signaling pathways and biological processes disrupted by the mutated genes were different.

Conclusions: Our study revealed spatial inter-tumoral heterogeneity and temporal clonal evolution in PanNENs, providing potential therapeutic targets for further prospective clinical trials.

Introduction

Pancreatic neuroendocrine neoplasms (panNENs) are heterogeneous tumors with distinct clinical syndromes and malignant potential. Unfortunately, the worldwide incidence and prevalence of panNENs have recently increased (1, 2). A subset of primary panNENs are aggressive and have a potential to metastasize. The 2019 World Health Organization (WHO) classification divides gastroenteropancreatic NENs into six groups—well-differentiated NETs G1, well-differentiated NETs G2, well-differentiated NETs G3, poorly-differentiated NECs—according to mitotic count and Ki-67 proliferation index, MiNEN and tumor-like lesions. Increased mitotic rate and high Ki-67 index are associated with a more aggressive clinical course and poorer prognosis (3, 4). Many studies have shown that distant metastases are related to poor overall survival in PanNENs (5, 6); however, the potential mechanisms of metastasis are poorly investigated and remain unclear.

Tumors progress and metastasize under Darwinian selection and evolution (7), and diverse genetic alterations create intra- and inter-tumoral heterogeneity. Difficulties with treatment, including drug resistance due to inter-tumoral heterogeneity and clonal evolution dynamics, are common. Some studies have shown concordance in mutational status between primary and distant metastases, where others do not (8, 9). Recent studies of breast cancer lung metastasis in mice showed that metastatic lesions were polyclonal, while omental metastases of high grade serous human ovarian cancer predominantly exhibited a single phylogenetic clade (10, 11). However, there are few studies investigating the clonal variation between primary and metastatic lesions of panNENs, or the different metastatic mechanisms among different grades. Studies of inter-tumoral heterogeneity in primary and metastatic lesions of panNENs are of great importance and could lead to more accurate clinical strategies for treatment.

The advent of next-generation sequencing (NGS) has increased understanding of tumor heterogeneity and clonal evolution (12, 13). In our study, we used a 450 oncogene sequencing panel to sequence pairs of primary tumors and hepatic metastases from synchronously metastasized panNENs from five patients. PyClone was applied to analyze the clonal populations and demonstrate the clonal distribution variation of primary and metastatic lesions. In addition, GO (Gene Ontology) and pathway enrichment analysis revealed distinct biological pathways that participated in primary and metastatic lesions in G2 neuroendocrine tumors (NETs) and neuroendocrine carcinomas (NECs).

Materials and Methods

Patients

Five panNEN patients with synchronous liver metastasis were included. Fresh tumor specimens were collected at the time of synchronous resection of primary tumors and liver metastases. Primary tissue samples, metastatic tissue samples, and adjacent normal tissue samples were collected. Our research was permitted by the Ethics Committee of the First Affiliated Hospital of Zhejiang University and the patients were informed of and gave consent to the use of tumor tissues for this research. All methods were performed in accordance with the declaration of Helsinki. Patient characteristics including patient demographics, pathologic TNM staging, histology stage, tumor thrombus, functional status, surgical approach, disease free survival (DFS), and patient status were obtained. TNM staging was adopted according to the 8th AJCC cancer staging system for neuroendocrine tumors of the pancreas. Grade was adopted according to the new WHO 2019 grading classifications.

Hematoxylin and Eosin and Immunohistochemistry

Hematoxylin and eosin (HE) staining and immunohistochemistry against Ki67 were carried out on each sample to confirm diagnosis and to determine the histological subtype. Briefly, the tumor samples were fixed in 4% neutral paraformaldehyde, dehydrated via a gradient ethanol, and embedded within paraffin blocks. Then, the histological sections (5 um) were prepared, deparaffinized, hydrated, and subjected to HE staining and immunohistochemistry (IHC). For IHC, the hydrated sections were first heated in antigen retrieval buffer, permeabilized with PBS containing 0.1% Triton X-100, and then blocked with 1% bovine serum albumin (BSA). Then, the sections were incubated with a Ki67 primary antibody at 4°C overnight and treated with 3% H₂O₂ to inhibit endogenous peroxidases. Then, sections were incubated with an HRP-linked secondary antibody for 1 h at room temperature; 3,30-diaminobenzidine was used as a chromogenic agent. Finally, the nuclei were stained with hematoxylin before dehydration and mounting.

Hybrid Selection and Sequencing

A custom hybridization capture panel including over 23,660 individually synthesized 5′-biotinylated DNA 120 bp oligonucleotides was used to target ~2.6 Mb of the human genome, including most exons of cancer-related genes and select introns of genes frequently rearranged in cancer (a 450 oncogene sequencing panel was used; gene names are listed in Supplementary Table 1). Hybridization capture employed xGen^® Lockdown^® Probes and Reagents (Integrated DNA Technologies, Version 3). Post-capture libraries were mixed, denatured, diluted, and then sequenced. For estimation of sequencing error rate, a PhiX spike-in was added as an external control to measure the percentage of reads with 0–4 mismatches, following the method described by Manley et al. (14). The average sequencing depth was 1000X for tissue-based deep sequencing.

Bioinformatics Pipeline for SNV and Short Indels

Alignment of raw reads to the human genome reference sequence (hg19) was done with the Burrows-Wheeler Aligner (BWA, v0.6.2), followed by PCR duplicate removal using the MarkDuplicates algorithm from Picard (version 1.47, http://picard.sourceforge.net/). Local realignment and base quality recalibration for single nucleotide variants (SNV) were performed using GATK (v3.1-1) and subsequently culled by MUTECT (v1.7).

Bioinformatics Pipeline for Copy Number Alternations

To identify CNA, aligned reads were first normalized within each bed by EXCAVATOR (version v2.2, http://sourceforge.net/projects/excavatortool/). Log ratio of read depths for each gene from tumor tissue and its matched normal blood control was then calculated. Tumor cellularity was estimated by allele frequencies of 4,088 sequenced SNPs (single nucleotide polymorphism), following the method in ASCAT.

Bioinformatics Pipeline for Gene Rearrangement

For detection of gene rearrangement, aligned reads with abnormal insert size of over 2,000 or zero bp were collected and used as discordant reads, i.e., paired-end reads that could not be closely mapped to a genome reference, with each read of paired-reads aligned to the same chromosomes or different chromosomes. Originally, the discordant reads with the distance <500 bp formed clusters were further assembled by fermi-lite (https://github.com/lh3/fermi-lite) to identify potential rearrangement breakpoints. The breakpoints were double-confirmed by BLAT and the resulted chimeric gene candidates were annotated.

Tumor Mutation Burden and Microsatellite Instability

Tumor mutation burden (TMB) was estimated by counting somatic mutations, including coding SNVs and indels per megabase of the sequence examined. Driver mutations and known germline alterations in dbSNP were not counted. MSI status was inferred based on MANTIS (15) score, and microsatellite regions were manually reviewed in Integrated Genomics Viewer (IGV) (16) for confirmation.

PyClone

A Bayesian clustering method, PyClone, was used to infer clonal population structures present in the tumor as previously described (17). Briefly, given the mutation allele frequencies for each sample, PyClone estimates cellular prevalence for each cluster in each sample. PyClone is freely available for academic use at https://github.com/Roth-Lab/pyclone.

GO and Reactome Pathway Enrichment Analysis

Related genes were selected using the Cytoscape GeneMANIA plugin (18). Genes from NEC and NET samples were combined separately and enrichment analysis for Biological Process was performed using the R software clusterProfiler (19) package with the pvalueCutoff set as 0.05. Reactome pathway enrichment analysis was performed using the ReactomePA package (20) with the p-value cutoff set as 0.05. Only the top 20 entries with a minimum adjusted p-value for primary and metastatic sites were included in the dotplot.

Results

Patient Demographics

Five panNEN patients with synchronous liver metastasis were included. All patients received synchronous resection of primary and liver metastases. The pathology of patients 1 and 2 was NEC. In the primary and metastatic lesions, the cells were poorly differentiated, composed of highly atypical neoplastic cells, with a mitotic rate > 20/10 high-power fields, and ki67 > 20%. Patient 3, patient 4, and patient 5 were G2 NET. The cells were arranged in nests and glandular tubes, were trabecular, had a mitotic rate of 2–20/10 high-power fields, and ki67 3–20%. All tumors were positive for Synaptophysin and Chromogranin A. Representative hematoxylin- and eosin-stained sections and IHC against ki67 in primary tumors and metastases are shown in Figure 1. TNM stage is also listed. All tumors were non-functioning. The tumors of patients 1 and 2 exhibited perineural invasion, and patient 5 was confirmed to have tumor thrombus. Among the surgical choices for these patients, four patients were treated with distal pancreatectomy and one by pancreaticoduodenectomy. Liver tumor resection was performed in all cases. Mean follow-up time was 51.7 months (46.9–54.6 months). Patient characteristics are listed in Table 1.

FIGURE 1

Figure 1. Representative hematoxylin and eosin and immunohistochemistry against ki67 in primary and metastatic lesions. In patient 1 and patient 2, HE staining of primary (pancreas) and metastatic lesions (liver) showed a mitotic count >20/10 HPF; the cells were poorly differentiated and composed of highly atypical neoplastic cells. Ki67 of patient 1 was 30–40%; patient 2 was 25–30%. In patients 3, 4, and 5, HE staining showed a mitotic count of 2–20/10 HPF. The cells were arranged in nests, glandular tubes, and were trabecular. Ki67 of patient 3 was 10–20%, patient 4 was 5–10%, and patient 5 was 8–10%. Magnification: 100x.

TABLE 1

Table 1. Characteristics of the five patients.

Sequencing

All primary and metastatic lesions of the five patients were sequenced (Figure 2, Supplementary Table 2) using a 450 oncogene sequencing panel. In patient 1, EPHA2, KRAS, SMARCB1, APC, and SPINK1 were altered in both primary and metastatic lesions, whereas RNF43 and AXIN2 were altered only in metastatic lesions. The VAF of EPHA2, KRAS, SMARCB1, APC all arose in metastatic lesion. In patient 2, TP53 and MEN1 were altered in both samples. Alteration of PIK3CA, RARA, CARD11, TET1, and PRSS1 were found only in the primary lesion while PRSS8 and PLA2G1B were found only in the metastatic lesion. The VAF of TP53 and MEN1 decreased in metastatic lesion. In patient 3, ETV1 was altered in the primary lesion and FGFR3 was altered in the metastatic lesion. In patient 4, MEN1, TSC2, and SIK1 were altered in both lesions and the VAF of TSC2 and MEN1 decreased in metastatic lesion. In patient 5, ATRX and MEN1 were altered in the primary and metastatic lesions. SMAD4 and TSC2 were altered in the primary lesion while FANCM was altered in the metastatic lesion. The VAF of ATRX arose from 0.53 to 0.57 in metastatic lesion while VAF of MEN1 decreased from 0.42 to 0.38. Mutated genes varied in the individuals, reflecting the heterogeneity of NETs and NECs. Compared with NETs, more mutations were revealed in NECs (more than five alterations in both primary lesions and metastatic lesions). And in NETs, gene alterations were less than four. After metastasis, the VAF elevated in one NEC patient, while the increasing VAF were hardly seen in NETs. In addition, in each patient, most of the mutated genes in the primary and metastatic lesions were common (except for patient 3). The TMB value varied from 1.6 to 6.4 mut/Mb in the five patients. All patients were classified as microsatellite stable (MSS).

FIGURE 2

Figure 2. Mutation characteristics and TMB. Mutation details of primary (p) and metastatic (m) samples of patients 1–5. Green blocks represent substitution/indel. Red blocks represent gene amplification. Blue blocks represent homozygous gene deletion. Purple blocks represent truncation.

PyClone

Based on variations identified using Bayesian clustering with PyClone, we identified three clusters in patients 1 and 2 (Figures 3A,B), one cluster in patient 3, and two clusters in patients 4 and 5 (Figures 3C–E); one cluster represents one subclone. Inter-tumoral heterogeneity was observed between the primary tumor samples and the metastases. The distribution of subclones changed, reflecting the evolution of tumors from primary to metastatic lesions. All tumor samples except for those of patient 3 showed evidence of subclonal structure.

FIGURE 3

Figure 3. PyClone analysis of the five patients. PyClone quantification of clonal populations in sequenced mutations of patients 1–5 (A–E, respectively). Cellular prevalence of clonal populations is the mean value of the cellular prevalence of mutations in the cluster. The predicted cellular frequencies represent the proportion of cancer cells in each set of clonal mutation.

NEC patients had three clonal populations in the primary lesions, while two clonal populations were identified in patient 1 and only one clonal population in patient 2. In addition, fewer clonal populations were found in G2 NETs and only one clonal population was identified in the metastatic lesions. Different clonal evolution patterns were identified between NETs and NECs.

GO Term and Reactome Analysis

We performed Gene Ontology (GO) enrichment analysis and Reactome analysis on sets of genes mutated in NETs and NECs. Mutated genes in the primary lesions of NET participated in sets of epigenetic changes and post-transcriptional modifications, such as histone methylation, chromatin modification, and macromolecule and protein methylation. In addition, the Reactome analysis showed that chromatin organization, SUMOylation, and the RUNX1 and NOTCH signaling pathways were most involved in the primary lesion. In the metastatic lesions of NETs, however, both Reactome and GO analysis showed that the mutated genes were prone to participation in histone and chromatin modification, DNA structure changes, DNA repair, and TP53 regulation. Histone and chromatin modification, DNA repair, and the regulation of TP53 and RUNX1 were involved in the primary and metastatic lesions in NET, whereas most other pathways were different (Figures 4A,B).

FIGURE 4

Figure 4. GO and pathway enrichment analysis of mutated genes involved in G2 NETs. GO enrichment analysis of biological process in mutated genes in NETs (A). Reactome pathway enrichment analysis in NETs (B). Dot size corresponds to number of genes; gene ratio is defined as percentage of genes in certain pathways compared with all genes in the samples.

In the primary lesions of NECs, GO analysis revealed that the mutated genes were most related to radiation response, cellular response to peptides and hormones, and transcription factor activity regulation. The Reactome analysis found that FGFR signaling, EGFR signaling, and NTRK2 signaling were largely involved in the NEC primary lesions. In the metastatic lesions of NECs, GO analysis showed that the mutated genes largely participated in organ development and differentiation processes. Reactome analysis revealed that mutated genes influenced DNA repair, chromatin modification, SUMOylation, WNT signaling, and the regulation of transcription and activity of TP53 and RUNX1. The common pathways between the primary and metastatic lesions in NECs were rare except for radiation response (Figures 5A,B).

FIGURE 5

Figure 5. GO and pathway enrichment analysis of mutated genes involved in NECs. GO enrichment analysis of biological processes in mutated genes in NECs (A). Reactome pathway enrichment analysis in NEC (B). Dot size corresponds to number of genes; gene ratio is defined as the percentage of genes in certain pathways compared with all genes in the samples.

These data revealed that in both NETs and NECs, signaling pathways and biological processes shared similarities and differences in the primary and metastatic lesions. However, the signaling pathways and biological process patterns were different between NETs and NECs. In the primary lesion, epigenetic changes and post-transcriptional modification participated in NETs, while FGFR signaling, EGFR signaling, and NTRK2 signaling were largely involved in NECs. In the metastatic lesions, although DNA repair and TP53 regulation were both involved, most of the signaling pathways and biological process disrupted by the mutated genes were different between NETs and NECs.

Discussion

Clonal evolution is defined by tumor heterogeneity over both space and time (21). Clonal evolutionary processes have been observed in many tumor types including pancreatic cancer (22), leukemia (23), and renal-cell carcinoma (24). Several large-scale genomic studies have characterized panNET genomes, including hundreds of somatic mutations and copy number variations, and reported that significantly mutated genes such as MEN1, DAXX, ATRX contribute to the mutagenic processes (25–27). Most panNECs harbor TP53 and RB1 alterations and lack neuroendocrine-related genetic changes (28). However, few researchers have studied the relationship between primary and metastatic lesions of individual panNENs from the perspective of clonal evolution. Our study performed NGS on five pairs of tumor lesions, revealing a comprehensive analysis of the course of tumor genomic evolution from primary lesion to metastatic lesion.

Tumors progress under Darwinian evolution, in which genetic variation alters molecular phenotypes in individual cells (29). Consequently, tumors often consist of multiple genes and different cell populations. These populations, known as clones, undergo selection in response to different tumor microenvironments or therapeutic interventions (30). Identifying dynamic clonal population structures can aid in predicting metastatic potential and chemotherapeutic resistance (17). PyClone is a Bayesian clustering method for grouping sets of somatic mutations into clonal clusters (17). PyClone was performed in the current study and identified three mutant clonal populations of primary lesions of panNEC and two clonal populations of panNET. Through the clone analysis, the sub-clones in each patient were identified, and the distribution of sub-clones reflected the evolutionary process of the tumors. Some clonal populations diminished, and some clonal populations expanded, which is one of the features of clonal evolution. Interestingly, metastatic lesions in NET patients tended to be from a single clone, which suggests monoclonal seeding from pancreatic lesion to liver and suggests a distinct clonal evolution mode from NEC. More importantly, the result should be confirmed by multiple lesion biopsies and NGS. It has been reported that panNETs exhibit a lower mutation burden (25), which may account for the small number of clonal populations in panNETs.

The failure of therapy and drug resistance is partially caused by intratumor heterogeneity, which provides diverse genetic material under evolutionary selection. In these five patients, the mutated genes varied in the individuals. The data from our study support the notion that pancreatic neuroendocrine carcinomas are fundamentally discordant from neuroendocrine tumors. Apart from their clinicopathological features, including histologic architecture, hormone production, and malignancy, the genetic profile differed. In patient 2, with the shortest overall survival, 45.9 months, TP53 was mutated in both the primary and metastatic lesions. It has been reported that in panNEC, ~70% of the tumors harbor TP53 mutations (28, 31, 32). Alternatively, in well-differentiated neuroendocrine neoplasms of the pancreas, TP53 mutations rarely occur. TP53 mutations, have been consistently associated with poor prognosis in cancers (33). Consistent with the literature, deleterious TP53 mutations were uncommon in panNET in our study. KRAS mutation was observed in one NEC patient and mutations in KRAS codon 12 were independently associated with a worse survival vs. wild-type KRAS (34). Somatic mutations of MEN1 occur in 30–44% of panNETs (25, 35, 36). In our study, patients 4 and 5 carried MEN1 mutations and patient 5 had a mutation of ATRX. DAXX and ATRX are mutually exclusive inactivating mutations, and no tumor with a mutation in ATRX had a mutation in DAXX. It has been reported that mutations in the MEN1 and the DAXX/ATRX genes are associated with prolonged survival compared with patients with tumors that lack these mutations (35). However, Fei Yuan et al. (37) and Marinoni et al. (38) reported that mutations of DAXX/ATRX were associated with a shortened survival. TSC2 mutations, which inhibit the mTOR signaling pathway and are mutually exclusive of mutations of PTEN, were revealed in patients 4 and 5. It has been reported that mTOR pathway genes (PTEN and TSC2) did not predict the pNET patients' survival (37). In addition, we observed high rates of common mutations in the primary and metastatc samples in these patients, except for patient 3. Recent studies on malignancies such as head and neck squamous cell carcinoma, non-small cell lung cancer, endometrial cancer, and small cell lung cancer all observed high rates of common mutations in primary and metastatic lesions (39–41).

Clonal evolution has two patterns, linear evolution or branched evolution. The linear model states that tumor cells acquire mutations over time, and that the strongest tumor becomes dominant. The branching evolution model states that parallel tumor cell clones acquire different genomic mutations over time (42, 43). In our study, both NEC and NET patients shared similar mutations in the primary and metastatic lesions, indicating the tumors were homogeneous at different locations and revealing a probable linear evolution model in PanNENs. This result should be confirmed by multiple biopsy analysis.

The mutated genes in the primary and metastatic lesions of NETs and NECs were analyzed via GO and Reactome analysis to explore the affected biological processes and pathways. The results revealed that the involved biological process and pathways were largely different among the primary and metastatic lesions of NETs and NECs. In NEC primary lesions, the EGFR, FGFR, and NTRK3 pathways were also involved. EGFR, FGFR, and NTRK3 targeted therapy may therefore be new options for NEC patients; however, clinical trials are required to validate this hypothesis. Further, mutated genes in the primary lesions of NETs exhibited sets of epigenetic changes, raising the intriguing possibility of epigenetic inhibitor treatment. Interestingly, though TP53 mutation was not detected in both primary and metastatic lesions in NETs, the Reactome analysis demonstrated that transcriptional regulation of TP53 was involved in NET primary lesions, whereas both the transcriptional regulation and activity regulation of TP53 were involved in NET metastasis. These results indicated that TP53 dysregulation may also play an important role in NET carcinogenesis and metastasis. In addition, we found that the DNA repair pathway participated in both NET and NEC metastasis. It has been reported that failed DNA repair can lead to carcinogenesis and tumor genome instability (44). And defects in DNA repair pathways may indicate a potential vulnerability to DNA-damaging therapies such as platinum (45). Additionally, oxaliplatin-based chemotherapy also shows a relatively high response rate in NETs (46). These findings not only provide better understanding of the mechanisms of NEC and NET carcinogenesis and metastasis, partially explaining the different clonal evolution patterns, but also give new insight into potential therapeutic approaches.

Immune checkpoint inhibitors have recently gained the attention of oncologists. PD-1/PD-L1 expression, tumor mutation burden, and DNA mismatch repair deficiency (dMMR) have been demonstrated as three potential biomarkers for the use of immune checkpoint inhibitors (47–50). Our results showed low TMB and MSS in both primary and metastatic lesions, which is consistent with other reports showing that neuroendocrine tumors have a relatively low mutation burden compared with other tumors (51) and may not benefit from immunotherapy.

This is the first study comparing primary and metastatic lesions of panNENs within the same person using NGS. Our study showed the mutation variation of primary and metastatic lesions and revealed clonal evolution dynamics in panNENs. We further detected different signaling pathways involved in the clonal evolution of NETs and NECs. However, our study had some limitations. The sample size was limited mostly due to the low incidence of PanNENs, and further studies including whole genome sequencing and other experimental studies are needed to verify the results. Our study does shed light on potential novel approaches for predicting and treating panNEN patients with synchronous liver metastasis.

Conclusions

Our study revealed spatial inter-tumoral heterogeneity and temporal clonal evolution in PanNENs, which provides potential therapeutic targets for further prospective clinical trials.

Data Availability Statement

According to national legislation/guidelines, specifically the Administrative Regulations of the People's Republic of China on Human Genetic Resources (http://www.gov.cn/zhengce/content/2019-06/10/content_5398829.htm, http://english.www.gov.cn/policies/latest_releases/2019/06/10/content_281476708945462.htm), no additional raw data is available at this time. Data of this project can be accessed after an approval application to the China National Genebank (CNGB, https://db.cngb.org/cnsa/). Please refer to https://db.cngb.org/, or email: Q05HQmRiJiN4MDAwNDA7Y25nYi5vcmc= for detailed application guidance. The accession code CNP0001685 should be included in the application.

Ethics Statement

The studies involving human participants were reviewed and approved by Research Ethics Committee of the First Affiliated Hospital, College of Medicine, Zhejiang University (2017-778). The patients/participants provided their written informed consent to participate in this study.

Author Contributions

ZT wrote the manuscript. LW, WS, and YZe analyzed the data. HZ and LL collected the clinical and pathological information from the cancer patients. YZh and CC designed the study. WF, WX, and PZ revised the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by National Natural Science Foundation of China (81472346), Natural Science Foundation of Zhejiang Province (LY15H030012). The funding source had no role in the design of the study, data collection, data analysis, or manuscript writing.

Conflict of Interest

WS and YZe were employed by the company OrigiMed.

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

The reviewer H-KW declared a shared affiliation, though no other collaboration, with several of the authors, ZT, LW, HZ, LL, YZh, CC, WX, WF, and PZ to the handling editor.

Supplementary Material

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

References

1. Yao JC, Hassan M, Phan A, Dagohoy C, Leary C, Mares JE, et al. One hundred years after “carcinoid”: epidemiology of and prognostic factors for neuroendocrine tumors in 35,825 cases in the United States. J Clin Oncol. (2008) 26:3063–72. doi: 10.1200/JCO.2007.15.4377

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Metz DC, Jensen RT. Gastrointestinal neuroendocrine tumors: pancreatic endocrine tumors. Gastroenterology. (2008) 135:1469–92. doi: 10.1053/j.gastro.2008.05.047

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Panzuto F, Boninsegna L, Fazio N, Campana D, Pia Brizzi M, Capurso G, et al. Metastatic and locally advanced pancreatic endocrine carcinomas: analysis of factors associated with disease progression. J Clin Oncol. (2011) 29:2372–7. doi: 10.1200/JCO.2010.33.0688

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Pape UF, Jann H, Muller-Nordhorn J, Bockelbrink A, Berndt U, Willich SN, et al. Prognostic relevance of a novel TNM classification system for upper gastroenteropancreatic neuroendocrine tumors. Cancer. (2008) 113:256–65. doi: 10.1002/cncr.23549

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Birnbaum DJ, Turrini O, Ewald J, Barbier L, Autret A, Hardwigsen J, et al. Pancreatic neuroendocrine tumor: a multivariate analysis of factors influencing survival. Eur J Surg Oncol. (2014) 40:1564–71. doi: 10.1016/j.ejso.2014.06.004

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Li Z, Du S, Feng W, Zhang W, Li G, Wei J, et al. Competing risks and cause-specific mortality in patients with pancreatic neuroendocrine tumors. Eur J Gastroenterol Hepatol. (2019) 31:749–55. doi: 10.1097/MEG.0000000000001350

PubMed Abstract | CrossRef Full Text | Google Scholar

7. McGranahan N, Swanton C. Clonal heterogeneity and tumor evolution: past, present, and the future. Cell. (2017) 168:613–28. doi: 10.1016/j.cell.2017.01.018

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Brannon AR, Vakiani E, Sylvester BE, Scott SN, McDermott G, Shah RH, et al. Comparative sequencing analysis reveals high genomic concordance between matched primary and metastatic colorectal cancer lesions. Genome Biol. (2014) 15:454. doi: 10.1186/s13059-014-0454-7

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Sebagh M, Allard MA, Bosselut N, Dao M, Vibert E, Lewin M, et al. Evidence of intermetastatic heterogeneity for pathological response and genetic mutations within colorectal liver metastases following preoperative chemotherapy. Oncotarget. (2016) 7:21591–600. doi: 10.18632/oncotarget.7809

PubMed Abstract | CrossRef Full Text | Google Scholar

10. McPherson A, Roth A, Laks E, Masud T, Bashashati A, Zhang AW, et al. Divergent modes of clonal spread and intraperitoneal mixing in high-grade serous ovarian cancer. Nat Genet. (2016) 48:758–67. doi: 10.1038/ng.3573

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Cheung KJ, Padmanaban V, Silvestri V, Schipper K, Cohen JD, Fairchild AN, et al. Polyclonal breast cancer metastases arise from collective dissemination of keratin 14-expressing tumor cell clusters. Proc Natl Acad Sci USA. (2016) 113:E854–63. doi: 10.1073/pnas.1508541113

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Zhao B, Hemann MT, Lauffenburger DA. Modeling tumor clonal evolution for drug combinations design. Trends Cancer. (2016) 2:144–58. doi: 10.1016/j.trecan.2016.02.001

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Ding L, Raphael BJ, Chen F, Wendl MC. Advances for studying clonal evolution in cancer. Cancer Lett. (2013) 340:212–9. doi: 10.1016/j.canlet.2012.12.028

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Manley LJ, Ma D, Levine SS. Monitoring error rates in illumina sequencing. J Biomol Tech. (2016) 27:125–8. doi: 10.7171/jbt.16-2704-002

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Kautto EA, Bonneville R, Miya J, Yu L, Krook MA, Reeser JW, et al. Performance evaluation for rapid detection of pan-cancer microsatellite instability with MANTIS. Oncotarget. (2017) 8:7452. doi: 10.18632/oncotarget.13918

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Robinson JT, Thorvaldsdóttir H, Winckler W, Guttman M, Lander ES, Getz G, et al. Integrative genomics viewer. Nat Biotechnol. (2011) 29:24–6. doi: 10.1038/nbt.1754

CrossRef Full Text | Google Scholar

17. Roth A, Khattra J, Yap D, Wan A, Laks E, Biele J, et al. PyClone: statistical inference of clonal population structure in cancer. Nat Methods. (2014) 11:396–8. doi: 10.1038/nmeth.2883

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Montojo J, Zuberi K, Rodriguez H, Bader GD, Morris Q. GeneMANIA: fast gene network construction and function prediction for Cytoscape. F1000Res. (2014) 3:153. doi: 10.12688/f1000research.4572.1

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. (2012) 16:284–7. doi: 10.1089/omi.2011.0118

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Yu G, He QY. ReactomePA: an R/Bioconductor package for reactome pathway analysis and visualization. Mol Biosyst. (2016) 12:477–9. doi: 10.1039/C5MB00663E

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Amirouchene-Angelozzi N, Swanton C, Bardelli A. Tumor evolution as a therapeutic target. Cancer Discov. (2017) 7:805–17. doi: 10.1158/2159-8290.CD-17-0343

CrossRef Full Text | Google Scholar

22. Yachida S, Jones S, Bozic I, Antal T, Leary R, Fu B, et al. Distant metastasis occurs late during the genetic evolution of pancreatic cancer. Nature. (2010) 467:1114–7. doi: 10.1038/nature09515

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Anderson K, Lutz C, van Delft FW, Bateman CM, Guo Y, Colman SM, et al. Genetic variegation of clonal architecture and propagating cells in leukaemia. Nature. (2011) 469:356–61. doi: 10.1038/nature09650

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Gerlinger M, Rowan AJ, Horswell S, Math M, Larkin J, Endesfelder D, et al. Intratumor heterogeneity and branched evolution revealed by multiregion sequencing. N Engl J Med. (2012) 366:883–92. doi: 10.1056/NEJMoa1113205

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Scarpa A, Chang DK, Nones K, Corbo V, Patch AM, Bailey P, et al. Whole-genome landscape of pancreatic neuroendocrine tumours. Nature. (2017) 543:65–71. doi: 10.1038/nature21063

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Pea A, Yu J, Marchionni L, Noe M, Luchini C, Pulvirenti A, et al. Genetic analysis of small well-differentiated pancreatic neuroendocrine tumors identifies subgroups with differing risks of liver metastases. Ann Surg. (2020) 271:566–73. doi: 10.1097/SLA.0000000000003022

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Chan CS, Laddha SV, Lewis PW, Koletsky MS, Robzyk K, Da Silva E, et al. ATRX, DAXX or MEN1 mutant pancreatic neuroendocrine tumors are a distinct alpha-cell signature subgroup. Nat Commun. (2018) 9:4158. doi: 10.1038/s41467-018-06498-2

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Konukiewitz B, Jesinghaus M, Steiger K, Schlitter AM, Kasajima A, Sipos B, et al. Pancreatic neuroendocrine carcinomas reveal a closer relationship to ductal adenocarcinomas than to neuroendocrine tumors G3. Hum Pathol. (2018) 77:70–9. doi: 10.1016/j.humpath.2018.03.018

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Gonzalgo ML, Jones PA. Mutagenic and epigenetic effects of DNA methylation. Mutat Res. (1997) 386:107–18. doi: 10.1016/S1383-5742(96)00047-6

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Aparicio S, Caldas C. The implications of clonal genome evolution for cancer medicine. N Engl J Med. (2013) 368:842–51. doi: 10.1056/NEJMra1204892

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Konukiewitz B, Schlitter AM, Jesinghaus M, Pfister D, Steiger K, Segler A, et al. Somatostatin receptor expression related to TP53 and RB1 alterations in pancreatic and extrapancreatic neuroendocrine neoplasms with a Ki67-index above 20. Modern Pathol. (2017) 30:587–98. doi: 10.1038/modpathol.2016.217

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Yachida S, Vakiani E, White CM, Zhong Y, Saunders T, Morgan R, et al. Small cell and large cell neuroendocrine carcinomas of the pancreas are genetically similar and distinct from well-differentiated pancreatic neuroendocrine tumors. Am J Surg Pathol. (2012) 36:173–84. doi: 10.1097/PAS.0b013e3182417d36

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Olivier M, Hollstein M, Hainaut P. TP53 mutations in human cancers: origins, consequences, and clinical use. Cold Spring Harbor Perspect Biol. (2010) 2:a001008. doi: 10.1101/cshperspect.a001008

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Hayama T, Hashiguchi Y, Okamoto K, Okada Y, Ono K, Shimada R, et al. G12V and G12C mutations in the gene KRAS are associated with a poorer prognosis in primary colorectal cancer. Int J Colorectal Dis. (2019) 34:1491–6. doi: 10.1007/s00384-019-03344-9

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Jiao Y, Shi C, Edil BH, de Wilde RF, Klimstra DS, Maitra A, et al. DAXX/ATRX, MEN1, and mTOR pathway genes are frequently altered in pancreatic neuroendocrine tumors. Science. (2011) 331:1199–203. doi: 10.1126/science.1200609

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Corbo V, Dalai I, Scardoni M, Barbi S, Beghelli S, Bersani S, et al. MEN1 in pancreatic endocrine tumors: analysis of gene and protein status in 169 sporadic neoplasms reveals alterations in the vast majority of cases. Endocr Relat Cancer. (2010) 17:771–83. doi: 10.1677/ERC-10-0028

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Yuan F, Shi M, Ji J, Shi H, Zhou C, Yu Y, et al. KRAS and DAXX/ATRX gene mutations are correlated with the clinicopathological features, advanced diseases, and poor prognosis in Chinese patients with pancreatic neuroendocrine tumors. Int J Biol Sci. (2014) 10:957–65. doi: 10.7150/ijbs.9773

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Marinoni I, Kurrer AS, Vassella E, Dettmer M, Rudolph T, Banz V, et al. Loss of DAXX and ATRX are associated with chromosome instability and reduced survival of patients with pancreatic neuroendocrine tumors. Gastroenterology. (2014) 146:453–60.e455. doi: 10.1053/j.gastro.2013.10.020

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Saber A, Hiltermann TJN, Kok K, Terpstra MM, de Lange K, Timens W, et al. Mutation patterns in small cell and non-small cell lung cancer patients suggest a different level of heterogeneity between primary and metastatic tumors. Carcinogenesis. (2017) 38:144–51. doi: 10.1093/carcin/bgw128

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Gibson WJ, Hoivik EA, Halle MK, Taylor-Weiner A, Cherniack AD, Berg A, et al. The genomic landscape and evolution of endometrial carcinoma progression and abdominopelvic metastasis. Nat Genet. (2016) 48:848–55. doi: 10.1038/ng.3602

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Hedberg ML, Goh G, Chiosea SI, Bauman JE, Freilino ML, Zeng Y, et al. Genetic landscape of metastatic and recurrent head and neck squamous cell carcinoma. J Clin Investig. (2016) 126:1606. doi: 10.1172/JCI86862

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Polyak K. Is breast tumor progression really linear? Clin Cancer Res. (2008) 14:339–341. doi: 10.1158/1078-0432.CCR-07-2188

CrossRef Full Text | Google Scholar

43. Yates LR, Campbell PJ. Evolution of the cancer genome. Nat Rev Genet. (2012) 13:795–806. doi: 10.1038/nrg3317

CrossRef Full Text | Google Scholar

44. Bouwman P, Jonkers J. The effects of deregulated DNA damage signalling on cancer chemotherapy response and resistance. Nat Rev Cancer. (2012) 12:587–98. doi: 10.1038/nrc3342

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Bever KM, Le DT. DNA repair defects and implications for immunotherapy. J Clin Investig. (2018) 128:4236–42. doi: 10.1172/JCI122010

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Kunz PL, Balise RR, Fehrenbacher L, Pan M, Venook AP, Fisher GA, et al. Oxaliplatin-fluoropyrimidine chemotherapy plus bevacizumab in advanced neuroendocrine tumors: an analysis of 2 phase II trials. Pancreas. (2016) 45:1394–400. doi: 10.1097/MPA.0000000000000659

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Rosenberg JE, Hoffman-Censits J, Powles T, van der Heijden MS, Balar AV, Necchi A, et al. Atezolizumab in patients with locally advanced and metastatic urothelial carcinoma who have progressed following treatment with platinum-based chemotherapy: a single-arm, multicentre, phase 2 trial. Lancet. (2016) 387:1909–20. doi: 10.1016/S0140-6736(16)00561-4

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Le DT, Uram JN, Wang H, Bartlett BR, Kemberling H, Eyring AD, et al. PD-1 blockade in tumors with mismatch-repair deficiency. N Engl J Med. (2015) 372:2509–20. doi: 10.1056/NEJMoa1500596

CrossRef Full Text | Google Scholar

49. Hodges TR, Ott M, Xiu J, Gatalica Z, Swensen J, Zhou S, et al. Mutational burden, immune checkpoint expression, and mismatch repair in glioma: implications for immune checkpoint immunotherapy. Neurooncology. (2017) 19:1047–57. doi: 10.1093/neuonc/nox026

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Teng F, Meng X, Kong L, Yu J. Progress and challenges of predictive biomarkers of anti PD-1/PD-L1 immunotherapy: a systematic review. Cancer Lett. (2018) 414:166–73. doi: 10.1016/j.canlet.2017.11.014

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Cunha LL, Marcello MA, Rocha-Santos V, Ward LS. Immunotherapy against endocrine malignancies: immune checkpoint inhibitors lead the way. Endocr Relat Cancer. (2017) 24:T261–81. doi: 10.1530/ERC-17-0222

PubMed Abstract | CrossRef Full Text | Google Scholar