Maternal DNMT3A-dependent de novo methylation of the zygotic paternal genome inhibits gene expression in the early embryo

De novo DNA methylation (DNAme) during mammalian spermatogenesis yields a densely methylated genome, with the exception of CpG islands (CGIs), which are hypomethylated in sperm. Following fertilization, the paternal genome undergoes widespread DNAme loss before the first S-phase. Paradoxically, recent mass spectrometry analysis revealed that a low level of de novo DNAme occurs exclusively on the zygotic paternal genome. However, the loci involved and impact on genic transcription was not addressed. Here, we employ allele-specific analysis of wholegenome bisulphite sequencing (WGBS) data and show that a number of genomic regions, including several dozen CGI promoters, are de novo methylated on the paternal genome in 2-cell embryos. A subset of these promoters maintains DNAme through development to the blastocyst stage. Consistent with zygotic paternal DNAme acquisition (PDA), many of these loci are hypermethylated in androgenetic blastocysts but hypomethylated in parthenogenetic blastocysts. Strikingly, PDA is lost following maternal deletion of Dnmt3a. Furthermore, a subset of promoters showing PDA which are normally transcribed from the paternal allele in blastocysts show premature transcription at the 4-cell stage in maternal Dnmt3a knockout embryos. These observations uncover an unexpected role for maternal DNMT3A activity in postfertilization epigenetic reprogramming and transcriptional silencing of the paternal genome.


ABSTRACT
De novo DNA methylation (DNAme) during mammalian spermatogenesis yields a densely methylated genome, with the exception of CpG islands (CGIs), which are hypomethylated in sperm. Following fertilization, the paternal genome undergoes widespread DNAme loss before the first S-phase. Paradoxically, recent mass spectrometry analysis revealed that a low level of de novo DNAme occurs exclusively on the zygotic paternal genome. However, the loci involved and impact on genic transcription was not addressed. Here, we employ allele-specific analysis of whole-

INTRODUCTION
Male germ cell development in mammals is characterized by widespread de novo DNA methylation (DNAme) and compaction of DNA via histone-to-protamine exchange 1

. As
DNAme is largely maintained throughout spermatogenesis, the genomes of mature spermatozoa harbor characteristically high levels of DNAme 2 . Following fertilization, the paternal genome undergoes another profound change in chromatin state, including replacement of protamines with histones and a global reduction in DNAme before the first S-phase 3,4 . Subsequently, DNAme levels on both parental genomes are progressively reduced with each DNA replication cycle 5 . While passive demethylation in the early embryo is likely explained by sequestration of the maintenance DNA methyltransferase DNMT1 in the cytoplasm [6][7][8] , the mechanism of active demethylation remains controversial. Indeed, while TET3-mediated oxidation followed by base excision repair has been implicated in this process, a TET-independent mechanism is also clearly involved [9][10][11][12][13][14][15] . Regardless, whole-genome bisulphite sequencing (WGBS) analyses reveal that DNAme levels on both parental genomes reach a low point in inner cell mass (ICM) cells of embryonic day 3.5 (E3.5) mouse blastocysts, followed by widespread de novo DNAme during post-implantation development [16][17][18] . This wave of genome-wide demethylation followed by remethylation is conserved in human embryonic development, albeit with slower kinetics [19][20][21] . Notably, disruption of the machinery required for the establishment or maintenance of DNAme result in infertility and/or embryonic lethality in mice, revealing the importance of DNAme homeostasis in early mammalian development 7,8,13,[22][23][24][25][26][27] .
Although the extent of histone post-translational modification (PTM) maintenance following fertilization is controversial 34,36,37 , the retention and/or rapid deposition of H3K4 methylation at CGI promoters may protect these regions against de novo DNAme 38 in the developing male germline as well as the early embryo. Indeed, most CGI promoters are enriched for H3K4me3 and remain hypomethylated on both parental genomes throughout early embryonic development and in adult tissues 16,39 . H3K4me2/3 also likely facilitate the initiation of transcription from the paternal genome during zygotic gene activation (ZGA) 28,34 , which marks the transition between oocyte and embryonic transcriptional programmes 20 . As in sperm, CGIs in oocytes are generally hypomethylated and harbor nucleosomes enriched for H3K4me3 36 and/or H3K27me3 40,41 . However, a subset of CGIs, are de novo methylated in growing oocytes by DNMT3A, which is highly expressed in oocytes [42][43][44] . Paradoxically, despite widespread DNAme loss from the paternal genome in mouse zygotes, maternal DNMT3A is also clearly detected in the paternal pronucleus at this stage 14,27 and ongoing de novo DNAme is required for maintaining paternal allele DNAme of the paternally imprinted H19 locus in early embryos 45 . Further evidence that the zygotic paternal genome is subject to de novo DNAme was demonstrated in a recent study employing immunofluorescence (IF) and ultrasensitive liquid chromatography/mass spectrometry 14 . However, the genomic regions subject to such paternal DNAme acquisition (PDA) in the early embryo and the relevance of this phenomenon to transcriptional regulation of the paternal genome has not been systematically addressed.
To determine which loci gain DNAme immediately following fertilization, we carried out an allele-specific analysis of WGBS data from 2-cell (2C) F1 hybrid embryos 16 and identified specific genomic regions, including CGI promoters, that show PDA. Corroborating these enigmatic findings, we observe PDA of an overlapping set of CGI promoters in androgenetic but not in parthenogenetic blastocysts. Allele-specific analysis of ChIP-seq data from 2C embryos reveals that PDA is accompanied by loss of H3K4me3 over the same regions, indicating that such DNAme may inhibit transcription from the paternal allele in early embryonic development. Indeed, we show that PDA is lost in the absence of maternal DNMT3A and a subset of hypomethylated genes are concomitantly upregulated specifically from the paternal allele in the 4C embryo. Taken together, these experiments reveal that beyond its role in maternal imprinting, DNMT3A methylates a subset of genes on the paternal genome in the zygote, inhibiting their expression in preimplantation development.

The paternal genome undergoes de novo DNAme following fertilization
To trace parent-specific DNAme levels following fertilization and throughout embryonic development with single-nucleotide resolution, we first processed publicly available WGBS data derived from primordial germ cells (PGCs), spermatozoa and MII oocytes 2,31,46 . We then applied our recently developed allele-specific pipeline MEA 47 to WGBS data generated from 2C (55X coverage) as well as 4C, ICM, E6.5 and E7.5 F1 hybrid embryos 16 . This integrated analysis yielded female/maternal and male/paternal DNAme profiles through fertilization and early development ( Fig. S1a-b). Consistent with previous IF data 3,4 , comparison of DNAme levels in mature gametes and 2C embryos reveals an overall decrease on the maternal and paternal genomes of 8% and 43%, respectively (Fig. S1c). Surprisingly however, coincident with global DNAme loss across the paternal genome, robust DNAme gain (defined as an increase of ≥30%) was detected at ~4% of all hypomethylated regions (≤20%) in sperm, totalling 1.4 Mbp of the mappable genome (Fig. 1a). De novo DNAme of the paternal genome is consistent with the rapid translocation of maternal DNMT3A into the zygotic paternal pronucleus ( Fig. 1c). Furthermore, a subset of TSSs showing PDA maintain such DNAme on the paternal genome to the blastocyst stage (Fig. 1d). While PDA is not restricted to CpG-rich regions, we focused our analyses on CGI promoters, as DNAme is reported to have the strongest impact on transcription of this class of promoters 48,49 . To generate a curated list of CpG-rich hypomethylated TSSs in sperm, we categorized promoters by CpG density (high, intermediate and low), as described for the human genome 49 . As expected, DNAme levels and CpG density are anti-correlated ( Fig. S2a). To minimize the potential confounding effects of strain-specific differences when comparing DNAme levels of parental genomes, we next determined the variation in DNAme levels in C57BL/6J versus DBA/2J sperm using published WGBS datasets 2,16 . Methylation profiles from these strains show a strong correlation (r 2 =0.98,  Table 1). Notably, lower levels of H3K4me3 in sperm are unlikely to explain their propensity to gain DNAme, as the distribution of H3K4me3 levels in sperm was not significantly different between these CGIs and those that showed no gain in DNAme (p=0.46, Fig. 1g). However, relative to CGIs that remain unmethylated, PDA CGI promoters show a significantly greater decrease of H3K4me3 levels over the paternal allele in 2C embryos (p=2.25E-5). Though we cannot discriminate between active demethylation of H3K4 and histone H3 turnover, these results reveal that de novo DNAme of CGIs on the paternal genome is generally accompanied by a reduction of H3K4me3. For example, Tuba3a, Gdap2 and Bpi gain DNAme across 28, 58 and 7 CpGs proximal to their TSSs (+/-300bp), respectively, coincident with loss of H3K4me3 on the paternal allele (Fig. 1h), while the control gene Ufc1 displays persistent DNA hypomethylation and H3K4me3 enrichment on the paternal genome both pre-and post-fertilization (Fig. 1i). Surprisingly, only 11 of the PDA genes are methylated (≥20% DNAme) in MII oocytes and none gain DNAme on the maternal genome in 2C embryos (Fig. S1e). In contrast, maternally imprinted genes, such as Impact and Snurf, show maternal allele-specific DNA hypermethylation in the ICM and paternal allele-specific enrichment of H3K4me3 throughout early embryonic development, as expected (Fig. S2c). Taken together, these analyses reveal that a specific subset of CGI promoters are de novo DNA methylated exclusively on the paternal genome following fertilization, concomitant with reduced H3K4me3.

DNAme at many PDA loci is maintained through the blastocyst stage
To determine whether DNAme at PDA sites persists through the wave of global DNAme erasure, we scored paternal DNAme levels at these loci using WGBS data from F1 hybrid ICM cells 16 . Of 63 PDA loci, 53 had sufficient allele-specific WGBS coverage to ascertain paternal DNAme levels in ICM cells. Relative to CGI promoters that remain hypomethylated following fertilization, those that show PDA retain higher DNAme on the paternal allele in ICM cells (p=4.97E-6, Fig. 2a), albeit at lower levels than observed in 2C embryos. In contrast, 19 CGI promoters that show PDA, including Gdap2 (Fig. 2b), are hypomethylated (mean <5%) in blastocysts (Supplemental Table 1   To confirm the persistence of paternal allele-specific DNAme and expand upon the number of loci at which PDA is likely occurring in the early embryo, we isolated genomic DNA from isogenic (C57BL/6NJcl) androgenetic blastocysts, which contain only paternally inherited genomes, and conducted WGBS. Compared to CGI promoters that remain hypomethylated in 2C embryos, those that show PDA exhibited a significantly greater level of DNAme in such bipaternal blastocysts (p=9.30E-8, Fig. 2a).
Furthermore, the majority of CGI promoters that show persistence of DNAme (≥20%) on the paternal allele in ICM cells, including 8 PDA genes, show ≥10% DNAme in androgenetic blastocysts (Fig. S2d). In contrast, analysis of our previously published WGBS data from parthenogenetic blastocysts 43 , in which both genomes are maternally derived, reveals that DNAme remains low at all of the CGI promoters showing zygotic PDA (Fig. 2a), with the exception of 6 of the 11 that are already hypermethylated in MII oocytes (Supplemental Table 1). Thus, loci showing PDA are de novo methylated in the early embryo exclusively when at least one paternal genome is present.
Given that androgenetic and parthenogenetic blastocysts are uniparental, we were able to extend our parental genome-specific DNAme analysis from 18 to 90% of all annotated autosomal CGI promoters. As expected, maternally imprinted CGI promoters, such as Peg10 and Zdbf2, which are hypermethylated exclusively on the maternal allele in normal embryos, remain hypermethylated in parthenogenetic blastocysts (Fig. 2c). However, only one CGI promoter (Ccdc114) shows a ≥20% DNAme gain in these cells relative to germinal vesicle (GV) oocytes. In contrast, 28 CGI promoters show such a gain in androgenetic blastocysts relative to sperm, confirming that the paternal genome is the preferred target for such post-fertilization de novo DNAme ( Fig. 2d and Supplemental Table 1). While 17 of these loci (including Thy1) do not harbor a genetic variant and could therefore not be assessed in the F1 hybrid dataset, six, including Tuba3a, Bpi, Them7, Shisa7, Syn3 and A230077H06Rik, were also identified as PDA genes in our allele-specific analysis of F1 hybrid embryos (Fig.   2e) and four (H1fnt, Dbx2, Tbx4 and Prss39) showed a gain in DNAme of 9-25%, just below our original threshold of >30%. These results indicate that PDA occurs independent of DBA/2J-specific variants and is thus a bona fide parent-of-origin effect at these loci.

Relationship between histone PTMs and PDA
To determine whether promoter regions showing PDA share a common DNA motif that may render them susceptible to de novo DNAme, we performed motif discovery using HOMER 50 and MEME 51 . However, we did not detect any common DNA motif around the TSSs (+/-300 bp) of genes that show PDA. Therefore, we focussed on the relationship between histone PTMs and paternal DNAme acquisition and/or subsequent DNAme maintenance at these sites. We analyzed published H3K4me2 34 , H3K4me3 30,36 , H3K9me3 17 , H3K27me3 30,40 and H3K36me3 37 ChIP-seq data from sperm, oocytes, zygote (1C), 2C and ICM cells using ChromHMM 52 . No single histone PTM or combination thereof in sperm predicted PDA at CGI promoters versus those that remain hypomethylated ( Fig. S3a-b). However, a striking enrichment of H3K9me3 was observed in the zygote and 2C stage at PDA sites (Fig. S3b) and those that show persistence of DNAme to the ICM show an even greater enrichment of H3K9me3 in the zygote and 2C embryo (Fig. S3c). Indeed, the paternal allele of many CGI promoters show a positive association between the levels of H3K9me3 at the 2C stage and DNAme in the ICM (Fig. S3d), including that of the PDA gene Tuba3a (Fig. S3d-e).  Surprisingly, only 17 of these hypomethylated genes, including Tuba3a and Gdap2, are transcribed at this stage ( Fig. 3b-c). To determine whether the inactive genes harbor histone marks associated with transcriptional repression, we integrated previously published H3K27me3 and H3K9me3 ChIP-seq datasets derived from oocytes and F1 hybrid embryos 17,40 . Consistent with a role for H3K9me3 in promoting DNAme maintenance, all 6 hypermethylated PDA loci are enriched for this histone PTM in oocytes ( Fig. S4a-b). In contrast, 23 of the 40 silent hypomethylated CGI promoters are embedded within H3K27me3-enriched domains (TSS +/-10kb, RPKM ≥0.5), 17 of which are maintained on the maternal allele at least to the late 2C stage (Fig. 3d). The TSS of Bpi for example is embedded within an extended H3K27me3 domain in oocytes which persists to the 2C stage (Fig. 3e). Notably, H3K27me3 was recently implicated in transcriptional silencing of genes subject to non-canonical maternal imprinting in the mouse 41,57 . These results reveal that at the 2C stage, the paternal and maternal alleles of a subset of PDA loci are distinctly marked by DNAme and H3K27me3, respectively.

Maternal DNMT3A is required for PDA
As DNMT3A is required for de novo DNAme in the female germline 31 and maternal DNMT3A persists in the zygote (Fig. S1d), we next wished to determine whether maternal DNMT3A is responsible for PDA. We crossed oocyte-specific homozygous Dnmt3a knock-out (matKO) C57BL/6J females with wild-type (WT) DBA/2J males and performed WGBS on F1 hybrids at the early-mid 2C stage (Fig. 4a, 5X mean sequencing coverage). As expected, global maternal DNAme levels in 2C matKO embryos mirror those in Dnmt3a matKO GV oocytes (Fig. S5a). Furthermore, while the paternal gametic differentially methylated regions (gDMRs) H19, Dlk1-Gtl2 and Rasgrf1 show no decrease, DNAme at maternal gDMRs is dramatically reduced in matKO relative to control embryos (Fig. 4b). Importantly, a positive correlation was observed when comparing paternal DNAme levels of CGI promoters with sufficient allele-specific coverage in our control dataset with published 2C WGBS data 16 , with the majority of PDA genes consistently showing relatively high DNAme on the paternal allele (Fig.   S5b). In the absence of maternal DNMT3A however, DNAme is lost on the paternal allele of all 23 CGI promoters showing PDA for which allelic methylation can be deduced (Fig. 4c). At the Tuba3a, Gdap2 and Snx1 promoters for example, PDA is lost in matKO embryos across all CpGs in the promoter region (Fig. 4d). While the remaining 40 PDA loci lack allele-specific coverage in our matKO samples, analysis of total DNAme levels reveals that 33 are clearly hypomethylated (<4%) (Supplemental Table 1), as shown for Bpi (Fig. 4d). Taken together, these results demonstrate that DNAme acquisition on the paternal allele immediately following fertilization is dependent upon maternal DNMT3A.

Loss of PDA results in ectopic expression from the paternal allele
Hypermethylation of CGI promoters is associated with transcriptional silencing 48,49 . To determine whether DNAme of genes showing PDA impacts their expression specifically from the paternal allele, we conducted strand-specific RNA-seq on early-mid 2C as well as 4C and blastocyst-stage F1 matKO embryos (Fig. 5a and Fig. S6a). A strong correlation was observed between matKO and WT 2C, 4C and blastocyst-stage embryos as well as previously published WT transcriptomes of the same developmental stages (Fig. S6b), indicating that maternal depletion of DNMT3A does not disrupt global transcriptional programming. This is consistent with previous observations that Dnmt3a matKO embryos develop normally until E9.5 26,58 . We next determined whether genes known to be regulated by maternally established genomic imprints show loss of imprinting in Dnmt3a matKO embryos. Consistent with their loss of DNAme (Fig. 4b), Mcts2, Plagl1, Snurf, Peg10 & Zdbf2 imprinted genes were significantly upregulated (≥2-fold change, Wald Chi-squared test, Benjamini-Hochberg adjusted p-value ≤0.1) from the maternal allele (Fig. S6c-e). Thus, DNMT3A expression in the oocyte is required for maternal allele-specific transcriptional silencing of a subset of genes in the preimplantation embryo.  A similar analysis of genes upregulated from the paternal allele in 2C matKO embryos yielded only Gdap2, a PDA gene, as significantly upregulated (Fig. 5b-c).  Table 1). Importantly, transcription of all upregulated PDA genes initiates from the aberrantly hypomethylated CGI promoter, excluding alternative promoter usage as an explanation for the increased expression from the paternal allele (Fig. 5c-d). Furthermore, no significant upregulation of PDA genes was observed from the maternal allele, consistent with DNAme-independent silencing at this stage (Fig. 3d), or when total (allele-agnostic) transcript levels were analyzed (Fig. S7ac). Thus, upregulation of PDA genes occurs exclusively from the paternal allele, raising the question, can expression data be exploited to identify additional PDA genes?
Due to the density of naturally occurring polymorphisms and depth of WGBS sequencing coverage, paternal DNAme dynamics could be measured in our earlier F1 hybrid analysis over only 4,434 of the 12,253 filtered autosomal CGI promoters (Fig.   1e). Of the additional 7,724 autosomal genes with CGI promoters that include at least 1 exonic genetic variant between parental strains and are expressed from the paternal allele in our RNA-seq data, 5 were significantly upregulated from the paternal genome in matKO 4C embryos, including Npm3, Eif3b, 2610524H06Rik, Crk & Pim1 (Fig. 5b-c).
Importantly, as for bona fide PDA loci, none of these genes show a change in expression from the maternal allele (Fig. S7b,d), nor are they upregulated from the paternal allele in matKO ICM (Fig. 5b-c). Furthermore, analysis of total DNAme levels clearly shows that these loci are methylated in WT but not in matKO 2C embryos, indicating that these CGI promoters are direct targets of maternal DNMT3A (Fig. 5e).
While the lack of genetic variants precludes the measurement of allele-specific DNAme over these regions, analysis of allele-agnostic DNAme levels in WT 2C embryos revealed that 4 of these 5 CGI promoters show a ≥18% gain in DNAme relative to their methylation state in sperm and oocyte (Supplemental Table 1). Thus, in addition to the PDA genes described above, these 4 loci are likely de novo methylated on the paternal allele in 2C embryos, and repressed by this mark in 4C embryos. Taken together, these results reveal that the absence of maternal DNMT3A in the early embryo leads to failure of de novo DNAme on the paternal allele of a subset of genes and, in turn, their ectopic expression.

DISCUSSION
Using allele-specific analyses of early embryos, we determined that 4% of all hypomethylated regions in sperm, including at least 63 CGI promoters, are de novo methylated on the paternal genome by the 2C stage. An Independent analysis of androgenetic blastocysts revealed 86 CGI promoters that show a ≥10% gain in DNAme relative to sperm, 15 of which overlap with the PDA genes identified in normal 2C embryos (Supplemental Table 1). These observations are particularly surprising, given that the hypermethylated sperm genome is globally demethylated in the mouse zygote. Isolation of androgenetic blastocyst. Diploid androgenones were prepared as previously described 75 . Oocyte and spermatozoa were isolated from B6D2F1/Jcl and C57BL/6NJcl mice (Clea Japan, Tokyo, Japan), respectively. Briefly, enucleated oocytes were in vitro fertilized and zygotes with two male pronuclei were cultured for 4 days in KSOM medium at 37 o C and 5% CO2 76 . 5 blastocyst pools were collected in triplicate for WGBS library construction.
Superovulation was induced using PMSG/hCG and MII oocytes were collected from oviducts. Dnmt3a 2lox ;Zp3-cre MII oocytes were artificially inseminated with DBA/2J spermatozoa. Cumulus cells were removed using hyaluronidase after insemination and embryos were cultured in KSOM at 37 o C and 5% CO2. Early-mid 2C embryos were DNMT3A immunofluorescence. Dnmt3a 2lox ;Zp3-cre females were mated with JF1 males and IF was performed on one-cell zygotes. DNMT3A was detected using the IMG-268 IMGENEX antibody, and DNA was counterstained using propidium iodide, as described previously 27 .
WGBS and RNA-seq library construction and sequencing. Lysates of androgenetic blastocysts were spiked with 0.1 ng lambda phage DNA and subjected to WGBS library construction according to the PBAT protocol for single-read sequencing 46  NGS data processing. Reads were trimmed using Trimmomatic v0.32 80 and processed for total and allele-specific alignments using MEA v1.0 47 using default parameters and the mm10 reference genome. 4 bases were removed from the 5' end of PBAT sequences. All publicly available NGS data (Supplemental Table 2) was reprocessed as above, with the exception of WGBS datasets from sperm 2 , oocytes 31 , parthenones 43 and primordial germ cells 46 , which were previously processed and filtered using identical parameters as in this study.
WGBS data analysis. DNAme levels over individual CpGs with ≥5x coverage (including allele-specific) were scored, with the exception of allele-specific alignments of WGBS datasets generated in this study, for which a ≥1x allele-specific read coverage cutoff was used to score methylation status. DNAme levels were calculated over CGI promoters using Bedops v2.4.27 81 and visualized using VisRseq v0.9.12 82 . Only CGI promoters that overlapped at least 2 informative CpGs separated by the maximum sequencing read length of the library were kept. Genome-wide 2-and 20kb bins were generated using Bedtools, and bins covered by at least 4 CpGs separated by over 1 read length in each dataset were used, and a random subset of 1,000 bins were visualized as parallel coordinate plots using VisRseq. Genome-wide 600bp bins with 100bp overlap were used to measure paternal allele DNAme dynamics in sperm and F1 hybrid 2C embryos. Bins covered by at least 2 CpGs separated by over 1 read length in each dataset were kept. Bins showing <20% DNAme in sperm and a ≥30% gain in 2C embryos were scored as showing PDA. Overlapping PDA regions were subsequently merged using Bedtools. NCBI RefSeq (default in VisRseq) and Ensembl Regulatory features (release 81) annotations were used to identify PDA regions that overlapped TSSs, gene bodies and enhancers.
ChIP-seq data analysis. Raw sequencing reads were reprocessed as described above into total and allele-specific genomic tracks. RPKM values were calculated over TSSs (+/-300bp) using VisRseq on the basis of normalized genomic tracks. ChromHMM v1.12 52 was employed to define distinct chromatin states (LearnModel, k=6) on the basis of filtered BAM files (BinarizeBam) using default parameters.
RNA-seq data processing. Raw sequencing reads were reprocessed as described above into total and allele-specific genomic tracks. Gene expression (RPKM) values over genic exons were calculated using VisRseq (NCBI Refseq). Paternal RPKM values for each gene showing PDA was then plotted on a parallel-coordinate plot using VisRseq. Oocyte CGI promoter expression (RPKM) was calculated over TSSs +/-300bp and normalized to total aligned reads. Correlograms were generated using Morpheus (https://software.broadinstitute.org/morpheus) on the basis of log2 transformed values.
For genome browser visualization, biological replicates were merged and total (alleleagnostic) and allele-specific genomic tracks were organized into UCSC Track Hubs as  Table 2 or the full list of data analyzed for this study.