Human Pancreatic β Cell lncRNAs Control Cell-Specific Regulatory Networks

Abstract

Recent studies have uncovered thousands of long non-coding RNAs (lncRNAs) in human pancreatic β cells. β cell lncRNAs are often cell type specific and exhibit dynamic regulation during differentiation or upon changing glucose concentrations. Although these features hint at a role of lncRNAs in β cell gene regulation and diabetes, the function of β cell lncRNAs remains largely unknown. In this study, we investigated the function of β cell-specific lncRNAs and transcription factors using transcript knockdowns and co-expression network analysis. This revealed lncRNAs that function in concert with transcription factors to regulate β cell-specific transcriptional networks. We further demonstrate that the lncRNA PLUTO affects local 3D chromatin structure and transcription of PDX1, encoding a key β cell transcription factor, and that both PLUTO and PDX1 are downregulated in islets from donors with type 2 diabetes or impaired glucose tolerance. These results implicate lncRNAs in the regulation of β cell-specific transcription factor networks.

Keywords: CRISPR interference; PDX1; PLUTO; chromatin; diabetes; lncRNAs; long noncoding RNAs; pancreatic islets; transcriptional networks; type 2 diabetes.

Graphical abstract

Figure 1

Knockdown of Selected β Cell lncRNAs Leads to Transcriptional Phenotypes (A) Schematic of the experimental plan. Lentivirally encoded amiRNAs were validated and transduced in duplicate (×2) or triplicate (×3) into ENDOC-βH1 cells as indicated and then analyzed with oligonucleotide expression arrays. (B) Differential gene expression analysis revealed genes that show significant up- or downregulation after knockdown of TFs or lncRNAs. For each TF or lncRNA, we combined all replicates transduced with the different target-specific amiRNAs and compared these with all replicates from five non-targeting controls. Differential expression was determined at p < 10⁻³ (ANOVA). (C) We compared gene expression data from all ten possible combinations of three versus two control non-targeting amiRNAs. Similarly, the two independent amiRNAs that target each TF or lncRNA were compared with all ten possible combinations of three control amiRNAs. For this analysis, we only considered the seven lncRNAs that were targeted by two independent amiRNAs. (D) Control comparisons result in a low number of differentially regulated genes (average 15 genes), whereas most TF and lncRNA comparisons yield higher numbers of differentially regulated genes. ^∗∗∗p < 10⁻⁴; ^∗∗p < 0.01; ns, not significant compared with control comparisons; Mann-Whitney test. (E) HI-LNC15 regulates its neighboring gene, NKX2.2, whereas HI-LNC12 knockdown (KD) does not affect its adjacent active gene, UNC5A (left). Further examples are shown in Figure S1G. RNAs were normalized to TBP mRNA and expressed relative to control amiRNAs; n = 3, error bars represent SEM; ^∗∗p < 0.01, ^∗p < 0.05 (Student’s t test).

Figure 2

Knockdown of lncRNAs Impairs Insulin Secretion (A) Examples of genes known to play a role in β cell function regulated by islet lncRNAs. (B) Glucose-stimulated insulin secretion was tested on T antigen-excised EndoC-βH3 cells after transduction with amiRNAs targeting the indicated lncRNAs or controls. Secreted or total insulin content was normalized to the number of cells per well and expressed as fold change over control amiRNA treatment at 2.8 mM glucose. Each bar represents an average from two independent amiRNA vectors and 12 separate wells from two independent experiments. Error bars represent SEM; ^∗∗∗p < 10⁻³, ^∗∗p < 0.01, ^∗p < 0.05 (Student’s t test).

Figure 3

Human Islet TFs and lncRNAs Regulate Common Genes (A) Heatmap displaying Pearson r values for all pairwise comparisons of fold changes in gene expression after knockdown of TFs and lncRNAs. Only genes significantly dysregulated at p < 10⁻³ under at least one condition were included in the analysis. (B) Unsupervised clustering analysis of fold change values after knockdown of five TFs and the five lncRNAs that displayed the strongest transcriptional changes. Only genes that were dysregulated at p <10⁻³ in at least one knockdown were selected. Blue represents downregulated and red represents upregulated genes. Controls represent control comparisons as described for Figure 1. (C) Examples of highly correlated transcriptional phenotypes. The plots show fold change values (Log2) after knockdown of the indicated pairs of genes. Only the top 100 most regulated genes for any of the two knockdowns were plotted. Pearson's correlation (r) and p values are displayed.

Figure 4

LncRNAs Regulate Enhancer Cluster Genes GSEA showed that genes that were downregulated upon knockdown of either islet TFs or lncRNAs were enriched in a set of 694 genes that is associated with human islet enhancer clusters (red dots) but not in ten control gene sets (black dots) that were expressed at similar levels as enhancer cluster genes.

Figure 5

Islet-Specific Coding and Noncoding RNAs Form Shared Co-expression Modules (A) Topological overlap matrix representing co-expression modules that were co-regulated across 64 human islet samples. Modules that were enriched in lncRNAs are marked with squares (hypergeometric test, p < 10⁻²). (B–D) Co-expression modules that showed enrichment in islet lncRNAs (B), islet enhancer cluster (EC)-associated genes (C), or a set of 94 islet-enriched TF genes (D). Five modules (M3, M7, M12, M18, and M20, marked in blue) out of seven modules that were enriched in lncRNAs were also enriched in ECs and TFs. (E) Module M3 was enriched in typical islet-specific biological process annotations. Right: examples of islet TFs and lncRNAs in module M3. (F) Correlation of the indicated lncRNAs and β cell-specific TF mRNAs across 64 islet samples. GAPDH is shown as a non-β cell reference. Pearson’s correlation values are displayed in the top left corner. The axes show expression values normalized across 64 islet samples. (G) Network diagram illustrating that TFs and lncRNAs often co-regulate the same genes, many of which are associated with enhancer clusters.

Figure 6

PLUTO Knockdown Decreases PDX1 mRNA (A) Downregulation of PLUTO (HI-LNC71) and PDX1 in islets from donors with T2D or IGT. Differential expression analysis was performed on control (n = 50) versus T2D (n = 10) or IGT (n = 15) samples. Boxplots represent expression normalized to the mean of control samples. Adjusted p values are shown. (B) Schematic of the human PDX1 locus and its associated enhancer cluster. A 4C-seq analysis was designed to identify regions interacting with the PDX1 promoter region in EndoC-βH1 cells. Red and orange vertical lines depict active and poised islet enhancers, respectively. F and R represent forward and reverse RNA-seq strands, respectively, and scales represent RPM. PLUTO (HI-LNC71) was generated from a de novo assembly of islet RNA-seq and differs from a transcript annotated in UCSC and RefSeq that originates from a PDX1 intronic region. (C) Downregulation of PLUTO or PDX1 using amiRNAs resulted in reduced PDX1 mRNA and protein levels. EndoC-βH1 cells were transduced with control (black), PLUTO (white), or PDX1 (turquoise) amiRNA vectors 80 hr prior to harvest. RNA levels were assessed by qPCR, normalized to TBP, and expressed as fold over control amiRNA samples (n = 4). For protein quantification, PDX1 levels were first normalized to the average of TBP and H3 levels and then compared with the control amiRNA sample. (D) Downregulation of PLUTO in human islet cells results in reduced PDX1 mRNA levels. Islet cells were dispersed and transduced with amiRNA vectors (n = 3) as in (B). (E) Downregulation of PLUTO in EndoC-βH3 cells using CRISPRi also decreases PDX1 mRNA. EndoC-βH3 cells were nucleofected with CRISPRi vectors 80 hr prior to harvest. RNA levels were assessed by qPCR and normalized to TBP and then to a control CRISPRi sample (n = 3). (F) PDX1 and PLUTO RNA levels were highly correlated in 64 human islet samples. (G) Knockdown of PDX1 and PLUTO resulted in differential expression of similar genes. Fold change value (Log2) of top 250 dysregulated genes following the PDX1 knockdown was plotted against the same genes following the PLUTO knockdown. (H) GSEA showed that genes that were downregulated upon knockdown of PDX1 and PLUTO were enriched in genes whose enhancers were bound by PDX1 (red) in islets but not in ten control gene sets (black) that were expressed at similar levels as PDX1-bound genes. (I) Knockdown of PDX1 and PLUTO resulted in differential expression of genes with similar biological process annotations. (J) Examples of known PDX1-regulated genes that are also co-regulated by PLUTO in parallel knockdown experiments. mRNA levels were assessed as in (B). Error bars denote SEM; ^∗∗∗p < 10⁻³, ^∗∗p < 0.01, ^∗p < 0.05 (Student’s t test).

Figure 7

PLUTO Regulates PDX1 Transcription and 3D Chromatin Structure (A) The mRNA stability of PDX1 was unaffected by PLUTO knockdown. PDX1 mRNA was measured in control and PLUTO amiRNA knockdown in EndoC-βH1 cells after Actinomycin D (ActD) treatment (n = 3). mRNA levels are presented as a percentage of levels observed at time = 0. (B) Knockdown of PLUTO was carried out as in Figure 6B, and this led to reduced PDX1 transcription, as assessed by qPCR analysis of intronic PDX1 RNA levels using hydrolysis probes. Values were normalized to TBP mRNA and expressed as fold over the control amiRNA sample (n = 4). (C) Schematic of selected epigenomic features of the PDX1 locus. (D) PLUTO is required for 3D contacts between the PDX1 promoter and distal enhancers. 3C analysis revealed that knockdown of PLUTO resulted in reduced contacts between the PDX1 promoter (anchor) and two enhancers (E1 and E2). Interaction signals were normalized to a control region on the PDX1 intron. CTL represents a negative control region that does not harbor interactions with the PDX1 promoter. Error bars denote ± SEM, and p values are from a Student’s t test. (E) PLUTO knockdown resulted in impaired 3D contacts between the PDX1 promoter and its adjacent enhancer cluster, causing reduced PDX1 transcriptional activity.