Article Text

Genome analysis identifies differences in the transcriptional targets of duodenal versus pancreatic neuroendocrine tumours
  1. Karen Rico1,
  2. Suzann Duan1,
  3. Ritu L Pandey2,
  4. Yuliang Chen2,
  5. Jayati T Chakrabarti2,
  6. Julie Starr3,
  7. Yana Zavros2,
  8. Tobias Else4,
  9. Bryson W Katona3,
  10. David C Metz3,
  11. Juanita L Merchant1
  1. 1Department of Medicine, University of Arizona Medical Center - University Campus, Tucson, Arizona, USA
  2. 2Department of Cellular and Molecular Medicine, University of Arizona Medical Center - University Campus, Tucson, Arizona, USA
  3. 3Department of Internal Medicine, Division of Gastroenterology, University of Pennsylvania Health System, Philadelphia, Pennsylvania, USA
  4. 4Department of Internal Medicine-Endocrinology, University of Michigan, Michigan Medicine, Ann Arbor, Michigan, USA
  1. Correspondence to Dr Juanita L Merchant; jmerchant{at}email.arizona.edu

Abstract

Objective Gastroenteropancreatic neuroendocrine tumours (GEP-NETs) encompass a diverse group of neoplasms that vary in their secretory products and in their location within the gastrointestinal tract. Their prevalence in the USA is increasing among all adult age groups.

Aim To identify the possible derivation of GEP-NETs using genome-wide analyses to distinguish small intestinal neuroendocrine tumours, specifically duodenal gastrinomas (DGASTs), from pancreatic neuroendocrine tumours.

Design Whole exome sequencing and RNA-sequencing were performed on surgically resected GEP-NETs (discovery cohort). RNA transcript profiles available in the Gene Expression Omnibus were analysed using R integrated software (validation cohort). Digital spatial profiling (DSP) was used to analyse paraffin-embedded GEP-NETs. Human duodenal organoids were treated with 5 or 10 ng/mL of tumor necrosis factor alpha (TNFα) prior to qPCR and western blot analysis of neuroendocrine cell specification genes.

Results Both the discovery and validation cohorts of small intestinal neuroendocrine tumours induced expression of mesenchymal and calcium signalling pathways coincident with a decrease in intestine-specific genes. In particular, calcium-related, smooth muscle and cytoskeletal genes increased in DGASTs, but did not correlate with MEN1 mutation status. Interleukin 17 (IL-17) and tumor necrosis factor alpha (TNFα) signalling pathways were elevated in the DGAST RNA-sequencing. However, DSP analysis confirmed a paucity of immune cells in DGASTs compared with the adjacent tumour-associated Brunner’s glands. Immunofluorescent analysis showed production of these proinflammatory cytokines and phosphorylated signal transducer and activator of transcription 3 (pSTAT3) by the tumours and stroma. Human duodenal organoids treated with TNFα induced neuroendocrine tumour genes, SYP, CHGA and NKX6.3.

Conclusions Stromal–epithelial interactions induce proinflammatory cytokines that promote Brunner’s gland reprogramming.

  • gastrin
  • neuroendocrine tumours
  • RNA expression
  • immunohistochemistry
  • inflammatory mediators
  • TNFalpha
  • organoids

Data availability statement

Data are available upon reasonable request.

http://creativecommons.org/licenses/by-nc/4.0/

This is an open access article distributed in accordance with the Creative Commons Attribution Non Commercial (CC BY-NC 4.0) license, which permits others to distribute, remix, adapt, build upon this work non-commercially, and license their derivative works on different terms, provided the original work is properly cited, appropriate credit is given, any changes made indicated, and the use is non-commercial. See: http://creativecommons.org/licenses/by-nc/4.0/.

Statistics from Altmetric.com

Summary box

What is already known about this subject?

  • Gastroenteropancreatic neuroendocrine tumours (GEP-NETs) are rare tumours that are typically characterised by location and secretory products.

  • The most malignant of these tumours are gastrinomas, which can develop sporadically or as part of multiple endocrine neoplasia type 1 (MEN1).

  • Of MEN1 gastrinomas 60% occur in the proximal duodenum.

  • Genomic analysis of GEP-NETs is more frequently reported for pancreas (pancreatic neuroendocrine tumours, PNETs) and not the duodenum (duodenal neuroendocrine tumours, DNETs).

What are the new findings?

  • Gene profile of a subset of DNETs expressing gastrin or duodenal gastrinomas (DGAST) compared with PNETs shows higher expression of calcium-dependent and stromal-related genes.

  • Despite a paucity of immune cells within the tumour, Interleukin 17 (IL-17) cytokines expressed by the tumour itself induce phosphorylated signal transducer and activator of transcription 3 (pSTAT3).

  • By contrast, tumor necrosis factor alpha (TNFα) was expressed in both the tumour and the stroma.

How might it impact on clinical practice in the foreseeable future?

  • DGAST tumours arise in the Brunner’s glands as scattered IL-17/TNFα/pSTAT3/SYP-expressing cells, suggesting that targeting these pathways might provide additional therapeutic options.

Introduction

Progress in genetic testing, advances in imaging and increased awareness are the major reasons for the increased incidence in gastroenteropancreatic neuroendocrine tumours (GEP-NETs).1–3 Although GEP-NETs comprise ~0.5% of all gastrointestinal (GI) cancers and are the second most prevalent GI tumour class after colon cancer,4–6 their rising prevalence has not yet provided further insights into their cellular origins. Nevertheless, the leading hypothesis is that this subcategory of neuroendocrine tumours (NETs) arises from enteroendocrine cells. Enteroendocrine cells display gene signatures common to both endocrine and neuronal cell types, fuelling controversy as to whether GEP-NETs develop from cells of endodermal or neural crest origin.7 8 A recent multiomics study of 25 different organoid lines derived from GEP-NETs suggests that the pathogenesis of these tumours involves the reprogramming of epithelial cell populations to exhibit a neuronal phenotype.9 Moreover the neuronal reprogramming is WNT and epidermal growth factor (EGF) independent, suggesting that activation of either or both extracellular signal regulated kinase 1/2 (ERK1/2) or AKT signalling occurs through alternative pathways. Yet despite extensive genomic analysis of NETs from the oesophagus, stomach, pancreas and colon, only one duodenal neuroendocrine carcinoma organoid was generated and only one organoid was from a gastrinoma identified exclusively in the gallbladder,9 underscoring the scarcity of this subset of GEP-NETs for evaluation. Of the GEP-NETs, the transcriptome of duodenal neuroendocrine tumours (DNETs) has been the least studied, and it has been assumed that they are similar to other small intestinal carcinoids or pancreatic neuroendocrine tumours (PNETs).9–12 In particular, many DNETs are duodenal gastrinomas (DGASTs) primarily due to their relationship with the multiple endocrine neoplasia type 1 (MEN1) syndrome.13

GEP-NETs are associated with MEN1, an autosomal dominant condition in which a non-functional, germline MEN1 allele renders neuroendocrine cells susceptible to transformation or loss of heterozygosity at the remaining wild-type MEN1 locus.14–16 MEN1 predisposes to neuroendocrine tumours or hyperplasia in the pituitary, parathyroid and foregut-derived tissues (thymus, lung, ECLomas), with the most predominant foregut tumours arising in the gastric antrum, proximal duodenum and pancreas.17 Compared with the more prevalent non-malignant pituitary and parathyroid gland tumours, GEP-NETs associated with MEN1 are at greater risk of metastasis, with non-functional and gastrin-secreting tumours (gastrinomas) being the most malignant.18 19 About 60% of MEN1-related gastrinomas occur in the submucosal Brunner’s glands (BGs) of the duodenum, are multiple and small (<1 cm), and frequently metastasise.13 19 20 In addition to their submucosal location, DGASTs also exhibit enteric glial markers.21

Anlauf et al22 initially reported that DGASTs arise within the BGs, which normally do not express GI hormones, including gastrin. Apparently, DGASTs do not arise from enteroendocrine cells of the small intestinal mucosa, suggesting reprogramming of the epithelial cell population in BGs. Therefore, the focus of the current study was to identify potential gene drivers of DNET tumours, specifically those considered gastrinomas (DGAST). RNA-sequencing (RNA-seq) of DGASTs and PNETs was performed to identify signalling pathways distinct to DGASTs. Genes more highly expressed in the DGASTs were validated through mining an online Gene Expression Omnibus (GEO) profile database and then using both immunohistochemistry (IHC) and digital spatial profiling (DSP) to compare expression in the BG versus DGAST tumours.

Materials and methods

Human samples

Five DGASTs and nine PNETs with corresponding blood samples were collected during surgery at the University of Pennsylvania (UPENN) and stored in liquid nitrogen (table 1). These de-identified biosamples were then transferred to the University of Arizona (UA) for genomic analysis. The UA Comprehensive Cancer Center (UACC) Biospecimen Repository (P30 CA23074) provided adjacent normal-appearing duodenum and pancreas obtained from a Whipple procedure for pancreatic adenocarcinoma. Additional normal pancreas (nPANC) and duodenum (nDUO) were obtained from a commercial source (BioChain, Newark, California). Repositories from the University of Michigan, UPENN and UACC provided formalin-fixed paraffin-embedded (FFPE) de-identified GEP-NETs and human duodenal organoids. Duodenal organoids were generated from normal duodenal biopsies obtained from the UA TARGHETS biorepository .

Table 1

DGAST and PNET patient demographics

Expression analysis of GEP-NET GEO data sets

RNA-seq counts from GSE98894 were obtained from the GEO database (http://www.ncbi.nlm.nih.gov/geo/). GSE98894 contained 194 samples comprising 113 PNETs and 81 small intestinal neuroendocrine tumours (SI-NETs).23 Out of these 113 PNETs, 83 were from primary tumours, and of the 81 SI-NETs 44 were from primary tumours. The remainder were from metastases. Differential analysis of PNET and SI-NET transcripts was performed using LIMMA.24 25 Genes differentially expressed with adjusted p values <0.05 were analysed for biological functions using the Kyoto Encyclopedia of Genes and Genomes (KEGG). Pearson’s correlation coefficient was used to identify genes that correlated with GASTRIN (GAST) from the 81 SI-NETs. Analysis of RNA-seq databases, correlations, Venn diagrams and heatmaps were performed using R (V.3.6.3).

RNA-seq and whole exome sequencing

RNA was extracted from four DGASTs, six PNETs, three nDUOs and two nPANCs using the QIAGEN RNeasy Mini Kit for tissue (Venlo, The Netherlands). RNA-seq was performed on Illumina platforms by Novogene and then analysed using fragments per kilobase of transcript per million mapped reads (FPKM) for Pearson’s correlation, principal component analysis and coexpression Venn diagrams (FPKM threshold set to 1). Differential expression analysis was conducted using the DESeq2 R package.26 Cluster analysis was carried out using log10(FPKM +1) of the differentially expressed genes (DEGs) in the four tissue groups. DEGs were then analysed with clusterProfiler for KEGG database enrichment.27 DNA was extracted from the buffy coat and tumour using QIAGEN PDNA Mini Kit (Cat #51304). Whole exome sequencing and analysis were performed by Novogene. All samples met quality control standards, which required Q30 above 80%. Sequencing results were compared with the reference CRCh37 human genome. The Genome Analysis Toolkit (GATK) was used to identify single nucleotide polymorphisms and small insertions and deletions from Binary Alignment Map (BAM) files. ANNOVAR was used to annotate variants.

Human duodenal organoids

Healthy nDUO for RNA-seq and whole exome sequencing (WES) was collected in Advanced DMEM/F12 (Gibco) during Whipple surgery. Duodenal organoid cultures were established from biopsies of four different individuals undergoing upper endoscopy by adapting the method of Spence and coworkers.28 29 Organoid lines were thawed from frozen stocks stored in liquid nitrogen and subcultured in complete organoid media consisting of a 1:1 (v/v) solution of L-Wnt3a/R-Spondin3/Noggin (WRN)- Conditioned Medium and Advanced DMEM/F12 (Thermo Fisher Scientific), 10 mM 4-(2-hydroxyethyl)-1-piperazineethanesulfonic acid (HEPES) (Invitrogen), 2 mM GlutaMAX (Invitrogen), 1 mM N-acetylcysteine, 1X N2 Supplement (Invitrogen), 1X B-27 Supplement Minus Vitamin A (Invitrogen), and 1X Penicillin-Streptomycin (Invitrogen). L-WRN Conditioned Medium was generated by subculturing L-Wnt3a/R-spondin 3/Noggin-(WRN)-expressing cells in AdvancedDMEM/F12 media supplemented with 20% fetal bovine serum and 1X Penicillin-Streptomycin. Conditioned Medium was collected every 2 days and sterile filtered. Complete organoid medium was supplemented with 100 ng/mL recombinant human EGF (R&D Systems), 2.5 µM CHIR99021, 10 µM Y-27632 and 500 nM A83-01 (all from Tocris, Minneapolis, Minnesota). Media were supplemented with 10 µM SB202190 and 10 nM nicotinamide (Tocris) during passaging and re-embedding steps to prevent anoikis. Organoids were used between passage numbers 3 and 7 and treated between days 10 and 14 postpassaging.

Human duodenal organoids were treated with 5 and 10 ng/mL of recombinant human TNFα (BioLegend, San Diego, California) diluted in complete organoid medium without A83-01 for 24 hours (for evaluating messengerRNA [mRNA] and Nkx6.3 protein expression) and 48 hours (for evaluating protein in nuclear and cytoplasmic fractions). Total RNA was extracted from duodenal organoids using the ReliaPrep miRNA Cell and Tissue Miniprep System (Promega, Madison, Wisconsin) following the manufacturer’s protocol. Copy DNA (cDNA) was synthesised using SuperScript VILO IV following treatment with ezDNAse to remove genomic DNA. Quantitative PCR (qPCR) was performed as described previously using 20 ng of cDNA. The qPCR programme consisted of initial 2 min denaturation at 95°C, followed by 40 cycles of 15 s at 95°C to denature and 1 min at 62°C to anneal. qPCR conditions were the same for all genes in online supplemental table 1. RNA input was normalised to Hypoxanthine phosphoribosyltransferase 1 (HPRT) mRNA and expressed as relative fold change using the 2−double deltaCt method.30

Immunohistochemistry

Tissues were fixed overnight in 10% formalin, de-paraffinised in xylene and then rehydrated in 70% ethanol. Antigen retrieval was performed in citrate buffer, pH 6.0, or Tris-EDTA, pH 9.0, for 30 min, and then blocked for 1 hour in 10% donkey serum diluted in Tris-buffered saline 0.1% Triton with 0.025% Tween-20 (TBST-TW). For immunofluorescence (IF), the primary antibodies were incubated in 1% bovine serum albumin (BSA) TBS-TW overnight at 4°C. Secondary antibodies (1:500) were incubated for 1 hour in 1% BSA with 0.025% Tween-20 (TBS-TW) at room temperature (RT) in the dark. For IHC, primary antibodies were incubated in 1% BSA Tris-buffered saline 0.1% Tween (TBS-T) overnight at 4°C. Slides were incubated in horseradish peroxidase (Cell Signaling, Danvers, Massachusetts) secondary antibody (1:300) for 60 min in 1% BSA TBS-T for 1 hour at RT. Slides were then incubated with diaminobenzidine (Vector, Burlingame, California) followed by counterstaining with haematoxylin (Leica, Buffalo Grove, Illinois). Dilutions of primary antibodies are listed in online supplemental table 2. Slides were mounted with ProLong Gold Antifade Mounting Medium with 4’, 6-diamidino-2-phenylindole (DAPI) (Invitrogen). Images were photographed on an Olympus BX53F microscope with an Olympus digital camera (Center Valley, Pennsylvania).

Digital spatial profiling

GEP-NETs were analysed using the NanoString Digital Spatial Profiler (DSP) (Seattle, Washington) to quantify a 40-plex panel of Ultraviolet (UV) photocleavable oligonucleotides-linked to neural-related antibodies after hybridising to 5 µm FFPE sections from four DGAST tumours (online supplemental figure 1). The tissue morphology was delineated by the IF detection of PanCK (epithelial cytokeratin, green), CD45+ (immune cells, red) proteins and SYTO13 (nuclei, blue) (NanoString GeoMx for Neural Profiling #GMX-PROCO-NCT-HNCP-12, #GMX-PROMOD-NCT-HADP-12 and #GMX-PROMOD-NCT-HPDP-12).

Figure 1

Transcriptome comparison of DGAST and PNET. (A) Schematic of DGAST and PNET RNA-seq comparison. (B) PCA in three dimensions of DGAST, PNET, nDUO and nPANC. (C) Pearson’s correlation coefficients (R2) for each normal tissue and GEP-NETs, R2 >0.8=biological replicate. (D) Coexpression Venn diagram of patient samples, with default threshold of FPKM set to 1 in each group. (E) Heatmap of DEGs in DGAST (n=4), PNET (n=6), nDUO (n=3) and nPANC (n=2). DGAST, duodenal gastrinoma; FPKM, fragments per kilobase of transcript per million mapped reads; GEP-NETs, gastroenteropancreatic neuroendocrine tumours; nDUO, normal duodenum; nPANC, normal pancreas; PC, principal component; PCA, principal component analysis; PNET, pancreatic neuroendocrine tumour; RNA-seq, RNA-sequencing.

Statistics

For the qPCR experiments comparing two groups and DSP region of interest (ROI) comparisons, the statistical analysis for significance was calculated using an unpaired t-test and a non-parametric Wilcoxon rank-sum test. For qPCR experiments, the significance was calculated using one-way analysis of variance with Tukey’s post-hoc test (GraphPad Prism). P<0.05 was considered significant. For Pearson’s correlation, an R2>0.8 was considered sufficiently significant to be biological replicates. The replicate of the experiments and sample size of the groups are described in the figure legends.

Results

RNA-seq in DGASTs versus PNETs

To compare the transcript differences of DGASTs and PNETs versus nDUO or nPANC, bulk RNA-seq analysis was performed using surgical resections of these GEP-NETs (figure 1A). Table 1 shows the demographics of the available specimens and table 2 shows MEN1 exon mutations identified by WES. A common exon 10 missense variant p.Ala541Thr (rs200921330) is observed in ~3.5% of Caucasian populations.31 Two neuroendocrine markers chromogranin A (CHGA) and synaptophysin (SYP) were elevated in the two tumour types compared with normal tissue (online supplemental figure 2A,B). Phenotypically, three of the four DGASTs exhibited a more restricted hormone expression profile than PNETs, with DGASTs expressing either gastrin or somatostatin transcripts (online supplemental figure 2C).

Figure 2

Differential expression of DGAST compared with nDUO. (A) Enrichment of KEGG pathways for genes expressed uniquely in DGAST tumour samples compared with nDUO and nPANC and PNET. (B) Volcano plot showing statistical significance versus fold change for DGAST versus nDUO DEGs. Significance of upregulated and downregulated genes was determined by adjusted p<0.05. (C) KEGG pathways of significantly upregulated genes based on enrichment analysis of 1342 overexpressed genes and (D) KEGG pathways of significantly downregulated genes based on enrichment analysis of 1291 underexpressed genes in DGAST compared with nDUO. Gene counts on the x-axis plotted as a function of bar colours that reflect the adjusted p value range indicated in the figure. (E) Heatmap of specific DEGs that function as nuclear regulatory factors for three nDUOs and four DGASTs. (F) Heatmap of specific DEGs that comprise signal transduction factors for three nDUOs and four DGASTs. Columns clustered by Z-transformation and Euclidean distance. DEGs, differentially expressed genes; DGAST, duodenal gastrinoma; KEGG, Kyoto Encyclopedia of Genes and Genomes; nDUO, normal duodenum; nPANC, normal pancreas; PNET, pancreatic neuroendocrine tumour. cAMP, cyclic adenosine monophosphate; ECM, extracellular matrix; GnRH, gonadotrophin-releasing hormone.

Table 2

Germline and somatic MEN1 mutations in subjects with DGAST and PNET

In the principal component analysis, FPKMs above 1 from DGASTs and PNETs were more related to each other than to their normal tissues of origin—nDUO and nPANC (figure 1B). Pearson’s correlation confirmed that DGASTs exhibited strong similarity to PNETs (R2>0.8) (figure 1C). Accordingly, these two GEP-NETs shared 1233 common coexpressed transcripts. However, 909 transcripts were distinctly expressed in DGAST compared with nDUO and PNET (figure 1D). Heatmap gene clustering of the DEGs between GEP-NETs showed some overlap between the DGAST and nDUO or PNET and nPANC. Moreover, focusing on those transcripts with adjusted p<0.05 facilitated determining distinct gene differences between the DGAST and PNET transcripts (figure 1E).

Among the 909 distinct DGAST transcripts, 203 encoded proteins. KEGG pathway analysis revealed significant enrichment of neuroactive-ligand receptor interaction and calcium signalling pathway genes (figure 2A and online supplemental figure 3). A volcano plot of the 2633 DEGs identified for DGASTs versus nDUO is shown (figure 2B). Upregulated KEGG pathways were neuroendocrine and calcium-related (figure 2C), while downregulated KEGG pathways in DGASTs relative to nDUO were cytokine and intestine-specific functions, such as nutrient digestion and absorption (figure 2D). In particular, heatmaps of upregulated genes encoding nuclear regulatory factors in the DGASTs highlighted genes associated with endocrine differentiation (HESX1, PAX6, ISL1, SMARCA1/SNF2L) distinct from downregulated signatures referable to intestinal differentiation (eg, CDX1, CDX2, GATA5, KLF5, SPDEF, TBX3) (figure 2E). Reduced Wnt signalling was also evident (TBX3, TCF7 and KLF5), as initially reported by Kawasaki et al9 using single cell RNA-sequencing (scRNA-seq) from NET-derived organoids. Immune-related receptor-signalling pathways, for example, IL-2, FLT3 and CD209, were downregulated in DGAST relative to nDUO (figure 2F). Moreover EGF receptor-related signalling, for example, ERBB2 and MAP2KK2 transcripts, was also reduced (figure 2F), as shown in the scRNA-seq analysis.9

Of the 203 DGAST DEGs compared with PNET, 163 were upregulated (figure 3A). These genes were analysed in a heatmap (figure 3B) and by KEGG pathway analysis (figure 3C). The most significant upregulated DEG genes were smooth muscle-specific genes such as DES, ACTG2, ACTA2, CNN1, FLNC, PLIN4 SYNM, TPM2 MYH11, MYLK, PYGM and LMOD1; stromal genes COL8A1, ASPN and ITGA5; and calcium-regulated genes (CAMK2A, CASQ2, SMTN, CAPN8 and TNC) (figure 3B,C). Transferrin (TF) and its receptor (TFR2) were also upregulated, consistent with expansion of the smooth muscle compartment (figure 3B). Interestingly, many of the growth factors upregulated were related to mesenchymal as opposed to epithelial cell growth, for example, ANGPTL1 (angiopoietin-like 1) CTGF (connective tissue) and SCRG1 (chondrogenesis) (figure 3B). The increase in stromal genes was likely related to the intermixing of tumour with S100B+ and smooth muscle actin (α-SMA)+ mesenchymal cells prevalent in DGASTs (online supplemental figure 1 and figure 3D–H). As reported previously21 and as shown here (figure 3G,H), S100B+ glial cells were often dispersed within the tumour. Receptors expressed in DGASTs relative to PNETs were stromal or neural-related, for example, P2RX2 (purinergic), TACR2 (tachykinin), CHRM1 (cholinergic), GALR2 (galanin), MST1R (macrophage stimulating 1 receptor), FGFR3 (fibroblast growth factor receptor), CD177 (neutrophils), HTR2B (serotonin) and OXTR (oxytocin). DGAST-related transcription factors upregulated relative to PNETs correlated with neuroendocrine differentiation (figure 3I). Signalling pathways overexpressed included tyrosine kinase receptors (PTK6, FGFR3) (figure 3J). Although the tissue landscape for NETs has been considered an ‘immune desert’,32 specific proinflammatory pathways were induced, for example, Angiopoietin-like1 (ANGPTL1), Interleukin 17 receptor E (IL17RE), Interleukin 17B (IL17B) and TNFAIP2 (TNFα-associated inducible protein2) (figure 3K).

Figure 3

DGAST shows distinct gene expression compared with PNET. (A) Volcano plot of 203 DEGs in DGAST versus PNET. (B) Heatmap of DGAST versus PNET genes clustered according to significance (p<0.05). Columns and rows were clustered by Z-transformation and Euclidean distance. (C) KEGG pathway enrichment analysis of genes upregulated in DGAST versus PNET. (D) H&E of a DGAST tumour (20×). (E) Magnified H&E of tumour clusters intercalated between thickened stroma. (F) IF with neuroendocrine marker synaptophysin (SYP, red). (G) IF of SYP+ tumour (labelled) stained with enteric glial marker S100B (white arrow, red) and mesenchymal marker α-smooth muscle actin (α-SMA, green). (H) IF of stromal markers showing S100B within DGAST tumour (white arrow). Inset shows IF for SYP staining the tumour. Nuclear stain: DAPI (blue). Genes from (B) organised into heatmaps of (I) transcription factors, (J) signalling pathways or (K) inflammation-related genes upregulated in DGAST versus PNET. DEGs, differentially expressed genes; DGAST, duodenal gastrinoma; IF, immunofluorescence; KEGG, Kyoto Encyclopedia of Genes and Genomes; PNET, pancreatic neuroendocrine tumour. DN-regulated, down-regulated. cAMP, cyclic adenosine monophosphate; cGMP-PKG, cyclic guanosine monophosphate-protein kinase G.

Validation using GEP-NETs

The GSE98894 public data set generated from biopsies of pancreatic and SI-NETs (online supplemental figure 4) was used to validate genes identified in the DGAST versus PNET analysis (figure 3). We compared 113 PNET with 81 SI-NET transcriptomes from primary and metastatic tumours originating in the pancreas or small intestine. Primary tumour samples comprised 83 PNETs and 44 SI-NETs. This revealed differential transcript clusters between PNETs and SI-NETs that were similar when compared across primary tumours or primary and metastatic samples combined. We compared the 203 significant DEGs in the ‘discovery cohort’ with the GEO transcripts and found over 60% of the genes (124 genes) in common between the two data sets (online supplemental figure 4A,B). Unsupervised clustering segregated GSE98894 transcripts into PNET and SI-NET groups. The genes downregulated in PNETs compared with DGASTs (figure 3B) were also suppressed in PNETs versus SI-NETs (online supplemental figure 4B). The PNET expression pattern showed greater uniformity compared with the SI-NETs due to greater heterogeneity in the intestinal tissue site. The same smooth muscle-specific genes (DES, ACTG2, CNN1, FLNC, PLIN4 SYNM, TPM2 MYH11, MYLK, PYGM, LMOD1) upregulated in SI-NETs (online supplemental figure 4A) were also observed in the DGAST tumours (figure 3B), of which some like ACTG2, MYLK and MYH11 are known IL-17 target genes.33–35 Additional analysis revealed a strong predilection for upregulated IL-17-RELA (NFκb) signalling in SI-NETs compared with PNETs (online supplemental figures 5A and 5B) . Like the discovery DEG data set, the transcripts of the GEO profile ‘validation cohort’ included cell signalling, cell adhesion, transcriptional regulation, neural development and cellular processes, which coincided with specific receptors and ion channels (online supplemental figure 4A). Of the 59 genes in the gastrin-positive SI-NETs (three) found within the GSE98894 cohort, HAP1, ARHGEF38 and SLC26A7 showed the most significant increase and were also identified in our four DGASTs (figure 3B). HAP1 (huntingtin-associated protein) is a neuronal protein involved in vesicle trafficking selectively expressed in antral G cells but not D cells, enterochromaffin (EC) or enterochromaffin-like (ECL) cells.36 Thus, both the ‘discovery’ data set and ‘validation’ cohort from the GEO profiles demonstrated a high degree of concurrence and underscored genes unique to DGASTs.

Figure 4

Digital spatial profiling of DGAST. Protein expression within the regions of interest was determined in abnormal tBG and the DGAST tumour for four patients: (A) SYP, (B) Ki67, (C) S100B, (D) CD45, (E) CD11b, (F) CD40 and (G) CD31. Shown are mean±SEM using Wilcoxon rank-sum test; *p<0.05, ***p<0.0002, ****p<0.0001. Immunofluorescence for IL-17b (white), SYP (red) and nuclei (DAPI, blue) in (H) nBG, (I) tBG and (J) DGAST tumours. Immunofluorescence for TNFα (white), SYP, (red) and nuclei (DAPI, blue) in (K) nBG, (L) tBG and (M) DGAST tumours. BG, Brunner’s gland; DGAST, duodenal gastrinoma; nBG, normal-appearing BG; nDUO, normal duodenum; SYP, synaptophysin; tBG, tumour-associated Brunner’s gland. AU, arbitrary units relative to housekeeping genes

DSP of DGASTs

Anlauf et al22 showed that MEN1-associated DGASTs arise within the BGs and not from enteroendocrine cells within the epithelium. Mesenchymal cell populations, for example, smooth muscle, fibroblasts, vasculature and neural crest-derived neuronal and enteric glial cells, encircled the tumour clusters (figure 3 and online supplemental figure 1). However, since RNA-seq does not provide spatial information, we further characterised DGAST tumours by DSP using a commercial panel of neural and immune antibodies. DSP allows quantitation of proteins in multiple ROIs per slide per patient. The tissue morphology in four patients with DGASTs, including the tumour-associated Brunner’s glands (tBGs), was defined using the cell markers PanCK (epithelial cells), CD45 (immune cells) and SYTO13 (nuclei) (online supplemental figure 6A,B). The epithelium strongly expressed PanCK, while inflammatory cells (CD45) were mainly observed in immune aggregates. Scattered cells within the tBGs were weakly positive for the neuronal marker SYP compared with the DGAST tumour (figure 4A). Moreover, GEP-NETs typically proliferate slowly as indicated by low Ki67 labelling.37 38 Accordingly, an 80% reduction in DGAST Ki67 compared with tBGs was observed (figure 4B). We further confirmed reduced expression of Ki67 in DGASTs compared with tBGs by immunohistochemical staining of tumours and adjacent BGs (online supplemental figure 6C,D). The glial cell marker S100B21 39 also decreased (figure 4C). As previously reported,32 40 41 we also found that DGASTs were nearly devoid of immune cells marked by CD45, CD11b and CD40 (figure 4D–F). Additionally, the vascular marker CD31 (platelet-endothelial cell adhesion molecule, PECAM) was also lower within the tumour clusters (figure 4G). Thus the DSP analysis showed that Ki67, stromal, immune and vascular gene transcripts were elevated in the tBGs relative to the DGASTs, while SYP expression was elevated in the tumour.

DGAST tumours express TNFα and IL-17b cytokines

The DSP suggested that the DGAST tumours lacked infiltrating immune cells. However, IL-17 and TNFα pathway components were identified in the discovery and validation cohorts. Therefore, we used IF to identify which cell populations expressed these proinflammatory cytokines. IL-17b was strongly expressed in the tBG, but not in normal BG (nBG) (figure 4H–J). IL-17 stimulates expression of TNFα.42 43 TNFα pathway components were expressed in the tumours and in small clusters in the stroma as well as the tumour (figure 4K–M). Like IL-17, nBGs did not express TNFα; instead the major source of TNFα occurred in the tBGs and tumours (figure 4K–M). Additionally, the expanded stromal compartment within the tBGs expressed both TNFα and SYP (figure 4K). Therefore we concluded that there was crosstalk between stromal cells and epithelial-derived BGs that contributes to their neoplastic reprogramming.

Since the DGAST tumours expressed both IL-17b and TNFα, we used IF to identify downstream signalling targets (figure 5). TNFAIP2 was identified in the ‘discovery’ cohort (figure 3K) and was found in epithelial and stromal cells of the tBG, but not the tumour (figure 5A,B). RelA/p65 is downstream of both IL-17 and TNFα signalling44 and localised in the epithelial cells of the tBG and tumour (figure 5C,D). However, pSTAT3 was strongly expressed in both the tBG stroma and tumour (figure 5E).

Figure 5

Immunofluorescent detection of TNFAIP2, NFκB and pSTAT3. Immunofluorescence analysis of known signalling targets in tBG versus the DGAST tumour for (A and B) TNFα (white), SYP (red) and TNFAIP2 (green); (C and D) p65 (white), SYP (red) and β-catenin (green); (E) tBG: pSTAT3 (white), SYP (red) and β-catenin (green); and (F) tumour: IL-17b (white), SYP (red) and pSTAT3 (green). DGAST, duodenal gastrinoma; SYP, synaptophysin; tBG, tumour-associated Brunner’s glands.

Stromal–epithelial interactions and cell specification in DGASTs

NKX6.3 is required to pattern the GI tract, particularly gastrin cell specification, and was a transcription factor identified in the DGAST ‘discovery’ cohort (figure 3D).45 Indeed, tumours with the highest GASTRIN transcript levels also showed elevated NKX6.3 (figure 6A). This correlated with IHC colocalisation of NKX6.3 protein to the nucleus and cytoplasm of tBG and DGASTs, respectively (figure 6B–D). The tBGs adjacent to the DGAST showed islands of NKX6.3 protein expression, without exhibiting the characteristic ‘rosette’ morphology of a NET or expressing gastrin (figure 6C,D), suggesting that NKX6.3 contributes to the regulation of cell specification prior to tumour formation and gastrin overexpression.

Figure 6

TNFα induces cell specification protein NXK6.3. (A) Gastrin and NKX6-3 mRNA in DGAST (n=4) compared with nDUO (n=5) determined by qPCR and analysed by Wilcoxon rank-sum (p=ns). (B) Representative H&E of DGAST 269. (C) IHC of NKX6.3 in DGAST 269 (dotted square) in the submucosa; adjacent BG (dotted circle). (D) IF of gastrin expression in DGAST 269. (E) IHC expression of NFκB-p65 in DGAST 269 (dotted square) in the submucosa; adjacent BG (dotted circle). (F) Representative light micrographs of duodenal organoids in culture on day 12 treated after treatment with BSA or TNFα (5 and 10 ng/mL) for 24 hours. (G) Western blot for NFκB-p65 and total STAT3. (H) NKX6.3 protein after TNFα treatment. (I) Western blot quantitation of two unique patient-derived organoids. (J) Time course of TNFα (10 ng/mL) induction of major signalling pathways NFκB, pSTAT3Y705, STAT3 and ERK1/2. Nuclear (N) marker: histone H3; cytoplasmic (C) marker: β-tubulin. (K) Immunofluorescent image of nuclear STAT3 after TNFα treatment. (L) SYP, (M) CHGA, (N) NKX6.3, (O) NGN3 and (P) NEUROD1 mRNA expression relative to HPRT was determined by qPCR using four duodenal organoid lines generated from endoscopic biopsies. Each dot represents an independent experimental replicate, with colour indicating a different patient-derived organoid line. Shown is mean±SEM. Significance determined by one-way analysis of variance with *p<0.05, **p<0.01, ns, not significant. BG, Brunner’s gland; DGAST, duodenal gastrinoma; IF, immunofluorescence; IHC, immunohistochemistry; nDUO, normal duodenum; SYP, synaptophysin.

In silico analysis of the NKX6.3 promoter revealed that (CCAAT-enhancer binding protein β (C/EBPβ), (Nuclear Factor κ B) (NFκB), GATA3 and ATF binding sites are present within the NKX6.3 5’UTR (5’ untranslated region) (GeneCards.org). IHC analysis showed that NFκB colocalised with NKX6.3 protein expression in the tBGs and DGASTs (figure 6C,E). TNFα is a potent inducer of NFκB, which we showed is coexpressed with TNFα expression in DGASTs. To test whether TNFα-treated duodenal organoids induced NKX6.3, organoids from four normal human duodenums were treated with TNFα (figure 6F). TNFα induced NFκb, STAT3 and NKX6.3 proteins (figure 6G–I). NFκb and STAT3 protein levels increased within 3 hours of treating human duodenal organoids with TNFα, while ERK activation peaked early and then decreased (figure 6J). In particular, TNFα induced phosphorylated STAT3 (p-STATY705) in the nucleus, which also regulates SYP expression. Collectively, SYP, CHGA and NKX6.3 transcripts increased, but NGN3 and NEUROD transcripts did not (figure 6L–P). Thus, we concluded that specific proinflammatory cytokines contribute to cell specification of nBG epithelial cells, but not necessarily through traditional developmental pathways for endocrine cells (online supplemental figure 7).

Discussion

The prevalence of GEP-NETs is about 2–5 per 100 000.1 However, due to their indolent behaviour, GEP-NETs are typically advanced by the time they are discovered.4 Most clinical studies have focused on therapies for the more prevalent PNETs and SI-NETs, for example, ileal carcinoids.23 Surgery, somatostatin receptor 2-targeted therapy, chemotherapy and mTOR (mammalian target of rapamycin) inhibition (everolimus) are examples of current treatment modalities.46 47 Epigenetic profiling of PNETs (in contrast to SI-NETs) has uncovered the presence of multiple mutations in genes involved in chromatin remodelling, that is, MEN1, DAXX and ATRX loci,12 although their discovery has not yet redirected current therapeutic options.48 Genetically engineered models (GEMMs) have provided a framework based on developmental paradigms for endocrine cells,49 while profiling of clonal human cell lines has generated lists of genetic targets to consider as novel therapeutic targets and biomarkers.9 However, these targets may not necessarily reflect actual tumour biology. Therefore, in this report, we used a combination of RNA-seq, DSP and organoid cultures to identify the molecular pathways that differ between human duodenal and pancreatic tissues and converge to generate NETs. Since DGASTs are more prone to metastasise,18 the analysis here identified novel insights into the pathogenesis of this GEP-NET subgroup.

Despite an intensive analysis of these tumours, we acknowledge that there are several important limitations to our study. Given the relatively rare occurrence of DNETs, we were only able to evaluate the genomics of a very small sample size, that is, four DGAST tumours (DGAST) and one non-functional DNET. Moreover, adjacent normal tissue from the same patient was not available. As a result, the comparison group was generated from normal duodenum and pancreas. In addition, the DGASTs were a mix of functional MEN1 mutations, that is, frameshift or nonsense. Extrapolation of these results to non-functional DNETs will require additional studies with collections occurring at multiple institutions. Even existing public databases as used here contain primarily SI-NETs without specific clinical information or designation of the site within the small bowel, limiting the ability to distinguish differences that may arise as a function of tumour location.

Nevertheless, our cohort of 10 GEP-NETs used in the RNA-seq analysis revealed that there were significant transcript differences between DGASTs and non-functional PNETs. Specifically, IL-17 and TNFα were two inflammatory cytokines strongly expressed in the DGASTs. In addition, stromal cells also expressed TNFα. IL-17 activates downstream targets in part by activating NFκb and pSTAT3. Moreover, STAT3 binds to the SYP promoter,42 50 51 suggesting that IL-17 expression might coincide with SYP expression. Additionally, IL-17 activates TNFα-related signalling.43

In DGASTs, the gene signatures for transcription factors involved in endocrine cell specification, for example, NKX2.5 and NKX6.3, were identified. Published cell specification results from GEMMs45 with in silico analysis of the NKX6.3 promoter (GeneCards.org) suggested regulation of this transcription factor by NFκb and ostensibly inflammatory cytokines. In Nkx6.3−/− mice, there is a significant decrease in gastrin cells (G cells) and an increase in somatostatin cells (D cells), indicating that Nkx6.3 functions to specify G cell and D cell fate.45 During development, G cells require the basic helix-loop-helix transcription factor neurogenin 3 (Ngn3).52 Other essential transcription factors implicated in G cell specification include Pax4, Pax6 and Pdx1.7 53 Indeed, PAX6 expression was elevated in DGASTs relative to nDUO. Since Anlauf et al22 indicated that the DGASTs arise within the BGs, we used IHC to determine whether there was overlap between cell specification proteins in normal versus aberrant BGs. We found that NKX6.3 protein was expressed in aberrant cells coalescing within the BGs and the DGAST but not in nBG tissue, suggesting that the cascade of genetic events proximal to tumour formation involved activation of NKX6.3 expression perhaps by chronic inflammation.

Progress in understanding GEP-NET pathogenesis is precluded in part by the paucity of tumour specimens. Here, we encountered similar limitations related to sample size in the discovery cohort. Further, it should be noted that our study reflects the broad heterogeneity of these tumours. By performing grouped analyses of functional and non-functional tumours, in addition to MEN1-associated and sporadic tumours, it remains difficult to parse out the specific transcriptional features associated with these subclassifications. Importantly, our findings to date are most relevant to informing the pathogenesis of DGASTs and benefited from prior studies establishing BGs as the site of precursor cells. Further investigations into other duodenal and SI-NETs are warranted to determine whether similar signalling mechanisms are relevant to those tumour types.

In summary, the studies here use multiple technologies to understand the similarities and differences in GEP-NET pathogenesis. An advantage to studying DGASTs is the presence of precursor lesions arising within BGs.22 Moreover, DGASTs employ some elements of tissue-specific endocrine cell specification pathways.45 Traditionally, clinicians separate PNETs from SI-NETs due to differences in their clinical behaviour and responsiveness to therapies. Moreover, DNETs are generally combined with PNETs based solely on anatomical proximity as opposed to different genomic signatures. However, the study here suggests that DNETs and PNETs should be further subdivided due to distinct natural histories that might impact treatment approaches. Although interactions between neuroendocrine and immune cells have been underappreciated, expression of inflammatory cytokines by precursor tBGs and DGASTs themselves suggests that cytokine suppression or immune checkpoint inhibitors might play a therapeutic role in those individuals predisposed to developing GEP-NETs.40 41 54 55

Data availability statement

Data are available upon reasonable request.

Ethics statements

Patient consent for publication

Ethics approval

Sample collection and studies were conducted with institutional review board (IRB) approval from each institution (University of Michigan, UM, #IRB #HUM00115310; UPENN (#814805);TARGHETS Biorepository at the University of Arizona (#1909985869) and UACC (#1808861863) and adhered to the principles of informed consent.

References

Supplementary materials

Footnotes

  • KR and SD are joint first authors.

  • Contributors KR: planning, conduct, wrote the initial draft. SD: planning and execution of organoid treatments, biochemistry studies and writing. RLP: performed bioinformatics and generated the graphs. YC: bioinformatics and generated the graphs. JC and YZ: prepared and validated the human organoids. TE: provided human samples. JS, BK and DCM: collected and validated the human samples. JLM: concept, experimental design, data review, manuscript generation, content guarantor and funding.

  • Funding The study was supported by PHS grants R01 DK45729-24 (JLM) and UACC P30 CA23074 (shared resource cores).

  • Competing interests None declared.

  • Provenance and peer review Not commissioned; externally peer reviewed.

  • Supplemental material This content has been supplied by the author(s). It has not been vetted by BMJ Publishing Group Limited (BMJ) and may not have been peer-reviewed. Any opinions or recommendations discussed are solely those of the author(s) and are not endorsed by BMJ. BMJ disclaims all liability and responsibility arising from any reliance placed on the content. Where the content includes any translated material, BMJ does not warrant the accuracy and reliability of the translations (including but not limited to local regulations, clinical guidelines, terminology, drug names and drug dosages), and is not responsible for any error and/or omissions arising from translation and adaptation or otherwise.

Request Permissions

If you wish to reuse any or all of this article please use the link below which will take you to the Copyright Clearance Center’s RightsLink service. You will be able to get a quick price and instant permission to reuse the content in many different ways.