DPPA2/4 and SUMO E3 ligase PIAS4 opposingly regulate zygotic transcriptional program

The molecular mechanism controlling the zygotic genome activation (ZGA) in mammals remains poorly understood. The 2-cell (2C)-like cells spontaneously emerging from cultures of mouse embryonic stem cells (ESCs) share some key transcriptional and epigenetic programs with 2C-stage embryos. By studying the transition of ESCs into 2C-like cells, we identified developmental pluripotency associated 2 and 4 (Dppa2/4) as important regulators controlling zygotic transcriptional program through directly up-regulating the expression of double homeobox (Dux). In addition, we found that DPPA2 protein is sumoylated and its activity is negatively regulated by small ubiquitin-like modifier (Sumo) E3 ligase protein inhibitor of activated STAT 4 (PIAS4). PIAS4 is down-regulated during ZGA process and during transitioning of ESCs into 2C-like cells. Depleting Pias4 or overexpressing Dppa2/4 is sufficient to activate 2C-like transcriptional program, whereas depleting Dppa2/4 or forced expression of Pias4 or Sumo2–Dppa2 inhibits 2C-like transcriptional program. Furthermore, ectopic expression of Pias4 or Sumo2–Dppa2 impairs early mouse embryo development. In summary, our study identifies key molecular rivals consisting of transcription factors and a Sumo2 E3 ligase that regulate zygotic transcriptional program upstream of Dux.


Introduction
Zygotic genome activation (ZGA) occurs predominantly at the 2-cell (2C) stage of mouse embryo and 4-to 8-cell stages in human embryo [1][2][3], which is essential for the development control passed from maternal to newly synthesized RNA and proteins. Any wrongdoing during ZGA may lead to termination of embryo development or have severe and long-term consequences for the life of an organism. However, the molecular regulation of ZGA in mammal is poorly understood. Recently, a rare subset of cells called 2C-like cells were found in mouse embryonic stem cell (ESC) cultures [4,5]. The 2C-like cells express high levels of ZGA transcripts including murine endogenous retrovirus (ERV)-L (MERVL) family of retroviruses and

Sumo2 and Sumo E3 ligase Pias4 inhibit zygotic transcriptional program in mouse ESCs
To study how zygotic transcriptional program is regulated, we first generated a 2C::tandem dimeric Tomato (tdTomato) reporter ESC line with transgenic expression of red fluorescence protein tdTomato under the control of MERVL 5 0 UTR [5]. As expected, around 0.5% 2C:: tdTomato-positive cells were typically present in ESC cultures (S1A Fig). Using this cell line, we confirmed that knocking down Sumo2, but not Sumo1, significantly increased the percentage of tdTomato-positive cells (S1A Fig). Importantly, knocking down Sumo2 also up-regulated the expression of MERVL and other classic 2C-specific genes including Dux, Cml2, zinc finger protein 352 (Zfp352), and Zscan4d (S1B Fig).
To figure out which Sumo E3 ligase is responsible for the repression of 2C-like state and zygotic transcriptional program, we knocked down nine known Sumo E3 ligases [31] individually with small interfering RNAs (siRNAs) in 2C::tdTomato ESCs. The flow cytometry results showed that only knocking down Pias4 significantly increased the fraction of 2C-like cells (Figs 1A and S1C). Consistently, quantitative PCR (qPCR) results showed that the expression of Dux and representative 2C-specific genes is up-regulated upon knocking down Pias4 ( Fig  1B). To confirm that the down-regulation of Pias4 is sufficient to globally promote the activation of zygotic transcriptional program, we performed RNA sequencing (RNA-seq) analysis on ESCs transfected with siRNAs against Pias4. The results showed that knocking down Pias4 globally increases the expression of 2C-specific ZGA transcripts (Fig 1C and S1 and S2 Tables). Importantly, a large fraction of Dux-induced genes or MERVL-long terminal repeat (LTR)driven genes were also up-regulated in Pias4-knockdown ESCs (Fig 1D and 1E). Consistently, Pias4 knocking down also up-regulated genes induced by miR-34a knockout [13], G9a knockout [5], or LINE1 [14] or CAF-1 (p150 and p60) [8] knockdown (S1D Fig.), all of which have been shown to promote the generation of 2C-like cells. Previously, using single-cell RNA-seq, Deng and colleagues analyzed gene expression profiles during different stages of mouse preimplantation development [32]. Interestingly, we found that genes induced upon knocking down Pias4 were also up-regulated in mouse embryos during the transitioning of zygote into the 2C stage, the time when ZGA occurs ( Fig 1F). These results suggest that Pias4 inhibits zygotic transcriptional program in mouse ESCs.
To exclude potential off-target effects of siRNAs, we generated Pias4-knockout ESCs (S1E Fig) using a clustered regularly interspaced short palindromic repeat (CRISPR)/CRISPR-associated protein 9 (Cas9) strategy and confirmed that knocking out Pias4 significantly increases the fraction of 2C-like cells as well as the expression of Dux and other 2C-specific genes (Figs 1G and 1H and S1F). Furthermore, we performed rescue experiment in Pias4−/− ESCs. We transfected Pias4−/− ESCs with wild-type or enzyme-dead mutant Pias4 (W356A) [33]. Western blotting analysis showed similar PIAS4 protein level between wild-type, Pias4, or Pias4 W356A rescued Pias4−/− ESCs (Fig 1I). qPCR results showed that wild-type but not enzymedead mutant Pias4 inhibited the expression of Dux and other 2C-specific genes ( Fig 1J). Altogether, these data indicate that Sumo2 and Pias4 inhibit ZGA transcripts in mouse ESCs.

PIAS4 protein is down-regulated in 2C-like cells and during ZGA
Since Pias4 is an inhibitor of 2C-like state and zygotic transcriptional program, we hypothesized that the expression of Pias4 is down-regulated in 2C-like cells and during ZGA. Unexpectedly, our single-cell qPCR results showed that Pias4 mRNA level was similar between 2C:: tdTomato-positive and negative cells (Fig 2A). The qPCR results are also confirmed by previous RNA-seq data (S2A Fig). Based on these results, we hypothesized that PIAS4 protein level is down-regulated in 2C-like cells. Indeed, immunofluorescence (IF) staining experiments showed that PIAS4 protein level is significantly lower in 2C::tdTomato-positive than negative cells (Fig 2B). Encouraged by these results, we further analyzed the Pias4 mRNA and protein level during early embryo development. Results from single-embryo quantitative reverse transcription PCR (RT-qPCR) showed that Pias4 mRNA was significantly decreased from zygote to the 2C embryo stage (Fig 2C), coincident with the up-regulation of ZGA transcripts including Dux, MERVL, and Zscan4d (S2B Fig). Consistent with qPCR results, IF staining showed that Pias4 protein was also significantly down-regulated during ZGA (Fig 2D). These results show that PIAS4 protein is down-regulated in 2C-like cells and during ZGA, consistent with its inhibitory role in zygotic transcriptional program. However, whereas both RNA and protein level of Pias4 are down-regulated during ZGA, only protein level of Pias4 is down-regulated in 2C-like cells, indicating different control mechanisms for Pias4 between early development in vivo and ESCs cultured in vitro.

Overexpression of Pias4 inhibits the zygotic transcriptional program and impairs early embryo development
Next, we tested whether Pias4 is sufficient to repress the zygotic transcriptional program by overexpressing Pias4 in ESCs (Fig 3A). qPCR results showed that around 2-fold overexpression (OE) of Pias4 could effectively decrease the expression of Dux and other 2C-specific genes ( Fig 3A). Consistent with the repression of zygotic transcriptional program, the percentage of 2C-like cells was also significantly decreased upon Pias4 OE (Fig 3B). To test whether Pias4 plays a similar inhibitory role during ZGA, we ectopically expressed Pias4 in mouse zygotes. In agreement with findings in ESCs, injecting in vitro transcribed green fluorescent protein (GFP)-Pias4 but not GFP or GFP-Pias4 W356A into zygotes induced embryonic arrest at the 2C stage (Fig 3C and 3D). Importantly, single-embryo RT-qPCR results showed that the up-regulation of Dux, MERVL, and Zscan4d were significantly inhibited in embryos injected with GFP-Pias4 at the early phase of the 2C stage ( Fig 3E). Moreover, Dux expression was only transiently activated at early 2C stage in GFP-injected embryos and quickly downregulated at the late 2C stage (Fig 3E). In contrast, the rapid down-regulation of Dux was not observed in Pias4-injected embryos (Fig 3E), suggesting that ZGA progression is severely interfered with by Pias4 OE. Together, these results demonstrate the potent function of Pias4 in inhibiting the transition of ESCs into 2C-like cells and ZGA processes during early embryonic development.

Knocking down Pias4 decreases the sumoylated protein form of multiple activators and inhibitors of zygotic transcriptional program
To figure out how PIAS4 inhibits zygotic transcriptional program, we decided to identify downstream substrate proteins of PIAS4. To this end, we performed pull-down of Sumo2 modified proteins followed by mass spectrometry (MS) analysis in control or Pias4 siRNAstransfected ESCs, based on an assumption that the abundance of Pias4 substrates in Sumo2 Left, representative dot plot; right, quantification of fraction of tdTomato-positive cells in 2C::tdTomato ESCs treated with control siRNAs or siRNAs against Pias4. Shown are mean ± SD, n = 6. The p-value was calculated by two-tailed Student's t test. (B) RT-qPCR of Dux and other 2C-specific genes in Pias4-knockdown ESCs. The β-actin gene was used as a control. For each gene, data were normalized to the mRNA level of ESCs transfected with control siRNAs. Shown are mean ± SD, n = 3. The p-value was calculated by two-tailed Student's t test. (C) GSEA for 2C-specific genes in ESCs transfected with control siRNAs or siRNAs against Pias4. For x-axis, genes were ranked based on the ratio of siNC versus si-Pias4 ESCs. (D) MA plots showing gene expression changes in Pias4-knockdown ESCs. Red dots indicate Dux-induced genes. Out of 123 Dux-induced genes, 75 were up-regulated in Pias4-knockdown ESCs. Fold enrichment and p-value are shown. The p-value was calculated by hypergeometric test. (E) MA plots showing gene expression changes in Pias4-knockdown ESCs. Red dots indicate MERVL-LTR-driven genes. Out of 161 MERVL-LTR-driven genes, 59 were up-regulated in Pias4-knockdown ESCs. Fold enrichment and p-value are shown. The p-value was calculated by hypergeometric test. (F) Expression of genes up-regulated in Pias4-knockdown ESCs in preimplantation mouse embryos. Center line, median; box limits, upper and lower quartiles; whiskers, 1.5× interquartile range. Preimplantation RNA-seq data are from [32]. (G) Fraction of 2C::tdTomato-positive cells in Pias4−/− ESCs. Left, representative fluorescence images. Right, flow cytometry analysis. Shown are mean ± SD, n = 3. The p-value was calculated by two-tailed Student's t test. Scale bars, 100 μm. (H) RT-qPCR of Dux and other 2C-specific genes in Pias4-knockout ESCs. The β-actin gene was used as a control. For each gene, data were normalized to the mRNA level of wild-type ESCs. Shown are mean ± SD, n = 3. The p-value was calculated by two-tailed Student's t test. (I) Western bots of PIAS4 protein level in wild-type ESCs and Pias4−/− ESCs rescued by empty, Pias4, or Pias4 W356A. Left, representative gel images; right, quantification of PIAS4 protein level, data were normalized to GAPDH and then to wild-type ESCs. Shown are mean ± SD, n = 3. (J) RT-qPCR of Dux and other 2C-specific genes in wild-type ESCs and Pias4−/− ESCs rescued by empty, Pias4, or Pias4 W356A. The β-actin gene was used as a control. For each gene, data were normalized to the mRNA level of wild-type ESCs. Shown are mean ± SD, n = 3. The p-value was calculated by one-way ANOVA followed by one-tailed Dunnett's test with empty versus wildtype and mutant Pias4 rescue. Source data for A, B, and G-J can be found in the supplemental data file (S1 Data pull-down should decrease upon Pias4 knockdown. We identified 97 potential candidate substrate proteins whose sumoylated form was decreased �1.25-fold in Pias4-knockout ESCs (S3A Fig and S3 Table). The decrease of sumoylated protein form might be caused by a decrease in sumoylation efficiency due to Pias4 knockout; alternatively, this could be due to Each dot represents one cell. n = 6. As expected, MERVL was significantly up-regulated in 2C:: tdTomato-positive cells. (B) IF staining of PIAS4 protein in 2C::tdTomato-positive and negative ESCs. Shown are representative images (left) and quantification of fluorescence intensity with mean ± SD (right). Scale bars, 20 μm. n = 77 tdTomato-positive cells and 1,478 tdTomato-negative cells. Each dot represents one cell. The p-value was calculated by two-tailed Student's t test. (C) Singleembryo RT-qPCR of Pias4 from zygote to the 8-cell stage. Shown is the Pias4 level normalized to zygote stage. Each dot represents one embryo. n = 6 for zygote, early 2C stage, and late 2C stage, and n = 4 for 4-cell and 8-cell stages. The p-value was calculated by two-tailed Student's t test. (D) IF staining of PIAS4 in mouse embryo from zygote to 8-cell stage. Top, representative images. Scale bars, 20 μm. Bottom, fluorescence intensity in nucleus. Shown are mean ± SD. Six to 10 embryos were analyzed at each stage. Each dot represents one nucleus. The p-value was calculated by one-way ANOVA followed by two-tailed Dunnett's test. Source data for A-D can be found in the supplemental data file (S1 Data

Dppa2 and Dppa4 are essential for activating zygotic transcriptional program in mouse ESCs
Next, we focused our study on Dppa2, a transcription factor expressed in germinal vesicle (GV) oocytes [35], whose dominant negative form arrests early mouse embryo development [36]. We transfected 2C::tdTomato reporter ESCs with siRNAs against Dppa2. The qPCR results showed that knocking down Dppa2 significantly decreased the expression of Dux and other 2C-specific genes (Fig 4A), consistent with the decrease in the percentage of 2C-like ESCs (Fig 4B). To confirm that Dppa2 is responsible for the up-regulation of zygotic transcriptional program upon Pias4 knockdown, we cotransfected siRNAs against Dppa2 and Pias4. The results show that knocking down Dppa2 can effectively prevent the up-regulation of Dux and other 2C-like ZGA genes upon Pias4 knockdown (S4A Fig). Since Dppa4 has been shown to play similar roles as Dppa2 in the activation of Dux and ZGA program [15,16], we also tested the effect of knocking down Dppa4 in ESCs. Indeed, qPCR results showed that knocking down Dppa4 decreased the expression of Dux and other 2C-specific genes at an extent similar to Dppa2 knockdown (S4B  Table). We found 50 2C-specific genes down-regulated �50% in double knockdown versus overexpression; Pias4, protein inhibitor of activated STAT 4; RT-qPCR, quantitative reverse transcription PCR; tdTomato, tandem dimeric Tomato; Zfp352, zinc finger protein 352; Zscan4d, zinc finger and SCAN domain containing 4D.
control. Among them, 31 genes were similarly down-regulated (differences within 20% range) in siDppa2, siDppa4, or siDppa2+siDppa4. Interestingly, six genes were further down-regulated �20% in double knockdown than any of the single knockdown. These results suggest that the function of Dppa2 and Dppa4 is largely cooperative in regulating the expression of 2C-specific genes, but another mode of action (e.g., redundant) likely exists. In addition, the overall profiles for knocking down Dppa2, Dppa4, or both were largely similar (S6A Fig). Intriguingly, Dppa2 and Dppa4 knockdown also down-regulated genes induced by other  [32]. Source data for A and B can be found in the supplemental data file (S1 Data). 2C, 2-cell; Dppa2, developmental pluripotency associated 2; Dux, double homeobox; ESC, embryonic stem cell; FC, fold change; FDR, false discovery rate; GSEA, gene set enrichment analysis; LTR, long terminal repeat; MERVL, murine endogenous retrovirus-L; NES, normalized enrichment score; RNA-seq, RNA sequencing; RT-qPCR, quantitative reverse transcription PCR; siNC, siRNA against negative control; siRNA, small interfering RNA; tdTomato, tandem dimeric Tomato; Zfp352, zinc finger protein 352; ZGA, zygotic genome activation; Zscan4, zinc finger and SCAN domain containing 4.

The sumoylation of DPPA2 by PIAS4 depends on K31 and K108 sites
Next, we verified that PIAS4 can catalyze the sumoylation of DPPA2 in mouse ESCs (Fig 5A). We overexpressed Flag-Dppa2 and hemagglutinin (HA)-Sumo2 in ESCs treated with siRNAs against negative control (siNCs) or Pias4. As expected, knocking down Pias4 decreased the fraction of sumoylated DPPA2 (Fig 5A). In addition, FLAG-DPPA2 and HA-SUMO2 were slightly but significantly increased in Pias4-knockdown ESCs (Fig 5B), suggesting that the decrease in the fraction of sumoylated DPPA2 is due to the decrease in sumoylation activity but not the expression of DPPA2 or HA-SUMO2 proteins in Pias4-knockdown ESCs. Furthermore, overexpressing PIAS4 increased the fraction of sumoylated DPPA2 in human embryonic kidney 293T (HEK293T) cells ( Fig 5C). To identify which lysine residues are important for the sumoylation of DPPA2 by PIAS4, we decided to map the potential sumoylation sites on DPPA2 protein. Using the in silico prediction programs SUMOplot (Abgent) and GPS-SUMO (The Cuckoo Workgroup), eight lysine residues (K31, K67, K108, K137, K159, K160, K226, and K291) were predicted as possible candidate sites for SUMO-conjugation in DPPA2 ( Fig 5D). Mutant proteins for all predicted lysine residues were then generated, and their effects on sumoylation of DPPA2 were examined in HEK293T cells. K31R and K108R mutations led to a substantial reduction in the intensity of the higher-molecular-mass sumoylated DPPA2 bands ( Fig 5E). These results suggest that the sumoylation of DPPA2 by PIAS4 depends on K31 and K108 sites.

Sumo2-Dppa2 inhibits the expression of 2C-specific genes and impairs early embryo development
To figure out the impact of sumoylation by Pias4 on the protein level of DPPA2, we performed western analysis of DPPA2 protein in control siRNA or ESCs treated with Pias4 siRNA. The results showed that the DPPA2 protein was only slightly increased upon Pias4 knockdown ( Fig 6A), indicating that Pias4 may regulate zygotic transcriptional program through regulating the sumoylation but not the protein level of DPPA2. In order to test the consequence of sumoylation of DPPA2 by PIAS4, a DPPA2 fused with SUMO2ΔGG at the N terminal was ectopically expressed in ESCs to mimic the sumoylated DPPA2 [37,38]. To avoid artifacts due to high expression level, SUMO2ΔGG-DPPA2 was expressed at a significantly lower level than endogenous DPPA2 protein (S7A and S7B Fig). Interestingly, we found that expressing SUMO2ΔGG-DPPA2 at this level can significantly decrease the expression of 2C-specific genes as well as the fraction of 2C-like cells (Fig 6B and 6C Since sumoylated DPPA2 has an inhibitory role for the zygotic transcriptional program, we reasoned that sumoylated DPPA2 must decrease in 2C::tdTomato-positive cells. To test this, we performed a proximity ligation assay (PLA) [39] using antibodies against SUMO2 and DPPA2. Indeed, the PLA assay results showed that the signal of colocalization of SUMO2 and DPPA2 is lower in 2C::tdTomato-positive cells (S7G Fig). Meanwhile, the protein level of DPPA2 is significantly increased in 2C::tdTomato-  Table of the predicted SUMOylation sites using in silico prediction programs SUMOplot and GPS-SUMO in DPPA2. Scores ware based on two criteria: direct amino acid match to SUMO consensus motif (ψKXE) or substitution of the consensus amino acid residues with amino acid residues exhibiting similar hydrophobicity. (E) Western blotting analysis of Flag IP samples in HEK293T cells transfected with flag-tagged lysine mutant, 3XHA-Sumo2, and Pias4. Quantification for relative sumoylation level is shown below. Data were normalized to HEK293T cells transfected with flag-tagged WT Dppa2, 3XHA-Sumo2, and Pias4. Shown are mean ± SD, n = 3. The p-value was calculated by one-way ANOVA followed by two-tailed Dunnett's test. Source data positive cells with no apparent alteration of localization (S7H and S7I Fig). Finally, we tested the inhibitory role of sumoylated DPPA2 for zygotic transcriptional program in mouse embryos. Compared to GFP and GFP-T2A-DPPA2, the expression of GFP-T2A-SUMO2ΔGG-DPPA2 in zygote significantly impaired the early embryo development (Fig 6E  and 6F), accompanied with insufficient up-regulation of 2C-specific genes at the early phase of the 2C stage ( Fig 6G). Therefore, sumoylated DPPA2 likely plays a dominant negative function to DPPA2. Taken together, these data showed that PIAS4 inhibits zygotic transcriptional program at least partially through sumoylation of DPPA2.

OE of Dppa2 activates zygotic transcriptional program
We then tested whether OE of Dppa2 is sufficient to activate zygotic transcriptional program. The flow cytometry analysis showed that around 4-fold OE of Dppa2 ( Fig 7A) significantly increased the fraction of 2C::tdTomato-positive cells (Fig 7B). In addition, Dppa2 OE further increased the fraction of 2C-like cells in Sumo2-or Pias4-knockdown ESCs (S8A Fig). Importantly, qPCR results showed that Dppa2 significantly activated the expression of Dux (around 7-fold) and other 2C-specific genes (Fig 7C). To find out the global impact of Dppa2 OE on zygotic transcriptional program, we performed RNA-seq in Dppa2-overexpressing ESCs. The results showed that Dppa2 OE globally up-regulated the expression of 2C-specific ZGA transcripts as well as Dux-induced genes and MERVL-LTR-driven genes (Fig 7D-7F and S5 and S6 Table). In addition, Dppa2 OE also up-regulated genes induced by CAF-1 knockdown, LINE1 knockdown, G9a knockout, and mir-34a knockout (S8B Fig). Moreover, genes induced by Dppa2 OE were also up-regulated in mouse embryos during the transitioning of zygote [32] into the 2C stage, when ZGA takes place ( Fig 7G). Interestingly, PCA analysis showed that Dppa2-OE or Pias4-knockdown ESCs clustered closely with other 2C-like cells including P150-, P60-, or LINE1-knockdown ESCs (S8C Fig). To gather more evidence that Dppa2 OE and Pias4−/− ESCs have properties of 2C-like cells, we performed chimera formation by aggregating ESCs with 8-cell embryos. The results showed that Dppa2 OE and Pias4−/− ESCs produced chimeric blastocysts with incorporation into both inner cell mass (ICM) and trophectoderm (TE) at a significantly higher frequency than wild-type ESCs ( Fig 7H). Together, these results indicate that overexpressing Dppa2 is sufficient to activate the expression of ZGA transcripts.
Next, we tested whether Dppa4 OE can activate the zygotic transcriptional program. Interestingly, qPCR results showed that Dppa4 alone cannot activate the expression of Dux and other 2C-specific genes. In addition, OE of Dppa4 slightly increased the expression of Dux and other 2C-specific genes in Dppa2-overexpressing ESCs (S8D Fig). However, the fraction of 2C-like cells was significantly increased in Dppa2/4-overexpressing versus Dppa2-overexpressing ESCs (S8E Fig). Together, these results support that OE of Dppa2 activates the zygotic transcriptional program and Dppa4 can enhance the function of Dppa2.

Dppa2 activates zygotic transcriptional program through directly activating Dux
To figure out how Dppa2 activates zygotic transcriptional program, we decided to identify the downstream targets of Dppa2. We searched through a previously published chromatin for A, B, and E can be found in the supplemental data file (S1 Data). BFP, blue fluorescent protein; Dppa2, developmental pluripotency associated 2; ESC, embryonic stem cell; GAPDH, glyceraldehyde 3-phosphate dehydrogenase; HA, hemagglutinin; HEK293T, human embryonic kidney 293T; IP, immunoprecipitation; PIAS4, protein inhibitor of activated STAT 4; siNC, siRNA against negative control; siRNA, small interfering RNA; SUMO, small ubiquitin-like modifier; WT, wild type.
https://doi.org/10.1371/journal.pbio.3000324.g005 DPPA2/4 and PIAS4 opposingly regulate zygotic transcriptional program ESCs treated with control and Pias4 siRNAs. Data are normalized to β-actin and wild-type siNC. Shown are presented as mean ± SD, n = 3. The pvalue was calculated by two-way ANOVA followed by two-tailed Dunnett's test. (E) Development of zygote injected with GFP, GFP-t2a-Dppa2, and GFP-t2a-Sumo2ΔGG-Dppa2 mRNAs. Shown is the percentage of normally developed embryos. n = 58, 76, and 61 embryos for GFP, GFP-t2a-SM2-Dppa2, and GFP-t2a-Dppa2 from four experiments, respectively. The p-value was calculated by Pearson's chi-squared test. (F) Development of zygotes that were injected with in vitro transcribed GFP, GFP-t2a-Dppa2, and GFP-t2a-Sumo2ΔGG-Dppa2 mRNAs. Shown are immunoprecipitation sequencing (ChIP-seq) data of Dppa2 in mouse ESCs [40] and found no enrichment of Dppa2 binding on 2C genes globally. Interestingly, though, we found that the gene body of Dux was bound by Dppa2 (Fig 8A), raising the possibility that Dppa2 regulates zygotic transcriptional program through directly activating Dux. We then confirmed the binding of Dppa2 to Dux gene body by ChIP-qPCR (Fig 8B). Consistently, Dppa2 OE significantly increased the expression of Dux (Fig 8C). Moreover, knocking down Dux significantly inhibited the effects of Dppa2 OE in activating 2C-specific genes ( Fig 8C) and increasing the fraction of 2C-like cells (Fig 8D). Taken together, these results suggest that Dppa2 activates zygotic transcriptional program through directly up-regulating Dux.

Discussion
The transcription factor DUX resides at the top of a transcriptional hierarchy to activate zygotic transcriptional program during the 2C stage [10][11][12]. However, DUX itself is only beginning to be expressed at the early phase of ZGA; therefore, other transcriptional factors activating the transcription of Dux must exist. Here, we confirmed recent findings [15,16] that two related transcription factors Dppa2 and Dppa4 activate zygotic transcription program through directly activating Dux (Fig 8E). Both Dppa2 and Dppa4 are required for the activation of ZGA transcripts in 2C-like cells. Whereas Dppa2 OE sufficiently activates zygotic transcriptional program, Dppa4 OE is not sufficient and only enhances the function of Dppa2. In addition, we identified a Sumo2 E3 ligase Pias4, which inhibits the zygotic transcriptional program, likely through sumoylating multiple inhibitors and activators in regulating zygotic transcription. We further showed that the function of Pias4 is at least partially through sumoylating DPPA2, which turns DPPA2 from an activator to a potent inhibitor of zygotic transcriptional program. Consistent with its essential function as a brake for zygotic transcriptional program, the PIAS4 protein level was down-regulated at the early phase of ZGA and in 2C-like ESCs. However, the exact mechanism in regulating Pias4 expression is different between ZGA and 2C-like ESCs. Future investigation is needed to clarify this discrepancy. Furthermore, Pias4 OE or sumoylated DPPA2 inhibited the emergence of 2C-like cells and impaired early embryo development in mouse. Together, our study identifies key molecular rivals upstream of Dux in regulating the transition of ESCs into 2C-like cells and zygotic transcriptional program.
https://doi.org/10.1371/journal.pbio.3000324.g006 DPPA2/4 and PIAS4 opposingly regulate zygotic transcriptional program  [31]. However, which E3 ligase is responsible for regulating ESC identity was unknown. In our study, we identified Pias4 as the specific E3 ligase repressing the transition of ESCs into 2Clike cells. We noticed that Pias4 deletion leads to decreased viability perinatally [46,47], consistent with its important function during embryogenesis. However, the survived Pias4−/− animals are largely normal and fertile. Therefore, alternative pathways may exist to repress zygotic transcriptional program. After fertilization, the maternal Pias4 mRNA is decreased at the early phase of the 2C stage. In contrast, during 2C-like state transition in ESCs, Pias4 mRNA level remains unchanged, but protein is dramatically down-regulated. The exact mechanisms underlying the regulation of Pias4 expression in these two different circumstances warrant future investigation.
Dppa2 and Dppa4 are closely linked genes tandemly locating on the same chromosome in both mouse and humans [35]. They both contain a SAP domain that is known as DNA binding domain [48]. In addition, Dppa4 was reported to directly bind to chromatin, particularly transcriptionally active regions [49]. More recently, Dppa2 and Dppa4 were reported to be upregulated during chemically induced reprogramming processes when 2C-like program is activated [50]. In another study, Dppa2 and Dppa4 were found as accelerators during reprogramming to pluripotency [51]. Mechanistically, they work as a heterodimer bound to chromatin to initiate global chromatin decompaction, which also happens to be an important event during 2C-like cell transition in ESCs and ZGA at the 2C stage [6,11,52]. These studies suggest that besides activating Dux expression, Dppa2 and Dppa4 may also facilitate chromatin remodeling to promote zygotic transcriptional program. Dppa2 and Dppa4 are expressed throughout the preimplantation development and in ESCs [35,53], but Dux is only activated at the 2C stage and in a minor population of ESCs. This controversy strongly suggests that additional factors must exist to restrain the activating function of Dppa2/4. Our study suggests that sumoylation  https://doi.org/10.1371/journal.pbio.3000324.g008 DPPA2/4 and PIAS4 opposingly regulate zygotic transcriptional program required for Dppa2/4, and the expression or activity of these factors is restrictive in the 2C stage and 2C-like cells. Furthermore, the accessible chromatin state of Dux gene may be a prerequisite for Dppa2/4 to function. It is worth noting that Dppa2 and Dppa4 single or double knockout did not affect early embryogenesis [54,55], possibly because of maternally deposited products. Dppa2-and Dppa4-knockout mice died early after birth, preventing further evaluation of fertility of these mice. Therefore, maternal knockout of Dppa2 and Dppa4 is required to provide definite evidence on their function in ZGA. However, injection of a Dppa2 lacking SAP domain [36] or a sumoylated Dppa2 (this study) impaired preimplantation development.
The dominant negative effect of mutant Dppa2 suggests that even though the redundant factors exist, they must share a similar working mechanism with Dppa2. Elucidating these mechanisms or identifying additional factors will contribute to the complete understanding of molecular mechanisms underlying ZGA control and may potentially lead to efficient means to produce a large quantity of 2C-like cells for research and therapeutic applications.

Ethics statement
Mice were bred and maintained under specific pathogen-free (SPF) conditions in the institutional animal facilities at Peking University and Tsinghua University. All animal protocols were approved by Institutional Animal Care and Use Committees (IACUC) of Peking University (IMM-WangYM-1) and Tsinghua University (16-NJ-1), both of which are accredited by the Association for Assessment and Accreditation of Laboratory Animal Care International (AAALAC).

Cell culture and construct of reporter cell lines
The ESC culture medium consisted of KnockOut DMEM (Gibco, Cat. # 10829081) with 15% FBS (Hyclone, Cat. # SH3007103), 1,000 U/ml mouse leukemia inhibitory factor (1,000 U/ml), 0.1 mM nonessential amino acids (Gibco, Cat. # 11140050), 0.1 mM β-mercaptoethanol, 1 mM L-glutamine (Gibco, Cat. # 25030081), and penicillin (100 U/ml) and streptomycin (100 μg/ml). For culture of ESC lines, the medium was changed daily, and cells were routinely passaged every other day. For generation of 2C::tdTomato reporter cell line, the MERVL-LTR-tdTomato reporter constructs made as previously described [5] were linearized and transfected into ESCs by electroporation. The cells were then selected with G418 for 7 d. Colonies containing tdTomato-positive cells were then picked and expanded. All cell lines were kept under constant drug selection with G418 to prevent transgene silencing. Mycoplasma detection tests were conducted routinely to ensure mycoplasma-free conditions throughout the study.

RNA extraction and RT-qPCR
Total RNA was extracted following standard Trizol protocol (Invitrogen, Cat. # 15596026). For RT-qPCR analysis, complementary DNA (cDNA) were generated from 500 ng total RNA using HiScript II Q RT SuperMix for qPCR (Vazyme, Cat. # R223). qPCR was performed using AceQ qPCR SYBR Green Mater Mix (Vazyme, Cat. # Q141) in 96-well dishes in three biological replicates on StepOnePlus Real-Time PCR System (Applied Biosystems) with standard protocols. The expression levels were plotted relative to b-actin. Primers for qPCR are listed in S7 Table.
For CRISPRi screening, a single gRNA was designed by http://crispr-era.stanford.edu/ for each candidate gene. The gRNA was cloned into a plasmid containing dCas9-krab. For screening, 0.05 million ESCs were plated in 12-well plate, and then 1 μg plasmid containing gRNA and dCas9-krab was transfected using lipofectamine 3000 (Invitrogen, Cat. # L3000) following the manufacturer's protocol. ESCs were selected with 300 μg/ml hygromycin 24 h post transfection for 4-5 d. Cells were then collected for flow cytometry analysis and RT-qPCR experiments.

siRNA transfection
ESCs were transfected with siRNA (S8 Table) using DharmaFECT 1 (Dharmacon, Cat. # T-2001) reagent following the manufacturer's instructions. NC is a siRNA sequence derived from Caenorhabditis elegans genome that does not target any mammalian genes, provided by GenePharma (Shanghai, China) as a control for siRNA transfection. Typically 2 days after transfection, cells were collected for flow cytometry, immunoblot, and RNA extraction.

Sumoylation IP-MS
His6-Sumo2 pull-down was performed as previously described [56]. Briefly, cells were washed in PBS and lysed in lysis buffer (6 M guanidinium-HCl, 10 mM Tris, 100 mM sodium phosphate buffer [pH 8.0], 5 mM β-mercaptoethanol, 1 mM imidazole). Ni2+ NTA magnetic agarose beads (50 μl, Qiagen) were added to cell lysis and incubated at 4˚C overnight. After incubation, the beads were washed once in lysis buffer, once in wash buffer (pH 8.0) (8 M urea, 10 mM Tris, 100 mM sodium phosphate buffer [pH 8.0], 0.1% Triton X-100, 5 mM βmercaptoethanol), and three times in wash buffer (pH 6.3) (8 M urea, 10 mM Tris, 100 mM sodium phosphate buffer [pH 6.3], 0.1% Triton X-100, 5 mM β-mercaptoethanol, 10 mM imidazole). Sumoylated proteins were eluted from the beads using elution buffer (100 mM sodium phosphate buffer [pH 6.8], 200 mM imidazole) for MS analysis. Immunoprecipitated proteins were analyzed by SDS-PAGE. For LC-MS/MS analysis, the eluted peptides were sprayed into a Velos Pro Orbitrap Elite mass spectrometer (Thermo Scientific, USA) equipped with a nano-ESI source. The mass spectrometer was operated in data-dependent mode with a full MS scan in FT mode at a resolution of 120,000 followed by collision-induced dissociation (CID) MS/MS scans on the 15 most abundant ions in the initial MS scan. To identify Pias4 substrate proteins, we compared the count of distinct sequences that have significant scores (prot_sequences_sig) for each protein in siNC versus si-Pias4 samples. Candidate proteins selected for further CRISPRi screening were chosen based on the following criteria: The candidate has to be detected in at least two of four siNC samples; the average #PSM for a candidate must be >1 in siNC samples; the ratio of average #PSM of siPias4 versus siNC is �1.25-fold at the 48-hr or 72-hr time point. In total, 97 candidates passed these criteria. For screening, housekeeping genes Actb and Ctc1 were not included; Mtco2 and Tmlhe were not included, as no appropriate CRISPRi gRNA can be designed based on their sequences; Sumo1 was not included, as we already showed that knocking down Sumo1 decreased the expression of Dux and other 2C-like genes (S1B Fig).

Flow cytometry
Cells were trypsinized and resuspended in ice-cold PBS containing 2% FBS. Sorting was performed using a BD FACSAria III. During sorting, cells were collected in culture medium and kept at 4˚C. Quantifying the population of 2C::tdTomato-positive cells was performed by BD LSRFortessa SORP. Data were analyzed using FlowJo software.

In vitro transcription
The cDNA encoding the desired genes was amplified and cloned under the control of a T7 promoter. After linearization by a restriction enzyme PmeI, the construct was purified with phenol-chloroform extraction and ethanol precipitation. mRNA was synthesized by in vitro transcription using a HiScribe T7 ARCA mRNA Kit (NEB, Cat. # E2065) according to the manufacturer's instructions. The synthesized mRNA was purified by lithium chloride precipitation and diluted with nuclease-free water. mRNA aliquots were stored at −80˚C until use.

Chimeric blastocyst assay
The zona pellucida of 8 cell-stage embryos were removed by a short exposure to acidic Tyrode's solution (Sigma, Cat. # T1788). The denuded embryos were placed into each concaved microwell created by a smooth depression using the aggregation needle. Four ESCs labeled by mRuby2 fluorescence were transferred into each concaved microwell and cocultured with the denuded embryo overnight.

IF staining and PLA
Cells were fixed with 4% paraformaldehyde for 20 min at room temperature. After the fixation, cells were permeabilized with 0.25% Triton X-100 for 20 min at room temperature and blocked with 3% FBS in PBS for 1 h at room temperature. Cells were then incubated with primary antibodies (1:100, anti-PIAS4, Proteintech, Cat. # 14242-1-AP) diluted in PBS with 3% FBS for 2 h. After washing three times in PBS, the cells were incubated with secondary antibody (1:200, anti-rabbit IgG Alexa fluor 488) for 1 h and followed by DAPI staining. For preimplantation embryos' IF staining, 0.1% Tween-20 was added in PBS. For PLA, the permeabilized cells were incubated in PBS containing two primary antibodies (1:700, anti-Dppa2, Abcam, Cat. # ab91318; 1:500, anti-Sumo2, Cytoskeleton, Cat. # ASM23) followed by Duolink in situ PLA (Sigma-Aldrich) anti-mouse (minus) and anti-rabbit (plus) probes and detection reagents Green according to manufacturer's instructions.

RNA-seq and bioinformatics analysis
Total RNA was enriched twice with poly-T oligo-attached magnetic beads and then subjected to the synthesis of double-stranded (ds) cDNA. The ds-cDNA was ligated to adaptors from NEB and sequenced by Illumina Genome Analyzer (Novogene, Tianjin, China). Sequencing reads were aligned to the mouse genome (mm10) with STAR (version 2.5.0) using the GEN-CODE transcript annotation as transcriptome guide. All programs were processed following default settings except for special annotation. The FPKM value generated by Cufflinks was used to quantify the expression level. The enrichment of selected gene sets was calculated by java GSEA Desktop Application. R 3.1.1 and Matlab were used for the generation of scatter plot and boxplot. The 2C-specific ZGA genes are genes activated during ZGA (the 2C stage) that are also enriched in 2C::tdTomato + cells from [5]. The Dux-induced gene set was from [11]. The list of MERVL-LTR-driven transcripts was from [6]. The list of genes induced by mir-34a knockout was from [13]. The list of genes induced by G9a knockout was from [5]. The list of genes induced by LINE1 knocking down was from [14]. The lists of genes induced by p150 or p60 knocking down were from [8].

Quality control of sequencing data
The quality of every library was determined using the fastqc tool (http://www.bioinformatics. babraham.ac.uk/projects/fastqc/). Reads were subsequently trimmed and adapters clipped using the fastq-mcf (https://github.com/ExpressionAnalysis/ea-utils/blob/wiki/FastqMcf.md). Only reads with none of the known high-throughput sequencing adapters, longer than 25 bps, with a mean quality score above 30 and maximum 1 N-call were kept.

Quantification and statistical analysis
The number of independent experimental replications, the definition of center, and precision measures are reported in the figure legends (n, mean ± SD). p < 0.05 is generally considered as statistically significant. Statistical analyses were performed using the GraphPad Prism v6 software. Statistical significance was assessed by two-tailed t test except when specified in the figure legends. For boxplots, upper and lower whiskers are defined as respectively Q3 + 1.5 × IQR and Q1 − 1.5 × IQR, with Q1 and Q3 being the first and third quartile of the plotted distribution and IQR the interquartile range; the p-value was determined by Wilcoxon signed rank test. For multiple comparison, the p-value was calculated by one-way or two-way ANOVA followed with Dunnett's test. Shown are mean ± SD, n = 3. The p-value was calculated by two-way ANOVA followed by two-tailed Dunnett's test. (B) RT-qPCR of Dppa2 and Dppa4 (left), Dux, and other 2C-specific genes (right) in ESCs treated with siRNAs against Dppa2 and Dppa4 individually or in combination. The β-actin gene was used as a control. For each gene, data were normalized to the mRNA level of wild-type ESCs. Shown are mean ± SD, n = 3. The p-value was calculated by two-way ANOVA followed by two-tailed Dunnett's test. Source data for A and B can be found in the supplemental data file (S1 Data). 2C, 2-cell; Dppa, developmental pluripotency associated; Dux, double homeobox; ESC, embryonic stem cell; Pias4, protein inhibitor of activated STAT 4; RT-qPCR, quantitative reverse transcription PCR; siRNA, small interfering RNA. Dppa2-overexpressing ESCs treated with Sumo2 or Pias4 siRNAs. Data are presented as mean ± SD, n = 3. The p-value was calculated by two-tailed Student's t test. (B) Box-and-whisker plots showing expression of genes up-regulated by mir-34a KO, G9a KO, LINE1 knockdown, and Caf-1 p150 or p60 subunit knockdown in cells overexpressing Dppa2. The p-value was determined by Wilcoxon signed rank test. (C) PCA mapped scatter plot: global protein coding genes (left) and repeat elements (right). Data for G9a KO from [5], P150 and P60 knockdown and 2C::EGFP+ cells from [8], Dux-overexpressing cells from [11], and LINE1-knockdown cells from [14]. Data were normalized to the control cell line of each study to exclude batch effects before PCA processing. (D) RT-qPCR of Dux and other 2C-specific genes in control, Dppa2-, Dppa4-, or Dppa2/4-overexpressing ESCs. The β-actin gene was used as a control. For each gene, data were normalized to the mRNA level of control ESCs. Shown are mean ± SD, n = 3. The p-value was calculated by one-way ANOVA followed by two-tailed Dunnett's test. (E) Fraction of 2C::tdTomato-positive cells in control, Dppa2-, Dppa4-, or Dppa2/4-overexpressing ESCs. Data are presented as mean ± SD, n = 3. The pvalue was calculated by two-tailed Student's t test. Source data for A, D, and E can be found in the supplemental data file (S1 Data). 2C, 2-cell; Dppa, developmental pluripotency associated; Dux, double homeobox; EGFP, enhanced green fluorescent protein; ESC, embryonic stem cell; KO, knockout; LINE1, long interspersed nuclear element; mir-34a, microRNA 34a; PCA, principle component analysis; Pias4, protein inhibitor of activated STAT 4; RT-qPCR, quantitative reverse transcription PCR; siRNA, small interfering RNA; Sumo2, small ubiquitin-like modifier 2; tdTomato; tandem dimeric Tomato. (TIF) S1 Table. Table. Sequences of ChIP-qPCR primers. Sequence information of forward and reverse primers for ChIP-qPCR analysis of various positions in Dux promoter and gene body. ChIP-qPCR, chromatin immunoprecipitation followed by quantitative polymerase chain reaction. ChIP-qPCR, chromatin immunoprecipitation-quantitative PCR; Dux, double homeobox. (XLSX) S1 Data. Excel spreadsheet containing, in separate sheets, the underlying numerical data and data analysis for figure panels. Source data for main and supplemental figures in separate sheets. Mean, SD, and p-value are shown where applicable. (XLSX)