Establishment and genetic analysis of Gtl2 transcription termination mouse

Gtl2 is located within the Dlk1-Dio3 imprinting domain at the end of mouse chromosome 12 (Fig. 1a). Previous analysis of the Gtl2 locus has revealed a differentially methylated region (Gtl2-DMR), a microRNA (miR-1906), two CTCF binding sites and some active open chromatin regions between exons 1–5 [12, 33]. In order to silence the Gtl2 gene without affecting the action and activity of other transcription factors, we identified a highly efficient sgRNA site at the end of the first exon, and used 96 bp of each upstream and downstream sgRNA as the left and right homology arms to form ssDNA with the inserted 3 × polyA (Fig. 1b and c). The combined ssDNA and RNP complex (sgRNA and Cas9 protein) were microinjected into the male pronucleus of the zygote (Fig. S1a), cultured in vitro until the blastocyst stage (Fig. S1b), and the embryos were transferred into two surrogate mice [25]. Among the six newborn mice (Fig. S1c), five founder mice with successful 3 × polyA knock-in were identified by PCR (Fig. 1b and S1d). Taken, we successfully generated a mutant mouse model that only silenced the expression of ncRNA in Dlk1-Dio3 domain without deleting other active elements.

Fig. 1

Generation and phenotypic analysis of founder mice. a Schematic representation of Dlk1-Dio3 domain in chromosome12. Maternally expressed genes are in red, paternally expressed genes are in blue and silenced genes are in gray, miRNAs are shown as short vertical lines, differentially methylated regions are shown as circles. Filled circle, methylated; open circle, unmethylated. b Schematic representation of mutant mouse production by Easi-CRISPR: preparation of RNP components and ssDNA, microinjection, generation of founder offspring, PCR genotyping and off-target detection and obtaining of founder mice. c Schematic representation of 3 × polyA termination sequence insertion sites. d Survival and death number of embryos from E12.5 to E16.5. e Representative pictures of embryo at E12.5, E14.5 and E16.5, abdominal hemorrhage are indicated by yellow arrows, dead embryos are in red boxes. Scale bar, 2 mm. f Bar chat showing embryo weight of E12.5, E14.5 and E16.5, each dot represents an embryo. Error bars, ± SEM, n.s., p > 0.05

Since Gtl2 is a maternally expressed gene, we maintained the generated mice by crossing wild type (Gtl2+/+) or Gtl2 transcription termination occurs on the paternal alleles (Gtl2+/polyA) of male mice with Gtl2+/+ or Gtl2+/polyA females (Fig. S1e). We selected 3–4 generations of mice for the parents of embryos obtained. Mice developed normally in Gtl2+/polyA mice, while embryonic death began at E13.5 and all mice died at E16.5 (Fig. 1e) when Gtl2 transcription termination occurred on the maternal allele (Gtl2polyA/+) and the homozygous biallelic (Gtl2polyA/polyA) (Fig. S1f). The mortality rate of Gtl2polyA/+ was 0%, 21.43%, 47.06%, 77.78% and 100% from E12.5 to E16.5, respectively, and then the mortality rate of Gtl2polyA/polyA was 0%, 0%, 58.06%, 100% and 100% form E12.5 to E16.5, respectively (Fig. 1d and S1g). So far, we have found that there is a slight difference in the time of death between Gtl2polyA/+ and Gtl2polyA/polyA. According to the representative images of embryos, we observed that there was no significant difference in the four-genotype embryos at E12.5. There was no obvious change in Gtl2+/polyA, however, the Gtl2polyA/+ and Gtl2polyA/polyA mice showed abdominal hemorrhage at E14.5 (Fig. 1e). There was no obvious change in Gtl2+/polyA at E16.5, whereas the Gtl2polyA/+ and Gtl2polyA/polyA mice were all dead at this stage (Fig. 1e). Gtl2polyA/+ and Gtl2polyA/polyA were slightly but not significantly heavier than Gtl2+/+, and Gtl2+/polyA was almost the same compared to Gtl2+/+, and the trends of embryo weight change were the same at E12.5, E14.5 and E16.5 (Fig. 1f). Together, these data show that the silencing of ncRNAs on Dlk1-Dio3 domain results in little effect on Gtl2+/polyA, and a slight increase in embryo weight, abdominal hemorrhage and embryonic death on Gtl2polyA/+ and Gtl2polyA/polyA.

Termination of Gtl2 transcription resulted in silencing of the maternally expressed RNAs in the Dlk1-Dio3 domain

We next attempt to elucidate the molecular mechanisms responsible for the above results. To this end, we performed whole transcriptome sequencing in four-genotype embryos at E12.5 and detected 29,158 lncRNAs and 2395 miRNAs (Tab. S2). After excluding genes with that FPKM is 0 in all the four genotypes, we screened out differential expression lncRNAs (DElncRNAs) (Tab. S3) and differential expression miRNAs (DEmiRNAs) (Tab. S4), fold change more than 2.0, and p-value less than 0.05.

Firstly, we identified 2272 upregulated and 2454 downregulated lncRNAs in Gtl2+/polyA, 2326 upregulated and 3878 downregulated lncRNAs in Gtl2polyA/+ and 2236 upregulated and 2803 downregulated lncRNAs in Gtl2polyA/polyA (Fig. 2a) (Tab. S3). RNA-seq data of all lncRNAs showed that Gtl2polyA/+and Gtl2polyA/polyA had a significant downward trend in Dlk1-Dio3 domain (Fig. 2b), mainly in Gtl2, Rtl1as (miR-127), Rian and Mirg (Fig. 2c). The results of qRT-PCR were consistent with the trend of RNA-seq (Fig. 2d), indicating that the termination of Gtl2 transcription led to the silencing of lncRNA in the domain.

Fig. 2

Termination of Gtl2 transcription resulted in silencing of the maternally expressed RNA in the Dlk1-Dio3 domain. a Volcano plot of DElncRNAs in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5. b lncRNA-seq of embryos at E12.5. c The expression of Gtl2, Rtl1as, Rian and Mirg in RNA-seq data. d qRT-PCR analysis of Gtl2, Rtl1as, Rian and Mirg, β-actin was used as an internal control, n = 6 Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01; ***, p < 0.001. e Volcano plot of DEmiRNAs in Gtl2polyA/+ and Gtl2polyA/polyA. f miRNA-seq of embryos at E12.5. g The expression of miRNAs in Dlk1-Dio3 domain in RNA-seq data. h qRT-PCR analysis of miRNAs in Gtl2, Rtl1as, Rian and Mirg locus, U6 was used as an internal control, n = 5. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01; ***, p < 0.001

The lncRNAs and miRNAs in the Dlk1-Dio3 domain constitute a large polycistronic transcription unit, and maternally expressed intergenic transcripts downstream of Gtl2 are transcribed in the same direction as Gtl2, and there is no typical promoter sequence present in this region [5]. Therefore, we guessed that miRNAs in the domain would also be severely affected by the termination of Gtl2 transcription. We identified 36 upregulated and 57 downregulated miRNAs in Gtl2+/polyA, 68 upregulated and 148 downregulated miRNAs in Gtl2polyA/+ and 76 upregulated and 133 downregulated miRNAs in Gtl2polyA/polyA (Fig. 2e) (Tab. S4). As most of the downregulated miRNAs were located in the Dlk1-Dio3 domain in Gtl2polyA/+ and Gtl2polyA/polyA (Fig. 2f), we clustered the miRNAs in the domain and found that the four genotypes could be clearly divided into two groups: the survival group including Gtl2+/+ and Gtl2+/polyA, and the death group including Gtl2polyA/+ and Gtl2polyA/polyA, among them, the expression of miRNAs in Gtl2+/polyA were almost unchanged, while the expression of miRNAs in Gtl2polyA/+ and Gtl2polyA/polyA decreased significantly (Fig. 2g). In Gtl2, Rtl1as, Rian and Mirg loci, one miRNA with high expression and one with low expression were selected for qRT-PCR verification, the results were identical to RNA-seq data, indicating that the termination of Gtl2 transcription led to the silencing of miRNAs in the domain (Fig. 2h). Taken together, these data indicated that our model successfully terminated the transcription of the maternally expressed RNAs in the Dlk1-Dio3 domain and the termination of embryonic development might be closely related to the silencing of these genes.

Silencing of maternally expressed RNAs in Dlk1-Dio3 domain results in hemorrhage in embryos

Because platelets undergo activities such as activation and aggregation when injured or bleeding [34], the platelet activation in Kyoto Encyclopedia of Genes and Genomes (KEGG) that appeared in Gtl2polyA/+ and Gtl2polyA/polyA, but not in Gtl2+/polyA, caught our attention (Fig. 3a and S2a) (Tab. S5). Gene Ontology (GO) enrichment analysis revealed that GO terms enriched in both Gtl2polyA/+ and Gtl2polyA/polyA, but not in Gtl2+/polyA, were mainly associated with positive regulation of apoptotic process, embryo development, angiogenesis and wound healing (Fig. 3b and S2b) (Tab. S6). These features are associated with the phenotype of bleeding.

Fig. 3

Termination of maternally expressed RNAs in Dlk1-Dio3 domain transcription leads to hemorrhage. a KEGG analysis of target genes of DElncRNAs in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5. b GO analysis of target genes of DElncRNAs in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5. c Venn diagram shows the intersection of DEmiRNAs in Gtl2polyA/+, DEmiRNAs in Gtl2polyA/polyA, and all miRNAs in Dlk1-Dio3 domain. d KEGG analysis of target genes of DEmiRNAs in Dlk1-Dio3 domain in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5. e Blood smear of peripheral blood at E12.5, yellow arrows point to avtivatee and aggregated platelets, Scale bar, 40 μm. f GO analysis of target genes of DEmiRNAs in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5

We hypothesized that miRNA changes in the death group may play more roles. Therefore, we found the common DEmiRNAs of Gtl2polyA/+ and Gtl2polyA/polyA, and overlapped with the miRNAs in the Dlk1-Dio3 domain to obtain 95 DEmiRNAs (Fig. 3c). We further analyzed the KEGG and GO of target genes of these 95 DEmiRNAs in Gtl2polyA/+ and Gtl2polyA/polyA (Tab. S4). KEGG enrichment analysis revealed platelet aggregation in both Gtl2polyA/+ and Gtl2polyA/polyA (Fig. 3d and S2c) (Tab. S7). In parallel, activated and aggregated platelets were also identified from cord blood smears of Gtl2polyA/+ and Gtl2polyA/polyA embryos (Fig. 3e). This suggests that bleeding or even injury did occur in Gtl2polyA/+ and Gtl2polyA/polyA embryos. In addition, VEGF signaling pathway also appeared in Gtl2polyA/+ (Fig. 3e), and combined with the related terms such as angiogenesis that appeared in GO of DElncRNAs (Fig. 3b), we speculated that the termination of maternally expressed RNA transcription may cause some effects on blood vessels. GO enrichment analysis revealed liver development in both Gtl2polyA/+ and Gtl2polyA/polyA, but not in Gtl2+/polyA (Fig. 3f and S2d) (Tab. S8), and the same term was also appeared in the GO enrichment of lncRNAs (Fig. 3b). Taken together, the DEmiRNAs suggested that the organ affected by ncRNA silencing was the liver.

Silencing of maternally expressed RNAs in Dlk1-Dio3 domain results in apoptosis

Finally, we identified 249 upregulated and 515 downregulated mRNAs in Gtl2+/polyA, 609 upregulated and 3147 downregulated mRNAs in Gtl2polyA/+ and 220 upregulated and 214 downregulated mRNAs in Gtl2polyA/polyA (Fig. 4a) (Tab. S9). We found the Dlk1-Dio3 region in the mRNA profile (Fig. 4b), RNA-seq showed that the expression levels of Dlk1 and Rtl1 in Gtl2polyA/+and Gtl2polyA/polyA were 2 times higher than those in Gtl2+/+, but there was no significant difference between Gtl2+/polyA and Gtl2+/+. The expression levels of Dio3 in Gtl2polyA/polyA were 2 times higher than those in Gtl2+/+, while the expression level in Gtl2+/polyA and Gtl2polyA/+ was slightly lower than that in Gtl2+/+, which may account for the lower expression of Dio3. (Fig. 4c). The results of qRT-PCR were consistent with the trend of RNA-seq (Fig. 4d), indicating that the termination of transcription of the maternally expressed RNAs in Dlk1-Dio3 domain activates the expression of the paternally expressed genes.

Fig. 4

The genetic changes within the Dlk1-Dio3 domain leading to apoptosis. a Volcano plot of DEmRNAs in Gtl2polyA/+ and Gtl2polyA/polyA embryos at E12.5. b mRNA-seq of embryos at E12.5. c The expression of Dlk1, Rtl1 and Dio3 in RNA-seq data. d qRT-PCR analysis of Dlk1, Rtl1, and Dio3, β-actin was used as an internal control, n = 6. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01; ***, p < 0.001. e KEGG analysis of DEmRNAs in Gtl2polyA/+ and Gtl2polyA/polyA. f GO analysis of target genes of DEmRNAs in Gtl2polyA/+ and Gtl2polyA/polyA. g Relative activity of Caspase 9, n = 3. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05

KEGG enrichment analysis revealed that apoptosis was enriched in DEmRNAs of both Gtl2polyA/+ and Gtl2polyA/polyA, but not of Gtl2+/polyA (Fig. 4e and S3a) (Tab. S10). GO terms enriched in DEmRNAs of both Gtl2polyA/+ and Gtl2polyA/polyA were mainly associated with apoptosis, wound healing, VEGF signaling pathway and platelet aggregation (Fig. 4f and S3b) (Tab. S11). Among them, platelet aggregation was the same as the experimental results of ncRNA enrichment, wound healing and VEGF signaling pathway appeared again, which had to attract our attention and needed further verification.

In addition, apoptosis occurs simultaneously with the above pathways, we speculated that the termination of maternally expressed RNA transcription may cause some effects on blood vessels and hemorrhage, and the relative activity of Caspase 9 in Gtl2polyA/+ and Gtl2polyA/polyA was significantly higher than that in Gtl2+/+, but there was no significant change in Gtl2+/polyA (Fig. 4g). Taken together, these data suggest that silencing of maternally expressed RNAs in Dlk1-Dio3 domain cause apoptosis, which may be responsible for vascular damage and hemorrhage.

Silencing of maternally expressed RNAs in Dlk1-Dio3 domain results in vascular damaged and hemorrhage in fetal livers

The GO functional enrichment of the maternally expressed RNA in the interval was concentrated in the liver (Fig. 3b and f), so we speculated that the above phenomenon may be the most obvious in the liver. To determine the effect of silencing the maternally expressed RNAs in the Dlk1-Dio3 domain on fetal liver, we obtained fetal livers from all four-genotypes at E12.5. It was found that there was no significant change in appearance of the three mutants compared with Gtl2+/+ (Fig. 5a). Immunofluorescence of Ki67 at E14.5 fetal livers showing that significant cell proliferation occurred in Gtl2polyA/+ and Gtl2polyA/polyA, but not in Gtl2+/polyA fetal livers (Fig. S4). The above results have shown that the increase in body and fetal liver weight may be due to cell proliferation which was consistent with the results of KEGG and GO enrichment by the ncRNAs. The trend of weight in four-genotype fetal livers was consistent with that of embryos (Fig. 5b). To further observe the changes at the cellular level in the fetal livers, histological analysis showed that there was no significant change between Gtl2+/polyA and Gtl2+/+ (Fig. 5C). The area of blood vessels in the fetal livers of Gtl2polyA/+ and Gtl2polyA/polyA increased, and some blood vessels were damaged (Fig. 5c). At the same time, a variety of cells, such as hepatocytes and endothelial cells, showed suspected apoptosis phenomena such as cell volume reduction, membrane rupture, karyopyknosis and cell fragmentation (Fig. 5c). Single cell suspension of intact liver was prepared for caspase 9 activity detection. The results showed that the relative activities of Gtl2polyA/+ and Gtl2polyA/polyA were more than 2 times higher than that of Gtl2+/+ (Fig. 5d). TUNEL enabled us to locate apoptosis in fetal livers, which also confirmed the previous speculation that cells in fetal liver did indeed undergo apoptosis (Fig. 5e). Relative fluorescence intensity of TUNEL showed that the apoptosis degree of Gtl2polyA/+ and Gtl2polyA/polyA were significantly higher than that of Gtl2+/+, which was consistent with the previous results (Fig. 5f). These results suggest that termination transcription of maternally expressed RNA in Dlk1-Dio3 domain leads to apoptosis in fetal livers.

Fig. 5

Termination of maternally expressed RNAs in Dlk1-Dio3 domain transcription leads to apoptosis and vascular damage in fetal liver. a Representative pictures of fetal liver at E12.5. Scale bar, 1 mm. b Bar chat showing weight of fetal livers at E12.5, each dot represents a fetal liver. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05. c Representative images of HE staining of fetal livers of each genotype at E12.5. Scale bar, 40 μm. d Relative activity of Caspase 9 in fetal liver, n =  3. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *p < 0.05. e TUNEL fluorescence staining of apoptotic cells in fetal livers, Scale bar, 40 μm. f Quantitative analysis of the immunohistochemical staining for TUNEL, n = 5. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01. g Representative images of immunohistochemistry staining for CD31, apoptosis of endothelial cells are shown as yellow arrows. Scale bar, 40 μm. h Quantitative analysis of the immunohistochemical staining for CD31, n = 5. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01. i Number of blood vessels per square millimeter in each genotype fetal livers, n = 5. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01; ***, p < 0.001

Next, combined with the previous histological observation of vascular rupture, we further analyzed the blood vessels of the fetal liver. Immunohistochemistry of CD31 was localized to liver endothelial cells and the results showed fragments of endothelial cell apoptosis, which damaged blood vessels and caused erythrocytes to flow out of the damaged blood vessels, resulting in liver hemorrhage, which was consistent with previous histological analysis (Fig. 5g). Quantitative analysis of vessel area showed that Gtl2+/polyA, Gtl2polyA/+ and Gtl2polyA/polyA were 0.98-fold, 1.26-fold and 1.35-fold larger than that of Gtl2+/+ (Fig. 5h). The number of blood vessels per square millimeter showed that Gtl2+/polyA, Gtl2polyA/+ and Gtl2polyA/polyA were 1.03-fold, 0.58-fold and 0.56-fold of Gtl2+/+, respectively (Fig. 5i). That is, the average vessel area of Gtl2polyA/+ and Gtl2polyA/polyA was significantly larger than that of Gtl2+/+, but the number of vessels was significantly less than that of Gtl2+/+, and the cause of vascular injury and hepatic hemorrhage was apoptosis of endothelial cells.

Silencing of maternally expressed RNAs in Dlk1-Dio3 domain had no significant effect on other organs

In order to investigate whether silencing of maternally expressed RNAs in Dlk1-Dio3 domain has an effect on tissues other than liver, we conducted a further investigation in embryos. Firstly, histological morphology showed that no significant difference, except liver, was observed in the mutant embryo as a whole (Fig. 6a). Secondly, as previous studies have confirmed the effect of maternally expressed RNA silencing in this domain on brain, lung and skeletal muscle [10, 11, 35], we mainly observed these organs (Fig. 6b). We divided the brain into three parts, telencephalon, diencephalon and metencephalon and the results showed that the brain of E12.5 was not fully developed, and there was no significant difference in telencephalon, diencephalon and metencephalon in the mutant (Fig. 6b). The tracheas of all mutants and Gtl2+/+ in lung were well developed, and there were some small blood vessels, and the alveoli were not yet formed at this stage (Fig. 6b). There were no significant differences in the bones of the legs, nor any defects in the skeletal muscles (Fig. 6b). It is worth mentioning that our currently published study reported that Gtl2 transcription termination caused developmental abnormalities in the heart caused by slow proliferation of epicardial cells (Fig. 6b) [36]. These results suggest that silencing of maternally expressed RNAs transcription in Dlk1-Dio3 at E12.5 does not have a significant effect on major organs except the liver, and the exploration of the effect on the heart needs to be further improved.

Fig. 6

Termination of maternally expressed RNAs in Dlk1-Dio3 domain transcription had no significant effect on other embryonic organs. a Representative images of HE staining of embryos at E12.5. Scale bar, 1 mm. b Representative images of HE staining of major organs. Scale bar, 80 μm. c qRT-PCR analysis of lncRNAs and mRNAs in Dlk1-Dio3 domain of major organs, β-actin was used as an internal control, n = 5. Error bars, ± SEM, Student’s t-test, n.s., p > 0.05; *, p < 0.05; **, p < 0.01; ***, p < 0.001. d In situ hybridization of Dlk1 and Gtl2. Scale bar, 1 mm

qRT-PCR detected the expression of genes in each organ within the range, and the expression trend of maternally expressed RNA in each organ was basically the same, all of them were significantly down-regulated to almost silence, and the expression of paternally expressed genes, Dlk1 and Rtl1, in brain, heart, liver, lung and leg increased about twice, which is basically the same as that of embryo (Fig. 6c). The changes of Dio3 in heart and lung were not very obvious, which may be the reason for the low expression of Dio3 within the two organs (Fig. 6c). We selected Dlk1 and Gtl2 as high expression genes from maternally expressed RNA and paternally expressed RNA, respectively, and in situ hybridization showed the same expression as qRT-PCR, and both Dlk1 and Gtl2 were highly expressed in liver, but not in other organs (Fig. 6d and S5a, b). Thus, the above results indicate that the changes of genes in major organs in the interval are basically the same as those in the whole embryo, and the effect on the liver is greater than that on the other organs, which needs further observation.

Silencing of Gtl2, but not Gtl2-DMR and IG-DMR, leads to Dlk1 activation

To explore whether the cause of the phenotype was caused by changes in the gene itself or by the methylation status of the Gtl2-DMR and IG-DMR, we found a SNP(A/G) site on Dlk1 that can distinguish between C57BL/6N(G) and DBA/2J(A). After obtaining the cDNA of the embryo of the female mouse of DBA/2J and the male mouse of C57BL/6N after mating, the cDNA fragment containing the SNP site was sequenced. The expression in Dlk1 of Gtl2+/+ and Gtl2+/polyA was mono-allelic, while the expression of Gtl2polyA/+ and Gtl2polyA/polyA was bi-allelic. At the same time, we also tested the offspring obtained by mating male DBA/2J mice with female C57BL/6N mice, and the same conclusion was obtained (Fig. 7a). This means that the silencing of maternally expressed RNAs in Dlk1-Dio3 cause the imprinting loss of Dlk1, although, paternally expressed genes were not edited in the domain, Dlk1, Rtl1 and Dio3 were still affected, which was also likely to have some impact on embryonic development.

Fig. 7

Loss of imprinting occurs in Dlk1, which was not caused by Gtl2-DMR and IG-DMR. a Biallelic expression of Dlk1 in embryos at E12.5, red arrow points to SNPs. b, c Methylation analysis of Gtl2-DMR (b) and IG-DMR(c). Methylated CpG dinucleotides are represented by filled circles, while unmethylated CpG dinucleotides are represented by open circles

Subsequently, we examined Gtl2-DMR and IG-DMR. The lengths of IG-DMR and Gtl2-DMR detected in this study were 497 bp and 385 bp, correspondingly. These sequences contained 19 CpG and 2 SNP sites, and 15 CpG and 1 SNP sites in C57BL/6N and DBA/2J hybrids, respectively [11, 34] (Fig. 7b and c). The results shown that both of two DMRs maintained the parental origin specific methylation status in E12.5 embryos (Fig. 7b and c). Taken together, the above data indicated that the variation of paternally expressed genes in the domain were caused by silencing of Gtl2 rather than by changes in methylation status of Gtl2-DMR and IG-DMR.