Genomic characterization of a well-differentiated grade 3 pancreatic neuroendocrine tumor.

Pancreatic neuroendocrine neoplasms (PanNENs) represent a minority of pancreatic neoplasms that exhibit variability in prognosis. Ongoing mutational analyses of PanNENs have found recurrent abnormalities in chromatin remodeling genes (e.g., DAXX and ATRX), and mTOR pathway genes (e.g., TSC2, PTEN PIK3CA, and MEN1), some of which have relevance to patients with related familial syndromes. Most recently, grade 3 PanNENs have been divided into two groups based on differentiation, creating a new group of well-differentiated grade 3 neuroendocrine tumors (PanNETs) that have had a limited whole-genome level characterization to date. In a patient with a metastatic well-differentiated grade 3 PanNET, our study utilized whole-genome sequencing of liver metastases for the comparative analysis and detection of single-nucleotide variants, insertions and deletions, structural variants, and copy-number variants, with their biologic relevance confirmed by RNA sequencing. We found that this tumor most notably exhibited a TSC1-disrupting fusion, showed a novel CHD7-BEND2 fusion, and lacked any somatic variants in ATRX, DAXX, and MEN1.

INTRODUCTION

Pancreatic neuroendocrine neoplasms (PanNENs) are neoplasms of the endocrine pancreas. These tumors account for 3%–5% of all pancreatic malignancies and in general have a better prognosis than typical pancreatic exocrine malignancies. Positive outcomes are largely dependent on complete surgical excision, as 5-yr survival approaches 55% when tumors are localized and resected, in contrast to 15% when they are not resectable (Edge et al. 2010; Dasari et al. 2017). Although PanNENs can closely approximate the morphology of neuroendocrine tumors of the gastrointestinal tract, their biologies differ, and it is, therefore, clinically advantageous to investigate alternative diagnostic and treatment approaches (Kulke et al. 2011).

PanNEN classification uses tumor morphology and grade to partition well-differentiated from poorly differentiated entities (Klimstra et al. 2010). Poorly differentiated pancreatic neuroendocrine carcinomas (PanNEC) are characterized by both a high mitotic count (>20 per 10 high power field) and a high proliferation index (Ki-67 > 20%). Well-differentiated pancreatic neuroendocrine tumors (PanNETs) tend to have lower mitotic rates and proliferation indices, classified into low-grade (G1) or intermediate-grade (G2) on the basis of their Ki-67 staining profiles (<3% and between 3% and 20% of tumor cells, respectively). It is also now recognized that a less-common subset of morphologically well-differentiated PanNETs has an incongruously high proliferation index (Ki-67 > 20%). Despite their high proliferative rate, these tumors appear to be separate biologic entities from the poorly differentiated NECs, so are now classified as high-grade (G3) PanNETs (Coriat et al. 2016; Tang et al. 2016; Lloyd et al. 2017).

Recurring molecular aberrations in PanNETs have been described in studies using whole-exome and whole-genome sequencing. Namely, mutations in DAXX and ATRX (chromatin remodeling genes) and TSC2, PTEN, PIK3CA, and MEN1 (mTOR pathway genes) have been observed (Jiao et al. 2011; Scarpa et al. 2017; Wong et al. 2018). Moreover, individuals exhibiting germline mutations in TSC1/2 (tuberous sclerosis), MEN1 (multiple endocrine neoplasia 1), VHL (von Hippel–Lindau), and NF1 (neurofibromatosis type 1) are at increased risk of developing PanNETs (Jensen et al. 2008).

Scarpa et al. (2017) provide a large-scale whole-genome landscape of PanNETs, which included five G3 PanNETs in their analysis, but the literature to date has an otherwise limited characterization of molecular aberrancy in these tumors. To this end, this report represents a complete genomic and transcriptomic characterization of a patient with a metastatic well-differentiated G3 PanNET.

RESULTS

Clinical Presentation

A previously healthy 35-yr-old male presented with weight loss and pain. Imaging demonstrated a mass in the pancreatic tail and liver metastasis. The patient was referred to the Personalized OncoGenomics Program (POG) for whole-genome and transcriptome analysis (NCT02155621) (Laskin et al. 2015). A core needle biopsy of one of the liver lesions revealed a G3 well-differentiated neuroendocrine tumor (Ki-67 = 30%; chromogranin- and CK7-positive; CK20-negative; octreotide scan–negative) (Fig. 1).

Figure 1.

Tumor histopathology and Ki-67 immunohistochemistry. (A) Core biopsy of a liver metastasis stained for hematoxylin and eosin showing a metastatic well-differentiated neuroendocrine tumor clinicoradiologically consistent with pancreatic origin (magnification at 10× [left] and 40× [right]). (B) Paired Ki-67 immunohistochemistry showing a proliferative index of 30%, classifying the tumor as a well-differentiated grade 3 (G3) neuroendocrine tumor (magnification at 10× [left] and 40× [right]).

Genomic Analysis

Whole-genome sequencing of liver metastases and blood was performed for comparative analysis and detection of somatic single-nucleotide variants (SNVs), insertion and deletions (indels), structural variants, and copy-number variants. RNA sequencing was performed to support variants detected by genomic sequencing and identify altered gene expression. Tumor content, as based on the tumor/normal sequencing ratio and loss of heterozygosity, was estimated to be 95% and a diploid model was used to describe the observed copy-number changes.

Single-Nucleotide Variants

Comparative analysis of whole-genome sequencing from tumor and normal samples revealed 4232 somatic SNVs. Somatic SNVs were compared to a catalog of previously described mutation signatures (http://cancer.sanger.ac.uk/cosmic/signatures) using the signIT package (https://github.com/eyzhao/SignIT/; EY Zhao, E Pleasance, and M Jones, in prep.). COSMIC mutation signature 5 was the most dominant mutation signature (Supplemental Table S1; Supplemental Fig. S1), with 30.8% (95% interval 17.7%–43.4%) of somatic SNVs attributed to this signature. Signature 5 was also detected as a dominant signature in the majority of ICGC PanNETs indicating a common etiology of mutagenesis across PanNETs (Scarpa et al. 2017). There were 15 somatic nonsynonymous mutations affecting protein-coding genes detected in the tumor biopsy, seven of which were subclonal (Table 1). Despite being recurrently mutated in PanNETs, we did not identify mutations affecting DAXX, ATRX, or MEN1 in this G3 PanNET sample. Similarly, none of the five G3 PanNETs from the ICGC study harbored mutations in DAXX, ATRX, or MEN1.

Table 1.

Somatic nonsynonymous single-nucleotide variants (SNVs) and indels

Gene	Chr	Position	Ref	Alt	HGVS codons	HGVS protein	Type	Genotype	Predicted Effect	dbSNP
HSPG2	1	22149892	G	A	c.13093C>T	p.Arg4365*	SNV	Heterozygous	Stop gain	–
CNR2	1	24202047	G	T	c.61C>A	p.Pro21Thr	SNV	Heterozygous	Missense	–
DNAJC13	3	132182654	A	T	c.1885A>T	p.Thr629Ser	SNV	Heterozygous	Missense	–
MUC7	4	71347240	T	C	c.779T>C	p.Val260Ala	SNV	Subclonal	Missense	rs145745951
RAPGEF2	4	160279276	A	C	c.4485A>C	p.Gln1495His	SNV	Heterozygous	Missense	–
FRG1	4	190876269	T	C	c.395T>C	p.Ile132Thr	SNV	Subclonal	Missense	rs149113074
DPCR1	6	30918448	T	C	c.2207T>C	p.Phe736Ser	SNV	Subclonal	Missense	–
PTPN6	12	7069329	C	T	c.1514C>T	p.Ala505Val	SNV	Subclonal	Missense	rs185835888
AHNAK2	14	105418798	T	C	c.2990A>G	p.Lys997Arg	SNV	Subclonal	Missense	rs200638785
TXNDC2	18	9887371	C	A	c.895C>A	p.Leu299Ile	SNV	Subclonal	Missense	rs200059339
ZNF57	19	2917901	A	G	c.1282A>G	p.Thr428Ala	SNV	Subclonal	Missense	–
GRIK5	19	42507823	C	T	c.2276G>A	p.Arg759Gln	SNV	Heterozygous	Missense	–
SYNJ1	21	34004062	A	G	c.3965T>C	p.Ile1322Thr	SNV	Heterozygous	Missense	–
PLCXD1	X	215957	C	CG	c.928dupG	p.Val310fs	Insertion	Heterozygous	Frameshift	–
USP7	16	8989553	CTTCTTGAT	CATAA	c.2858_2865delATCAAGAAinsTTAT	p.Gln954fs	Indel	Heterozygous	Frameshift

Copy-Number Alterations

The genome was subject to loss of one copy each of Chromosome 9 and 13 and broad deletions affecting Chromosomes 1, 2, and 6 (Fig. 2A). Combinatorial analysis of copy-number aberrations (CNAs) and SNVs did not identify loss of heterozygosity of any of the altered protein-coding genes. Homozygous deletion of CDKN2A and CDKN2B on Chromosome 9 and a single-copy loss of RB1 were the main copy alterations of note. Comparison of RB1 expression in this tumor sample with RB1 expression across a compendium of all TCGA data sets as well as the average expression across all tissue types in the Illumina BodyMap reference RNA-seq data set was performed. Expression of RB1 was in the 6th percentile compared to TCGA compendium and down-regulated 3.48-fold compared to the Illumina BodyMap, indicating low expression of RB1 in this sample.

Figure 2.

Copy-number aberration and TSC1–TMEM71 structural rearrangement in the G3 PanNET. (A) Circos plot illustrating copy-number alterations and loss of heterozygosity (LOH) observed in the tumor sample. (B) Genomic sequencing read alignment of the genomic regions affected by the rearrangement. Green reads that aligned to TSC1 and purple reads that aligned to TMEM71 indicate split reads that support the rearrangement. (C) Illustration of TSC1–TMEM71 rearrangement generated by MAVIS (Reisle et al. 2018). The TMEM71 and TSC1 breakpoints are indicated at the chromosome, gene, and transcript level (B1 and B2, respectively, upper box). Protein-coding sequence associated with the respective transcripts are indicated by the black line and the amino acids included in the fusion product are indicated. Pfam (http://pfam.xfam.org) protein domains are indicated in the track below the transcript: Hamartin (PF04388). The predicted fusion transcript and protein-coding sequence are shown in the lower box (TSC1 exons are colored green; TMEM71 exons are colored blue).

TSC1–TMEM71 Rearrangement

Structural variants were identified using de novo assembly followed by variant detection. A novel somatic translocation between Chromosomes 8 and 9 (t(9;8)(q34.13;q24.22)) was detected in the genomic assembly (Fig. 2B) and supported by the transcriptomic assembly. The rearrangement resulted in an in-frame fusion of TSC1 exons 1–8 (encoding amino acids 1–246) with TMEM71 exons 8–10 (encoding amino acids 232–277) (Fig. 2C). TSC1 encodes the tumor suppressor, hamartin. Hamartin heterodimerizes with the GTPase, tuberin, encoded by the TSC2 gene. The TSC1/2 protein complex suppresses cell growth, largely by inhibiting the small G-protein Rheb, a crucial activator of the mTORC1 pathway (Castro et al. 2003; Inoki et al. 2003). The predicted TSC1–TMEM71 fusion protein identified in this tumor lacks protein sequence critical for TSC2 binding and is thus predicted to disrupt TSC1 function (Huang and Manning 2008). Deletion of exon 9 alone has been shown to disrupt TSC1–TSC2 dimer formation and downstream inhibition of mTORC1 kinase activity (Santiago Lima et al. 2014), supporting this fusion as indeed a loss of function variant. Moreover, loss of heterozygosity as a result of whole Chromosome 9 loss is predicted to render the tumor deficient in TSC1 activity (Fig. 2B). There is a significant reduction in the aligned RNA-seq read coverage of TSC1 exons 9–23, downstream from the breakpoint, compared to upstream of the breakpoint, in support of the hypothesis that this tumor lacks a functional full-length copy of TSC1 (Supplemental Fig. S2).

CHD7–BEND2 Rearrangement

A second novel somatic translocation between Chromosome 8 and the X chromosome was detected in the genomic assembly (t(8;X)(q12.2;p22.13)) (Fig. 3A), which was also supported by the transcriptomic assembly. The rearrangement results in a novel in-frame fusion between CHD7 exons 1–2 (encoding amino acids 1–555) and BEND2 exons 5–14 (encoding amino acids 165–800) (Fig. 3B). CHD7 encodes the chromodomain helicase DNA-binding protein 7, a chromatin remodeling enzyme involved in differentiation and transcription regulation (Schnetz et al. 2010; Feng et al. 2013). Germline loss of function mutations or deletions of CHD7 are found in patients with CHARGE syndrome, a disorder characterized by dysmorphic features and congenital anomalies in multiple organs. CHD7 remodeling activity is required for neural crest cell gene expression networks, linking a lack of CHD7 activity with CHARGE syndrome features (Bajpai et al. 2010). The translocation observed in this G3 PanNET disrupts the majority of the protein-coding region of CHD7, resulting in loss of the chromodomains and helicase domain (Fig. 3B).

Figure 3.

CHD7–BEND2 structural rearrangement and increased BEND2 expression. (A) Genomic sequencing read alignment of the genomic regions affected by the rearrangement. Gold reads that aligned to CHD7 and purple reads that aligned to BEND2 indicate reads that support the rearrangement. (B) Illustration of CHD7–BEND2 rearrangement. The CHD7 and BEND2 breakpoints are indicated at the chromosome, gene, and transcript level (B1 and B2, respectively, upper box). The protein-coding sequences associated with the respective transcripts are indicated by the black line, and the amino acids included in the fusion product are indicated. The predicted fusion transcript and protein-coding sequence are shown in the lower box (CHD7 exons are colored blue; BEND2 exons are colored green). Pfam protein domains are indicated by the tracks below the transcript: SNF2_N (F00176), helicase_C (PF00271), Chromo (PF00385), BRK (PF07533), BEN (PF10523). (C) Genomic sequencing reads (top) and RNA sequencing reads (bottom) aligned to the BEND2 gene. Genomic reads that are affected by the rearrangement are colored purple. The BEND2 fusion breakpoint is indicated by the dashed line. RNA-seq coverage downstream from the breakpoint (exons 5–14) is increased compared to upstream of the breakpoint (exons 1–4). (D) Boxplot of BEND2 RPKM (log10 transformed) across a panel of 40 TCGA data sets. BEND2 log10(RPKM) for this patient sample is indicated by the red horizontal line.

The CHD7 fusion partner BEND2 encodes BEN domain-containing 2, a relatively uncharacterized protein. Alignment of RNA-seq reads across the BEND2 genomic loci revealed an increase in expression of BEND2 downstream from the breakpoints, indicating the fusion of BEND2 with CHD7 enhances transcription of BEND2 (Fig. 3C). Furthermore, BEND2 was found to be an extreme high expression outlier in this tumor compared to all TCGA disease types (Fig. 3D). BEND2 protein is solely characterized by the presence of two BEN domains, which have been demonstrated to bind DNA and are typically found in proteins involved in transcription and chromatin regulation (Abhiman et al. 2008; Sathyan et al. 2011; Dai et al. 2013, 2015). Interestingly, a fusion between EWSR1 amino-terminal transactivation domain and BEND2 was identified in two patients in the ICGC study, one of which was a grade 3 PanNET (Scarpa et al. 2017). Moreover, the EWSR1–BEND2 rearrangement also similarly resulted in a fusion transcript containing exons 9–14 of BEND2, indicating a common breakpoint among these BEND2 rearrangements. In addition, a recurring fusion between MN1 and BEND2 was identified in three patients with central nervous system primitive neuroectodermal tumors (Sturm et al. 2016). This rearrangement also resulted in the juxtaposition of a transactivation domain donated by MN1 with BEND2. The cumulative evidence leads to the hypothesis that activation of BEND2, either by fusion with a transactivation domain or potentially by increased expression, may promote oncogenic activity.

Gene Expression Analysis

RNA-seq data were compiled as gene-specific RPKM values for all seven PanNEN samples (six PanNETs and one PanNEC) sequenced as part of the POG initiative. As these samples represented several different tumor grades (two grade 1, three grade 2, and two grade 3), we evaluated the impact of tumor grade on gene expression. A mixed-effects model fit to the principal components showed that tumor grade accounted for 5% of the variability in the transcriptomic data. Based on this observation, elastic net regression was used to identify genes that could characterize tumor grade in these samples (Supplemental Fig. S3). This approach revealed genes associated with PanNEN grade. Along with the MKI67, which encodes Ki-67, high expression of RIPPLY3, CIT, YBX2, and KIF19 was also associated with the G3 PanNEN signature. RIPPLY3 is expressed in islet β cells of the pancreas and may play a role in regulating proliferation and differentiation of β cells (Osipovich et al. 2014). Last, hierarchical clustering of the samples based on highly variable genes revealed that tumors with similar grade did not cluster together (Supplemental Fig. S4). These experiments suggested that tumor grade did not significantly bias global gene expression in these samples, but revealed potential gene signatures associated with tumor grade that merit further study.

Somatostatin receptors (SSTRs) are targeted by the peptide hormone somatostatin, eliciting an antiproliferative effect (Florio 2008; Rinke et al. 2009; Caplin et al. 2014). Somatostatin analogs (SSAs), such as octreotide, have played a role in the clinical setting for both tumor localization and therapeutic intervention. Given the clinical relevance of SSTRs, we investigated the expression of SSTR-encoding genes in this patient sample. There are five SSTR genes in the human genome that have been demonstrated to be heterogeneously expressed across gastro-entero-pancreatic neuroendocrine tumors (John et al. 1996; Papotti et al. 2002). To assess the relative SSTR expression in this patient sample, RNA expression of SSTR1–5 was compared among this patient sample and six additional metastatic PanNENs (Wong et al. 2018). In agreement with a previous report (Papotti et al. 2002), SSTR1 and SSTR2 were most frequently expressed across the cohort followed by SSTR3 and SSTR5 (Supplemental Fig. S5). Normalized SSTR gene expression was consistently found to be low in this patient sample compared to other PanNEN samples. This observation is in line with SSA scintillation imaging, which was negative for this patient, providing a rationale for not pursuing SSA-based therapy (Fig. 4).

Figure 4.

Low expression of SSTR genes in the G3 PanNET. SSTR1, SSTR2, SSTR3, and SSTR5 RPKM normalized to the mean of the cohort were plotted for seven metastatic PanNET samples. Relative SSTR gene expression for this patient sample is highlighted in red.

DISCUSSION

Single-Nucleotide and Copy-Number Variants

Here, we describe the molecular characterization of a patient with the recently recognized distinct class of PanNEN, PanNET G3. Similar to other reports that have investigated the PanNET genomic landscape across the broad spectrum of classifications, the mutational burden of this PanNET G3 sample was relatively low compared to pancreatic ductal adenocarcinoma. Only 15 nonsynonymous SNVs and indels were identified, many of which were subclonal. Despite a high rate of occurrence in PanNETs (Jiao et al. 2011; Yachida et al. 2012; Scarpa et al. 2017), mutations or copy losses were not identified in DAXX, ATRX, or MEN1. High-grade poorly differentiated PanNECs were found to have recurrent RB1 mutations and reduced protein expression, a characteristic not frequently observed in G1 and G2 PanNETs (Yachida et al. 2012; Bergsland et al. 2016). A single-copy loss and low expression of RB1 and homozygous loss of the upstream cell cycle–negative regulators CDKN2A and CDKN2B, which have also been observed to be mutated in a subset of PanNECs (Bergsland et al. 2016), may also contribute to the high-grade nature of this G3 PanNET.

Structural Rearrangements

Loss-of-function mutations in TSC1 and TSC2 occur in a subset of PanNETs (Jiao et al. 2011; Scarpa et al. 2017). Moreover, patients with the hereditary disorder tuberous sclerosis have an increased incidence of low-grade, well-differentiated PanNETs, further supporting a role for these tumor suppressors in PanNET oncogenesis (Larson et al. 2012). Mutations in other components of the PI3K/mTOR signaling pathway including PTEN, DEPDC5, and PIK3CA have also been observed (Jiao et al. 2011; Yachida et al. 2012; Wong et al. 2018). A novel somatic structural variant disrupting the TSC1 gene was identified in this tumor sample. The structural variant resulted in an in-frame fusion with the transmembrane protein-encoding gene, TMEM71, resulting in loss of key residues required for TSC1/2 protein complex formation. Gene fusions are more typically associated with oncogene activation as opposed to loss of function of tumor suppressors, rendering this a relatively uncommon mechanism of gene inactivation.

The mTOR inhibitor, everolimus, has been shown to significantly increase progression-free survival in phase II and III trials and is approved for use in the PanNET setting (Yao et al. 2008, 2010, 2011, 2016). Given the incidence of mTOR pathway alterations in PanNETs, there is further rationale for treatment of patients predicted to have an aberrant mTOR pathway with mTOR inhibitors. Continued investigation will be required to determine if TSC1/TSC2-deficient PanNETs have a differential response to mTOR inhibitors compared to the other PanNETs and how this impacts response to everolimus in grade 3 well-differentiated PanNETs, which were not included in RADIANT 3 (Yao et al. 2011). The majority of reported inactivating TSC1/2 variants detected in PanNETs are SNVs or indels. The aberrant TSC1 fusion identified in this tumor highlights a potentially underappreciated mechanism of tumor-suppressor inactivation by structural rearrangement. Incorporation of structural variant detection for tumor suppressors including TSC1 and other genes associated with potentially actionable pathways may further increase the proportion of patients that could benefit from targeted therapies.

A novel rearrangement between CHD7 and BEND2 disrupting the majority of CHD7 and driving expression of BEND2 was identified in this patient sample. It is unclear what the biological impact of such a fusion would have on the tumor, if any; however, recurrent fusions affecting BEND2 in this disease type indicate it may contribute to tumor biology (Scarpa et al. 2017). Although the function of the BEND2 protein has not been studied, the presence of two BEN domains suggests it may play a role in chromatin regulation and gene expression. Studies of other proteins that encode BEN domains have revealed a direct and motif-specific interaction of the BEN domain with DNA (Dai et al. 2013, 2015). Furthermore, another BEN domain–containing protein, BEND3, specifically binds heterochromatin and is proposed to contribute to transcriptional repression (Sathyan et al. 2011). Another homolog, BEND5, has been shown to similarly specifically localize to neural cortex and also stimulate transcription repression (Dai et al. 2013). Further investigation into the function of BEND2 will be required to assess the relevance of BEND2 rearrangements found in a subset of PanNETs.

SSTR Expression

SSA tumor localization was negative in this patient, indicating a lack of SSTR expression on tumor cells. In agreement with the pathology assessment, gene expression analysis of the five SSTR genes revealed low expression compared to a panel of previously reported metastatic PanNETs (Wong et al. 2018). Given the application of SSAs in the therapeutic setting, reduced SSTR expression as determined by lack of SSA imaging and gene expression analysis may indicate inherent resistance to SSA-based treatment.

CONCLUSION

This report represents the complete genomic and transcriptomic description of a well-differentiated grade 3 PanNET, which to date has had limited molecular characterization. Whole-genome sequencing revealed a lack of variants in the frequently mutated genes ATRX, DAXX, and MEN1 and identified novel structural rearrangements. Discrimination of this class of tumor in future studies will be critical for understanding the molecular mechanisms that drive their phenotypic distinction.

METHODS

Sample Collection and Processing

Following informed consent, patients underwent image-guided metastatic biopsies as part of the POG program of British Columbia (NCT02155621). Eight core biopsies using 18-22G biopsy needles were embedded in optimal cutting temperature (OCT) compound for subsequent sectioning and sequencing library preparation. Tumor sections were reviewed by a pathologist to confirm the diagnosis, evaluate tumor content and cellularity, and select areas most suitable for DNA and RNA extraction. Peripheral venous blood samples were obtained at the time of biopsy and leukocytes isolated for use as a germline reference. DNA and RNA were extracted for genomic and transcriptomic library construction, which have been previously described in detail (Jones et al. 2017).

Sequencing and Bioinformatics

Paired-end reads were generated on an Illumina HiSeq2500 sequencer and aligned to the human reference genome (GSCh37; available from http://www.bcgsc.ca/downloads/genomes/9606/hg19/1000genomes/bwa_ind/genome) by the BWA aligner (v0.7.6) (Li and Durbin 2010). Sequencing read lengths were 125 bp. Somatic single-nucleotide variants (SNVs) and small insertions/deletions were processed using SAMtools (Li et al. 2009), MutationSeq (v4.3.5) (Ding et al. 2012), Strelka (v1.0.6) (Saunders et al. 2012), and Trans-ABySS (v1.4.10) (Simpson et al. 2009). Regions of CNA were determined using CNAseq (v0.0.8) (Jones et al. 2010) and LOH by APOLLOH (v0.1.2) (Ha et al. 2012). Tumor content and ploidy models were estimated from sequencing data through analysis of the CNA ratios and allelic frequencies of each chromosome. This was then compared to in-house theoretical models for different ploidy at various tumor contents. Tumor content and sequencing coverage were as shown in Table 2.

Table 2.

Sequencing coverage table

Library ID	RNA/DNA	Aligned Reads	Coverage
P02594 (Matched Normal)	DNA (PCR-free WGS)	993322416	37×
P02600 (Tumor)	DNA (PCR-free WGS)	2386018119	94×
P02603 (Tumor)	RNA (strand-specific)	371042140

Structural variants were detected by de novo assembly of tumor reads using ABySS (v1.3.4) (Birol et al. 2009) and Trans-ABySS (Robertson et al. 2010), followed by variant discovery using DELLY (v0.7.3) (Rausch et al. 2012). Structural variants were visualized using MAVIS (Reisle et al. 2018). Genomic and transcriptomic sequencing read alignments were visualized using Integrative Genomics Viewer (http://software.broadinstitute.org/software/igv/). Reads were colored by insert size and chromosome of the paired read to reveal reads that supported the structural rearrangements.

Mutation Signature Analysis

The mutation signature profile was determined by classifying genomic SNVs into 96 classes based on variant and 3′/5′ mutation context to obtain a mutation catalog vector as described by Alexandrov et al. (2013). To determine the best fit to a consensus set of 30 mutation signatures (available at http://cancer.sanger.ac.uk/cosmic/signatures), we used Markov chain Monte Carlo (MCMC) to determine the full Bayesian solution posterior of the exposure vector e in the multinomial mixture model m = Se using the signIT package (https://github.com/eyzhao/SignIT/) (Zhao et al. 2017). S is a known 96 × 30 matrix of consensus signatures, m is the 96-element mutation catalog, and the members of e denote each signature's relative contribution to the overall mutation burden. The MCMC posterior sample was examined to model the variance inherent to m, the distributions of which are shown in Supplemental Figure S1. Where confidence intervals are reported in the text, 95% simulated interval was reported using the 2.5 percentile and 97.5 percentile values for each signature exposure. This allows estimation of statistical confidence, with the caveat that it relies upon the underlying assumption that the 30-signature reference set accurately reflects the space of possible signatures.

Gene Expression Analysis

RNA-seq reads were processed as previously described in Sheffield et al. (2015). Gene expression was determined as the number of reads per kilobase of transcript per million mapped reads (RPKM). For SSTR gene expression analysis, SSTR1-5 RPKM was compared across a cohort of seven PanNENs that were sequenced as part of the POG project (Wong et al. 2018). Normalized SSTR1, SSTR2, SSTR3, and SSTR5 expression was calculated by dividing RPKM by the mean RPKM for the cohort. For pan-TCGA BEND2 expression comparison, log10 transformed BEND2 RPKM from TCGA-sequenced tumor samples representing 40 disease and disease subtype data sets was plotted as a boxplot. BEND2 RPKM, as determined using the GAF 2.0 gene model (geneontology.org), was plotted as a red horizontal line.

The effect of tumor grade on the transcriptome was evaluated using the PVCA package in R (http://bioconductor.org/packages/release/bioc/html/pvca.html). Expression signatures for tumor grade were determined using elastic net regression. An elastic net model was fit to the whole transcriptome data, with the predictive outcome being tumor grade. Highly predictive features were quantified using weight analysis of the fitted model and clustered using Euclidean distance and complete hierarchical clustering.

Clustering based on gene expression (RPKM) was performed on the top 1% variable genes selected from seven PanNENs of varying grade that were sequenced as part of the POG initiative and plotted as a heatmap. Rows (genes) and columns (samples) were clustered using the average and Manhattan methods, respectively. Tumor histology, Ki-67 staining percentage, grade, sex, and tumor content as determined by bioinformatics analysis (biofx_tc) were indicated in the horizontal tracks.

ADDITIONAL INFORMATION

Data Deposition and Access

The whole-genome sequencing and RNA-seq data for this case are available as .bam files from the European Genome-phenome Archive (EGA; www.ebi.ac.uk/ega/home) as part of the study EGAS00001001159, accession ID EGAD00001004715. The TSC1-disrupting and CHD7–BEND2 gene fusions were submitted to ClinVar (https://www.ncbi.nlm.nih.gov/clinvar/) and can be found under accession numbers SCV000897763 and SCV000897764, respectively.

Ethics Statement

Patients provided written informed consent for metastatic biopsies, sequencing, and publication of results as part of the Personalized OncoGenomics Program of British Columbia (NCT 02155621, University of British Columbia Clinical Research Ethics Board approval no. H14-00681).

Acknowledgments

We gratefully acknowledge the participation of our patients and families, the POG team, and the generous support of the BC Cancer Foundation. We also acknowledge contributions toward equipment and infrastructure from Canada Foundation for Innovation and the BC Knowledge Development Fund. The results published here are in part based on data generated by the following projects and obtained from dbGaP (http://www.ncbi.nlm.nih.gov/gap): The Cancer Genome Atlas managed by the NCI and NHGRI (http://cancergenome.nih.gov) and the Genotype-Tissue Expression (GTEx) Project supported by the Common Fund of the Office of the Director of the National Institutes of Health (https://commonfund.nih.gov/GTEx).

Competing Interest Statement

The authors have declared no competing interest.

Referees

Daniel M. Halperin

Anonymous

Funding

Funding, under the form of unrestricted philanthropic donations, was generously provided by the BC Cancer Foundation (NRF10229) and Genome British Columbia (project B20POG). We also acknowledge contributions toward equipment and infrastructure from Genome Canada and Genome British Columbia (projects 202SEQ, 212SEQ, 12002), Canada Foundation for Innovation (projects 20070, 30198, 30981, 33408), CIHR Foundation (FDN-143288) programs, and the BC Knowledge Development Fund. Additionally, M.L.T. is supported by the University of British Columbia Clinician Investigator Program.