LINE-1 retrotransposons mediate cis-acting transcriptional control in human pluripotent stem cells and regulate early brain development

Abstract

Long interspersed nuclear element 1 (L1) retrotransposons represent a vast source of genetic variability. However, mechanistic analysis of whether and how L1s contribute to human developmental programs is lacking, in part due to the challenges associated with specific profiling and manipulation of human L1 expression. Here, we show that thousands of hominoid-specific L1 integrants are expressed in human induced pluripotent stem cells and cerebral organoids. The activity levels of individual L1 promoters vary widely and correlate with an active epigenetic state. Efficient on-target CRISPR interference (CRISPRi) silencing of L1s revealed nearly a hundred co-opted L1-derived chimeric transcripts, and L1 silencing resulted in changes in neural differentiation programs and reduced cerebral organoid size. Together, these data implicate L1s and L1-derived transcripts in hominoid-specific CNS developmental processes.

Keywords: LINE-1; cerebral organoids; epigenetics; evolution; hiPSCs; neurodevelopment; non-coding genome; pluripotency; transposable elements.

Graphical abstract

Figure 1

L1s are highly expressed in hiPSCs (A) Left: schematic of the experimental design (top) to profile L1 expression in hiPSCs. Bioinformatic rationale (bottom) to quantify L1 expression by using a unique mapping approach to quantify stranded expression of unique L1 loci. Middle: overview of a full-length L1 element with an intact 5′ UTR and two long open reading frames (ORF1 and ORF2). The zoom-in shows the bidirectional promoter encoding the L1 body in sense and ORF0 in antisense. Arrows indicate the sense and antisense promoters. Right: phylogenetic tree showing the evolutionary age of different L1 subfamilies. (B) Expression (reads per kilobase per million mapped reads [RPKM]) of full-length (>6 kb), evolutionarily young L1s in two hiPSC lines. n = 3 technical replicates. (C) Number of L1s belonging to primate-specific L1 subfamilies (L1HS to L1PA4) expressed in hiPSCs. (D) Normalized read counts in antisense (red) and sense (blue) and of full-length L1HS (top) and L1PA2 (bottom) per sample (n = 3 technical replicates, t test). Bars represent mean normalized expression, error bars show mean ± standard deviation. (E) Genome browser tracks showing normalized expression (RPKM) of an L1HS element in two hiPSC lines (dark blue, forward transcription; purple, reverse transcription). (F) Western blots (WBs) of ORF1p in two hiPSC lines (top) and actin-β (bottom) as a loading control protein. (G) Immunocytochemistry of the pluripotency marker NANOG (red) and the L1-derived protein ORF1p (green) in hiPSCs. DAPI nuclear staining (blue) is included in the overlay.

Figure 2

The epigenetic profile over L1 loci correlates with their expression (A and B) Schematics of the CUT&RUN and DNA methylation analysis workflow. (C) H3K4me3 peaks over the promoter of young (L1HS-PA4) >6 kb L1s (n = 2 biological replicates) (RPKM normalized). TSS, transcription start site. Profile plot at the top shows the summed signal. (D) Violin plots showing DNA methylation levels (ONT DNA sequencing data) over the promoter of full-length (>6 kb), evolutionarily young (L1HS–L1PA4) L1 elements in two hiPSC lines, including H3K4me3 peak-called L1s. Boxplot centers correspond to the median, hinges correspond to the first and third quartile, and whiskers stretch from the first and third quartile to +1.5 interquartile range (IQR). (E) Scatterplots showing a negative correlation between normalized expression of individual full-length (>6 kb) L1s (x axis; normalized counts by gene sizeFactors as calculated by DESeq2 [i.e., median of ratios]) and the percentage of methylated CpG sites at their promoters (y axis). R² and p value of fitted linear models between methylation and expression per L1 subfamily are shown. (F) Genome browser tracks showing, from top to bottom, H3K4me3 mark at the L1 promoter (RPKM normalized) and transcription of the element (dark blue, forward transcription; purple, reverse transcription, RPKM normalized). ONT DNA reads in the region, with black dots indicating methylated CpG sites, and methylation coverage of the L1 elements are shown at the bottom (light blue, hiPSC1; dark blue, hiPSC2).

Figure 3

CRISPRi-based silencing of L1s in hiPSCs does not affect human pluripotency (A) Schematic of the gRNA target sites within the full-length L1s. The gRNAs were designed to target the 5′ UTR of evolutionarily young, >6 kb L1s. (B) Schematic of the CRISPRi workflow and downstream analyses. (C) CUT&RUN analysis of dCas9 (gRNA1, control signal) (left) and bulk RNA-seq data showing the normalized expression of uniquely mapped, full-length (>6 kb), evolutionarily young (L1HS to L1PA4) L1s in control vs. L1-CRISPRi hiPSCs (all tracks RPKM normalized). (D) H3K4me3 peaks over the promoter of young (L1HS–PA4) >6 kb L1s in control vs. L1-CRISPRi hiPSCs. Profile plot at the top shows the summed signal. (E) Genome browser tracks illustrating dCas9 enrichment (gRNA1, control signal), H3K4me3 loss over L1 elements’ promoter, and the loss of expression in control vs. L1-CRISPRi hiPSCs (lilac, dCas9 CUT&RUN signal; green, H3K4me3 CUT&RUN signal; dark blue, forward transcription; purple, reverse transcription) (all tracks RPKM normalized). (F) Expression of L1 families analyzed using TEtranscripts in control vs. L1-CRISPRi hiPSCs (heatmap showing normalized expression, scaled by row). (G) Western blots (WBs) of ORF1p (top) and actin-β (bottom) in control vs. L1-CRISPRi hiPSCs. (H) Immunostaining of the pluripotency marker NANOG (red) and the L1-derived protein ORF1p (green) in control vs. L1-CRISPRi hiPSCs. DAPI nuclear staining is in blue. (I) Mass spectrometry (MS) data showing changes (p_adj) in ORF1p levels upon L1-CRISPRi in hiPSCs. See STAR Methods for statistical analysis details. Bars show protein levels (× 1,000), and error bars correspond to ± standard deviation. (J) Left: heatmap showing log2 normalized expression (RNA-seq) of pluripotency and differentiation markers in control vs. L1-CRISPRi hiPSCs. Right: heatmap of MS data showing expression of pluripotency markers in control vs. L1-CRISPRi hiPSCs (heatmap showing log2 normalized expression).

Figure 4

L1s drive the expression of protein-coding genes and long non-coding RNAs in hiPSCs (A) Scatterplot showing mean gene expression in L1-CRISPRi (y axis) and control (x axis) hiPSCs and summary of the differential expression analysis (DEA) (bulk RNA-seq, n = 3 replicates/condition). Blue dots, significantly downregulated genes; red dots, significantly upregulated genes; gray dots, non-significant (DESeq2: Wald test; p_adj < 0.05; log2[fold change] > 1). (B) Top: distribution between protein-coding and non-coding downregulated genes. Bottom: number of downregulated genes split by their distance from the nearest full-length L1. In purple, number of genes that overlap with a full-length L1 (intragenic). (C) Heatmap showing all the normalized expression of downregulated protein-coding genes that overlap with a full-length L1 (n = 3 replicates/condition, heatmap showing log2 normalized expression). (D) Schematic showing the expression of a canonical transcriptional isoform from the canonical gene promoter vs. expression of an L1-driven alternative transcriptional isoform. (E) Genome browser tracks showing the expression of L1-driven alternative transcriptional isoforms at the ELAPOR2 (top) and PPP1R1C (bottom) loci in hiPSCs and its silencing upon L1-CRISPRi. Left shows normalized transcription (RPKM) in L1-CRISPRi and control hiPSCs and ONT long-read direct RNA from control hiPSCs showing the different isoforms. Right shows a zoom-in to the L1 elements showing, from top to bottom, bulk RNA-seq data (dark blue, forward transcription; purple, reverse transcription), H3K4me3 CUT&RUN (green, control vs. L1-CRISPRi hiPSCs), and ONT direct RNA reads from control hiPSCs. (F) Normalized exon expression between control and L1-CRISPRi hiPSCs of the first three and last three exons of ELAPOR2 and all the exons of PPP1R1C (p_adj DESeq2: Wald test, n = 3 replicates per condition). Intronic L1 position is indicated in blue. Bars represent mean normalized expression, and error bars show mean ± standard error.

Figure 5

Expression of evolutionarily young full-length L1s in cerebral organoids (A) Schematic of the workflow for the generation of unguided cerebral organoids and downstream analyses. (B) Bright-field images of differentiating cerebral organoids at different time points and immunohistochemistry on day 15 organoids for ZO1 (red) and PAX6 (green). DAPI is included in blue in the overlay. (C) Left: uniform manifold approximation and projection (UMAP) showing clusters found in day 15 unguided cerebral organoids. Middle, top: UMAP displaying the different cell types found in day 15 cerebral organoids. Middle, bottom: bar plot showing the percentage cell type composition of the cerebral organoids at day 15 of differentiation. Right: dot plot displaying selected neuronal and NPC markers used to characterize the cell clusters (dot size shows the percentage of cells expressing the gene, color indicates average expression in each cell type). (D) Expression of uniquely mapped evolutionarily young full-length L1s in hiPSCs and day 15 cerebral organoids. (E) Violin plots of the methylation status over the promoter of L1HS-L1PA4 >6 kb L1s in hiPSCs vs. day 15 cerebral organoids. Zoom-in plots indicate mean methylation levels (red dot) per condition per subfamily. Boxplot centers correspond to median, hinges correspond to the first and third quartile, and whiskers stretch from the first and third quartile to +1.5 IQR. (F) Genome browser tracks showing normalized transcription (RPKM) of L1s in hiPSCs and day 15 cerebral organoids (dark blue, forward transcription; purple, reverse transcription) and ONT DNA reads across day 15 organoids and hiPSCs, black dots indicating methylated CpGs, and methylation coverage of the L1 elements at the bottom (green, organoids; blue, hiPSCs).

Figure 6

L1s regulate early neural differentiation in cerebral organoids (A) Workflow for the differentiation of L1-CRISPRi and control hiPSCs into day 15 unguided cerebral organoids and downstream analyses. (B) Differentiation of the unguided cerebral organoids (hiPSC1 control vs. L1-CRISPRi). Bright-field images included from days 4, 9, and 14 of differentiation. Right: immunostaining showing neural rosettes and expression of the tight junction marker ZO1 (red) and the neural progenitor cell marker PAX6 (green). Nuclear staining DAPI is included in the overlay (blue). (C) Normalized expression (RPKM) of uniquely mapped evolutionarily young (L1HS–L1PA4), >6 kb L1s in control vs. L1-CRISPRi unguided cerebral organoids at day 15 of differentiation. (D) Left: UMAP showing clustering of day 15 organoids (control + L1-CRISPRi). Middle: UMAP colored by cell-cycle score (S + G2M). Right: UMAPs showing the identified cell types in control and L1-CRISPRi day 15 cerebral organoids. (E) Dot plot showing the expression of selected markers used for the characterization of the cell clusters (dot size shows the percentage of cells expressing the gene, color indicates average expression in cell type per condition). (F) Bar plot of the cell type distribution across batches and cell lines of control vs. L1-CRISPRi organoids. (G) Left: heatmap showing commonly downregulated genes (L1-CRISPRi vs. control organoids) in all cell lines and guides (genes selected based on hiPSC1 gRNA1 avg log2FC < −0.25, hiPSC1 gRNA2, and hiPSC2 gRNA1 avg log2FC < 0) in the snRNA-seq data. The full list of downregulated genes can be found in Table S1. Right: heatmap with all the upregulated genes found in all cell lines and guides (genes selected based on hiPSC1 gRNA1 avg log2FC > 0.25, hiPSC1 gRNA2, and hiPSC2 gRNA1 avg log2FC > 0). The average expression of the upregulated genes in NPCs across conditions is also shown. (H) Violin plots of selected genes upregulated in control (gray) compared to L1-CRISPRi (purple) organoids in cluster 0 and cluster 1 (Seurat::FindMarkers, Wilcoxon test; number of nuclei per condition as in Table S2). (I) Gene set enrichment analysis of upregulated genes from clusters 0 and 1 (control vs. L1-CRISPRi organoids). (J) Growth curves showing the size measurement distribution of control vs. L1-CRISPRi cerebral organoids from day 2 to day 15 of differentiation. Control is gray, L1-CRISPRi is purple. Area was measured as square micrometers and was assessed with Fiji ImageJ. Student’s t test was performed to compare the two conditions at each time point.