Characterization of enhancer activity in early human neurodevelopment using Massively parallel reporter assay (MPRA) and forebrain organoids

Abstract

Regulation of gene expression through enhancers is one of the major processes shaping the structure and function of the human brain during development. High-throughput assays have predicted thousands of enhancers involved in neurodevelopment, and confirming their activity through orthogonal functional assays is crucial. Here, we utilized Massively Parallel Reporter Assays (MPRAs) in stem cells and forebrain organoids to evaluate the activity of ~7,000 gene-linked enhancers previously identified in human fetal tissues and brain organoids. We used a Gaussian mixture model to evaluate the contribution of background noise in the measured activity signal to confirm the activity of ~35% of the tested enhancers, with most showing temporal-specific activity, suggesting their evolving role in neurodevelopment. The temporal specificity was further supported by the correlation of activity with gene expression. Our findings provide a valuable gene regulatory resource to the scientific community.

Fig 1.. Experimental design and overall lentiMPRA results.

A. ~7,000 candidate enhancer regions were selected among the 96,000 gene-linked enhancers identified by Amiri et al [23]. The lentiMPRA library was synthesized on a custom array and cloned into a lentiMPRA vector, packaged into lentivirus, and introduced into 3 iPSC lines before organoid differentiation. B. At 3 different timepoints, DNA/RNA sequencing was used for estimating the enhancer activity from the ratio of corresponding barcodes. Enhancers were considered active by MPRA if their activity was significantly above the background (blue curve) derived from a gaussian mixture model of the activity of negative controls. C. The Venn diagram shows counts of MPRA-active enhancers across time points. BC, barcode; LTR, long terminal repeat; mP, minimal promoter; TD, terminal differentiation.

Fig 2.. Characteristics of active enhancers.

A. Heatmap of enhancer activity (RNA/DNA) of all candidate enhancers tested in MPRA experiment. Clustering was performed using Ward variance minimization algorithm. Right bar annotates MPRA-active (red) or inactive (blue) enhancers. B. Upset plot showing the number of active (red) and inactive (blue) enhancers overlapping external datasets. In the box, bars represent the fraction of overlapping active and inactive enhancers relative to the total number of active and inactive enhancers. P-values were calculated using Fisher’s exact test (*p-value < 0.05; **p-value < 0.01). C. Cumulative plot of number of transcription factor binding sites identified by FIMO in positive controls, negative controls, and tested enhancers. No statistical significantly difference was observed between active and inactive enhancers. D. Cumulative plot of the expression (RPKM) of genes predicted to be regulated by tested enhancers in TD0 and TD30 organoids in the endogenous genomic context by Amiri et al. [23]. Genes regulated by MPRA-active enhancers have significantly higher expression than inactive enhancers. The background curve represents the expression of all ~96,375 gene-linked enhancers from Amiri et al [23].

Fig 3.. MPRA active enhancers correlate with differences in expression of linked genes across different time points.

A. Scatter plots of change in MPRA activity for an MPRA-active enhancer (y-axis) and change of expression for the gene(s) linked to the enhancer in the endogenous genomic context. Activity and expression are compared for TD0 and iPSC. Each dot represents a gene-enhancer pair. iPSC or TD0 active enhancers are shown as red circles in the upper panel; inactive enhancers are represented in the bottom panel as blue circles. B. More correlated enhancer-gene pairs are observed in MPRA-active enhancers compared to inactive enhancers. The x-axis represents the cutoff for the change in enhancer activity and gene expression. The y-axis represents the fraction of enhancer-gene pairs passing cut off for the change in the enhancer activity and gene expression. The dashed line represents cutoff at log2FC=0.8, which is the same as marked in panel (A). C. Example of an iPSC-only active enhancer (MPRA activity shown in the heatmap) and the expression of the linked gene (G3BP2), with ZNF263 as the predicted binding transcription factor. D. Example of a TD0 active enhancer (MPRA activity shown in the heatmap) and the expression of the linked gene (CUX1) with ASCL1 as predicted binding transcription factor. For comparison, a mildly active enhancer at TD0 for this gene is identified downstream of the active enhancer. TSS, Transcription Start Site.