Comprehensive molecular portrait reveals genetic diversity and distinct molecular subtypes of small intestinal neuroendocrine tumors

Abstract

Small intestinal neuroendocrine tumors (siNETs) are rare bowel tumors arising from malignant enteroendocrine cells, which normally regulate digestion throughout the intestine. Though infrequent, their incidence is rising through better diagnosis, fostering research into their origin and treatment. To date, siNETs are considered to be a single entity and are clinically treated as such. Here, by performing a multi-omics analysis of siNETs, we unveil four distinct molecular groups with strong clinical relevance and provide a resource to study their origin and clinical features. Transcriptomic, genetic and DNA methylation profiles identify two groups linked to distinct enteroendocrine differentiation patterns, another with a strong immune phenotype, and the last with mesenchymal properties. This latter subtype displays the worst prognosis and resistance to treatments in line with infiltration of cancer-associated fibroblasts. These data provide insights into the origin and diversity of these rare diseases, in the hope of improving clinical research into their management.

Fig. 1. Clinical description of siNET pathology.

A valuable resource for the study of siNETs. a Schematic representation of the siNET cohort. SNP (Single Nucleotide Polymorphism); FFPE (Formalin-Fixed Paraffin-Embedded). Created with BioRender.com. b Representative hematoxylin and eosin (H&E) and Chromogranin-A immunohistochemistry (IHC) staining for pT1, pT2 and pT3 statuses. All 122 cases have undergone H&E staining and Chromogranin-A IHC for diagnosis. Scale bar 1 mm. c Table of the clinical characteristics of the siNET cohort. WHO (World Health Organization); G (grade); ENETS (European Neuroendocrine Tumor Society); UICC (Unio Internationalis Contra Cancrum); LN (Lymph Node); SRI (Somatostatin Receptor Imaging); CI (Confidence Interval); NR (Not Reached) d left panels: representation of four out twenty clinical variables used for the multiple correspondence analysis (MCA) in the first two planes (functionality, uncontrolled carcinoid syndrome, Chromogranin A (CgA) dosage at surgery time and number of metastatic sites variables) (see Supplementary Fig.1). Right panels: Kaplan-Meier, log-rank test and Cox proportional hazards regression model methods were used to study overall (OS) and progression-free (PFS) survival for the four selected clinical features. The shaded areas represent 95% confidence intervals around the curves; p-values are indicated below the graphs.

Fig. 2. Transcriptional landscape of siNETs defines four molecular features and the mesenchymal character predicts poor prognosis.

Exploring the transcriptional landscape of siNETs by RNA-sequencing. a Determination of the optimal number of gene clusters in the siNETs cohort with consensus clustering based on the 10% most variant genes among the discovery cohort of siNETs. CDF (Cumulative Distribution Function) b Leiden clustering analysis (w/ 4 clusters) projected in UMAP-2D plan. c Pathway enrichment analysis of the 4 siNETs cluster described in (a, b) (epithelial in green, vesicular in red, immune in blue and mesenchymal in gold). Adjusted two-tailed p-values (hypergeometric test) are shown on each bar (BH: Benjamini-Hochberg). d Hierarchical clustering (clustering method: Ward’s; distance: Euclidean): genes (row) are split according to the Leiden gene expression clusters: epithelial (green), vesicular (red), immune (blue), mesenchymal (gold). Recapitulative ssGSEA score of each cluster is indicated as bottom annotation (see Supplementary Fig. 2). TPM (Transcripts Per Million). e Kaplan-Meier, log-rank test and Cox proportional hazards regression model methods were used to study overall (OS) survival for epithelial, vesicular, immune, mesenchymal statuses respectively in green, red, blue and gold. The shaded areas represent 95% confidence intervals around the curves; the median survival rate is represented by dotted line; p-values are indicated below the graphs.

Fig. 3. siNETS subtypes arise from different types of neuroendocrine differentiation.

siNETs tumors are defined by different precursors and/or types of differentiation. a Unsupervised clustering of ssGSEA scores for fetal small intestine single-cell RNASeq-based cell type signatures in the full siNETs cohort (n = 206). b Unsupervised clustering of gene expression of transcription factors, described as essential for intestinal differentiation, across the siNETs cohort. ATOH1 (Atonal BHLH Transcription Factor 1); ARX (Aristaless Related Homeobox); DCLK1 (Doublecortin Like Kinase 1); LMX1A (LIM Homeobox Transcription Factor 1 Alpha); ssGSEA (Single Sample Gene Set Enrichment Analysis). c, d Unsupervised clusterings of gene expression of endocrine hormones or converting enzymes: full list (c) or short list (d) across the entire siNETs cohort. CCK (Cholecystokinin); GAST (Gastrin); GCG (Pro-glucagon: Glucagon-like peptide 1 et 2, Oxyntomodulin); GHRL (Ghrelin); GIP (Gastric inhibitory polypeptide); HDC (Histidine decarboxylase); LEP (Leptin); MLN (Motilin); NTS (Neurotensin); NUCB2 (Nesfatin-1); PCSK1 (prohormone convertase 1 for glucagon-like peptide 1 et 2, oxyntomodulin); PYY (Peptide YY); SCTR (Secretin); SI (Small Intestine); SST (Somatostatin); TPH1 (Tryptophan Hydroxylase 1, biosynthesis of serotonin). Clustering method: Ward’s; distance: Spearman. Recapitulative ssGSEA score of each gene cluster (epithelial, vesicular, mesenchymal and immune) is indicated as bottom annotation. LN (Lymph Node); MN (Mesenteric Node); T (Tumor); TPM (Transcripts Per Million).

Fig. 4. Genomic landscape mutations of siNETs.

Human siNETs are poorly mutated without clear genetic drivers. OncoPrint of DNA somatic alterations from 24 tumors sequenced in whole genome. The percentage on the right indicates the mutation frequency of each gene across samples. OncoPrint is split in 3 panels: frequently altered genes across WGS cohort (top), frequently altered genes among TCGA cancer pathways (middle), frequently altered whole chromosomes (bottom). Mutation counts (substitutions and small indels) and COSMIC mutational signatures distribution per sample were inferred from whole genome data. Phylogenetic relationships between primary and metastasis samples are indicated with a solid line above the tissue annotation track. Tumor purity, ploidy and fraction of genome altered (FGA) were estimated from FACETS analysis. LOH (Loss Of Heterozygosity); MN (Mesenteric Node); T (Tumor).

Fig. 5. Chromosomal rearrangements delineate distinct profiles and their association with siNETs survival.

siNETs can be defined by specific chromosomic rearrangements. a Frequency plots of DNA gains (yellow) and losses (blue) in the discovery cohort. RB1 (Retinoblastoma 1). b Unsupervised hierarchical clustering of copy number alterations in the ileum discovery cohort. A 20260 genes x 99 samples matrix encoded with −1 (loss), 0 (no alteration), +1 (gain) was used (see Methods). Tumor purity, ploidy and fraction of genome altered (FGA) were estimated from ASCAT analysis. Recapitulative ssGSEA score of each gene cluster (epithelial, vesicular, mesenchymal and immune) is indicated as bottom annotation. CNA (Copy Number Alterations). Clustering method: Ward’s; distance: binary. c, d Kaplan-Meier, log-rank test and Cox proportional hazards regression model methods were used to study overall (OS) and progression-free (PFS) survival for chr18.del and chr4.10.14.gains statuses. chr4.10.14.gains was set to YES if at least one of chromosomes 4, 10 or 14 was gained, NO otherwise. The shaded areas represent 95% confidence intervals around the curves; p-values are indicated below the graphs.

Fig. 6. The DNA methylation landscape of siNETs.

Methylation refines the siNET groups. a Unsupervised k-means consensus clustering was found to be optimal with 4 groups for the top 1,000 most variable CpGs methylation probes (in row: A, B, C and D) and samples (in column, Sample.Methyl.Clust: epithelial-enriched, hypomethylated, unifocal-enriched and FGA-enriched). Clustering method: Ward’s; distance: Euclidean. Recapitulative ssGSEA score of each gene cluster (epithelial, vesicular, mesenchymal and immune), chr18.del and chr4.10.14.gains statuses together with methylation subtypes are indicated as bottom annotations. b Boxplots of the methylation levels (mean beta value) of the 1000 most variants probes (left) and each of the four probe clusters (right: A-D) across the four k-means based methylation clusters (hypomethylated, n = 41; epithelial-enriched, n = 18; unifocal-enriched, n = 26 and FGA-enriched, n = 22). Two-tailed exact p-values were determined by Kruskal-Wallis test. Boxplots: center line = median, box range 25th–75th percentile, minimum/maximum denoted by whiskers. c Barplots for co-occurrence analysis of methylation clusters vs expression clusters, FGA levels, chr18 deletion, chr4.10.14 gains and tumor type (unifocal or multifocal). Significance was determined by two-tailed Fisher’s exact test.

Fig. 7. Multifocal and unifocal forms belong to two different pathologies.

Omics can define unifocal vs multifocal forms of siNETs a Venn-diagram of somatic point mutations and small indels on the whole-genome scale for available WGS samples of multifocal patients P101, P17, P41 and P68 (patient P68 was analyzed using a tumor as a paired sample). b Volcano plot for multifocal (n = 45) vs unifocal (n = 57) differential gene expression (DE) analysis; p-values are indicated on the y-axis (two-tailed Wald test with Benjamini-Hochberg correction). c Hierarchical clustering based on the 9 highly significantly differentially expressed coding genes ( | log2(FoldChange)|>1 & pval.adj < 1e-2). ATP8B3 (ATPase Phospholipid Transporting 8B3); CAMKK1 (Calcium/Calmodulin Dependent Protein Kinase Kinase 1); HOXC10 (Homeobox C10); LIPF (Lipase F, Gastric Type); NCR1 (Natural Cytotoxicity Triggering Receptor 1); NWD2 (NACHT And WD Repeat Domain Containing 2); PIK3C2G (Phosphatidylinositol-4-Phosphate 3-Kinase Catalytic Subunit Type 2 Gamma); PLA2G2C (Phospholipase A2 Group IIC); SEMA3E (Semaphorin 3E). d Volcano plot for multifocal (n = 44) vs unifocal (n = 54) differential methylation (DM) analysis; p-values are indicated on the y-axis (two-tailed T-test with Benjamini-Hochberg correction). e Hierarchical clustering based on the 836 significantly differentially methylated probes (delta(Beta)>10%, pval.adj < 5e-2). f Multifocal (n = 45) vs unifocal (n = 57) boxplots of gene expression for endocrine hormones or converting enzymes (log2TPM). GAST (Gastrin); GCG (Pro-glucagon: Glucagon-like peptide 1 et 2, Oxyntomodulin); GIP (Gastric inhibitory polypeptide); HDC (Histidine decarboxylase); NTS (Neurotensin); SST (Somatostatin). Boxplots: center line = median, box range 25th–75th percentile, minimum/maximum denoted by whiskers. Significance was determined by Mann-Whitney U tests. Clustering method: Ward’s; distance: Spearman. Recapitulative ssGSEA score of each gene cluster (epithelial, vesicular, mesenchymal and immune), chr18.del and chr4.10.14.gains statuses together with methylation subtypes are indicated as bottom annotations. FGA (Fraction genome altered); TPM (Transcripts Per Million).

Fig. 8. Cancer-associated fibroblasts induce resistance to treatments, enhance cancer cell proliferation and decrease survival.

CAFs favor severe forms of siNETs. a Unsupervised clustering of MCPcounter scores estimating the abundance of immune and stromal cells infiltrate across the full siNETs cohort (n = 206). Clustering method: Ward’s; distance: Spearman. X-cell derived Immune and Stromal Scores, recapitulative ssGSEA score of each gene cluster (epithelial, vesicular, mesenchymal and immune), chr18.del and chr4.10.14.gains statuses together with methylation subtypes are indicated as bottom annotations. LN (Lymph Node); MN (Mesenteric Node); T (Tumor). Clustering method: Ward’s; distance: Spearman. b Representative anti-alpha-smooth muscle actin (αSMA) immunohistochemistry (IHC) staining illustrating CAF infiltration of siNETs (low CAF score vs high CAF score). Twenty-one tumors have undergone αSMA IHC. Scale bar 100 µm. c Prognostic value for OS (Overall Survival) and PFS (Progression Free Survival) of the CAFs score in siNETs tumors. A cut-off in EPIC-related CAF scores distribution was applied (score=0.025) to label samples either as low (below cut-off) or high (above cut-off) for CAF infiltration. Kaplan-Meier, log-rank test and Cox proportional hazards regression model methods were used to study overall (left) and progression-free (right) survival. The shaded areas represent 95% confidence intervals around the curves; the median survival rate is represented by dotted line; p-values are indicated below the graphs. d Schematic representation of purification, culture and production of conditioned media from patients’ siNETs associated fibroblasts and transfer to the GOT-1 cell line. Created with BioRender.com. e Analysis of 2-months proliferation of GOT-1 cells treated with CAF conditioned medium (p = 0.0286). Two-tailed exact p-value was calculated using a Mann-Whitney U test. Error bars represent means ± SEM, n = 4. Source data are provided as a Source Data file. RLU (Relative Light Units).