Integrated Single-Cell analysis Reveals molecular correlation of maternal-to-zygotic transition between human and pig embryo.

Skip to main content An official website of the United States government Here's how you know Here's how you know Official websites use .gov A .gov website belongs to an official government organization in the United States. Secure .gov websites use HTTPS A lock ( Lock Locked padlock icon ) or https:// means you've safely connected to the .gov website. Share sensitive information only on official, secure websites. Search Log in Dashboard Publications Account settings Log out Search… Search NCBI Primary site navigation Search Logged in as: Dashboard Publications Account settings Log in Search PMC Full-Text Archive Search in PMC Journal List User Guide PERMALINK Copy As a library, NLM provides access to scientific literature. Inclusion in an NLM database does not imply endorsement of, or agreement with, the contents by NLM or the National Institutes of Health. Learn more: PMC Disclaimer | PMC Copyright Notice J Transl Med . 2026 Mar 4;24:500. doi: 10.1186/s12967-026-07904-w Search in PMC Search in PubMed View in NLM Catalog Add to search Show available content in en fr Integrated Single-Cell analysis Reveals molecular correlation of maternal-to-zygotic transition between human and pig embryo Jianlin Fan Jianlin Fan 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China Find articles by Jianlin Fan 1, # , Qinjian Li Qinjian Li 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Qinjian Li 1, 2, # , Zhongyu Yuan Zhongyu Yuan 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Zhongyu Yuan 1, 2, # , Xueqing Liu Xueqing Liu 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Xueqing Liu 2 , Chaoqian Jiang Chaoqian Jiang 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Chaoqian Jiang 2 , Xiaokang Xu Xiaokang Xu 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Xiaokang Xu 1, 2 , Zhonghua Liu Zhonghua Liu 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Zhonghua Liu 1, 2, ✉ , Yanshuang Mu Yanshuang Mu 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China Find articles by Yanshuang Mu 1, 2, ✉ Author information Article notes Copyright and License information 1 Key Laboratory of Animal Cellular and Genetic Engineering of Heilongjiang Province, Northeast Agricultural University, Harbin, China 2 College of Life Science, Northeast Agricultural University, Harbin, China ✉ Corresponding author. # Contributed equally. Received 2025 Oct 4; Accepted 2026 Feb 13; Collection date 2026. © The Author(s) 2026 Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/ . PMC Copyright notice PMCID: PMC13069747  PMID: 41782053 Abstract Background Early embryonic development begins with a transcriptionally quiescent phase that is regulated by maternally deposited products. Developmental control subsequently shifts to the zygotic genome through zygotic genome activation (ZGA) during the maternal-to-zygotic transition. This transition represents a critical developmental window that governs embryonic timing and cell fate determination. Notably, the characteristics and temporal dynamics of ZGA occur in distinct waves and display substantial divergence across species. Methods Single-cell RNA sequencing data from 139 human and porcine samples, spanning developmental stages from the oocyte to the morula, were integrated using Seurat v5 with reciprocal principal component analysis. Downstream analyses included orthologous gene mapping, dimensionality reduction, differential gene expression analysis, functional enrichment analysis, and pseudotime trajectory inference using Monocle3. Transcription factor (TF) regulatory networks were reconstructed using SCENIC (v1.3.1). Results Our analyses revealed both conserved and species-specific features during embryonic development from the oocyte to the morula stage in humans and pigs. The core TFs associated with major ZGA events differed markedly between the two species, indicating species-specific regulatory programs. In contrast, chromatin remodeling–related genes were highly enriched among maternal transcripts and exhibited functional conservation across species. Furthermore, MT genes showed consistently higher average expression levels in human embryos compared with porcine embryos across all examined stages, whereas ATP8 expression demonstrated pronounced heterogeneity among pig embryos. Conclusion In conclusion, early embryonic development in humans and pigs is governed by conserved regulatory mechanisms underlying key biological processes, such as chromatin remodeling and major ZGA events. However, notable differences in metabolic activity and signaling pathways suggest the presence of species-specific adaptive strategies during early development. Supplementary Information The online version contains supplementary material available at 10.1186/s12967-026-07904-w. Graphical Abstract Introduction Embryonic development begins with a single cell that contains the oocyte cytoplasm and the combined parental genomes, which merge to form a unified nucleus. During the earliest stages, the embryo is transcriptionally inactive and relies exclusively on maternally derived RNAs and proteins to support initial developmental processes. As development proceeds, regulatory control shifts from maternal factors to the newly formed zygotic genome during the maternal-to-zygotic transition (MZT) [ 1 , 2 ]. This transition encompasses two fundamental processes that reprogram the embryo and reestablish totipotency [ 3 ]: the selective degradation of maternal RNAs and proteins that are no longer required [ 4 , 5 ], and the activation of previously silent zygotic genes through zygotic genome activation (ZGA). ZGA initiates transcription from the embryonic genome and involves key transcription factors (TFs), such as members of the TPRX family in humans and OBOX and DUX proteins in mice, which together establish species-specific gene regulatory networks [ 6 , 7 ]. Mammalian ZGA occurs in two sequential waves: minor ZGA and major ZGA. Although minor ZGA involves relatively low levels of transcriptional activity, it is essential for priming the activation of genes expressed during major ZGA [ 8 ]. The timing and molecular features of these activation waves vary across species. In mice, minor ZGA initiates at the one-cell stage and is followed by major ZGA at the two-cell stage [ 9 ]. In contrast, pigs and humans display more comparable ZGA dynamics: in pigs, major ZGA takes place between the four- and eight-cell stages [ 10 – 12 ], whereas in humans it predominantly occurs at the eight-cell stage [ 9 ]. During ZGA, a wide range of fundamental biological processes are initiated, including the clearance of maternal regulatory programs, activation of the embryonic genome, establishment of stage-specific transcriptional regulatory networks, extensive epigenetic remodeling, RNA turnover, and metabolic reprogramming [ 6 ]. Prior to ZGA, the zygotic genome remains largely transcriptionally inactive; following activation, embryonic development proceeds under the control of the embryo’s own genetic information [ 6 , 7 , 13 ]. This developmental window is therefore critical, as it ultimately determines both the tempo of embryonic progression and subsequent cell fate decisions [ 6 , 7 , 13 , 14 ]. Pigs offer significant advantages for advancing human reproductive research due to their short generation intervals, large multiparous litters, high efficiency of genome editing, and their established utility as donors for human xenotransplantation [ 15 – 18 ]. Porcine oocytes are comparable to human oocytes in both size and proteomic expression profiles, and the two species share similar mechanisms of early embryo production [ 19 – 21 ]. However, despite these similarities, post-fertilization embryonic development differs substantially between pigs and humans [ 7 , 13 , 14 ]. Cross-species comparisons of single-cell transcriptomes from early embryos provide a powerful approach for elucidating the molecular foundations underlying both conserved and species-specific developmental processes [ 7 , 13 ]. Characterizing these similarities and differences yields critical insights into the regulatory mechanisms governing early mammalian development [ 7 , 13 , 22 ]. Moreover, such knowledge strengthens the utility of pigs as translational models in human medicine, supports the optimization of xenogeneic organ donation strategies, and facilitates the generation of human–pig chimeras [ 23 ]. In this study, we integrated single-cell RNA sequencing datasets from human and porcine embryos covering developmental stages from the oocyte to the morula in order to systematically characterize cross-species transcriptomic dynamics. This comparative analysis aimed to identify both conserved regulatory factors and species-specific gene expression programs operating during early embryogenesis. By delineating shared and divergent molecular mechanisms, our findings contribute to a deeper understanding of early mammalian development and offer valuable insights relevant to reproductive biology, developmental biology, and emerging biotechnological applications. Results Sample selection and primary data Characterization To investigate the relationship between humans and pigs during the MZT in early embryogenesis, we integrated single-cell RNA sequencing data from a total of 139 samples encompassing developmental stages from the oocyte to the morula in both species. Orthologous genes were identified using the Ensembl database, and these genes consistently exhibited high and comparable expression proportions across developmental stages in both species, exceeding approximately 75–85% (Fig. 1 A). Within the curated gene sets, orthologous genes accounted for 66–76% of expressed genes in human embryos and 88–93% in porcine embryos (Fig. 1 B). Based on these results, orthologous genes were selected for subsequent comparative analyses. For data preprocessing, transcript abundance was normalized using log-transformed transcripts per million (TPM) values. The normalized data were then scaled using the ScaleData function in Seurat with default parameters (mean centering and unit variance scaling). Dimensionality reduction was performed by retaining the first 50 principal components, which were used for downstream clustering, visualization, and reconstruction of developmental trajectories. Fig. 1. Open in a new tab Comparative characteristic analysis of single-cell RNA sequencing data for human and pig embryos from oocyte-morula stage. ( A ) column charts show homologous gene expression (in TPM) proportions in the original gene expression matrices for human and pig embryos, from oocyte to morula stage. ( B ) column charts compare homologous genes’ proportions in original vs. selected matrices (TPM > 0/1 at oocyte or 8-cell stages) for each species. ‘o’ indicates original matrices; ‘0’/‘1’ indicates selected matrices. ( C ) 3D UMAP of gene expression after dimensionality reduction for 139 samples, using robust principal component analysis (RPCA) integration, with blue pseudo-time trajectory showed for human and red for pig. ( D ) box plots compare gene amounts between species at each stage. ( E ) box plots compare total expression between species at each stage. ( F ) box plots compare expression percentages of MT encoded genes between species at each stage. ( G ) expressions trajectories of POU5F1 , SOX2 , KLF4 , NANOG , MYC across stages in both species, respectively. ( H ) heatmap depicts top 5 upregulated markers’ expression at each stage for species Dimensionality reduction based on the full gene expression matrix revealed pronounced interspecies differences during early developmental stages up to the 4-cell stage, which progressively diminished by the 8-cell and morula stages (Supplementary Fig. S1 A). Overall developmental trajectories were broadly similar between humans and pigs, showing relatively limited transcriptional changes prior to the 8-cell stage. Notably, in pigs, transcriptional differences between the 2-cell and 4-cell stages were more pronounced, resulting in a developmental trajectory that more closely aligned the 4-cell stage with the 8-cell stage compared with that observed in humans (Supplementary Fig. S1 A). Although substantial transcriptional divergence was evident between the 4-cell and 8-cell stages in both species, these differences became less pronounced during the transition from the 8-cell stage to the morula (Supplementary Fig. S1 A). Cross-species data integration for developmental stage alignment revealed that the earliest and latest stages exhibited the highest degree of correspondence between humans and pigs (Fig. 1 C). As development progressed, the transcriptional distance between the human 2-cell and 4-cell stages increased, reflecting a pattern similar to that observed in pigs. Notably, the human 4-cell stage emerged as a distinct transcriptional state along the pseudotime trajectory. Interspecies divergence became increasingly apparent as embryos advanced from the oocyte to the 4-cell stage. During the transition from the human 4-cell to 8-cell (H8) stage, embryos displayed transcriptomic similarity to both the porcine 4-cell (P4) and 8-cell (P8) stages, indicating partial developmental overlap between species (Fig. 1 C). To account for intragroup heterogeneity and improve stage specificity, we refined the dataset by excluding samples with low internal consistency or ambiguous stage assignment based on pseudotime positioning and correlation analyses. This filtering reduced the total number of samples to 116 for subsequent analyses (Fig. 1 Cand Supplementary Fig. S2 A–C). We further characterized embryonic development by examining both gene expression levels and gene amount, defined here as the number of actively transcribed genes at each developmental stage. Across both species, these metrics exhibited broadly similar developmental trends. From the oocyte to the morula stage, gene amount initially increased, followed by a marked decline and stabilization at relatively low levels, accompanied by increased intragroup variability. Notably, this downregulation began at the 4-cell stage in human embryos but occurred earlier, at the 2-cell stage, in porcine embryos. In pigs, gene amount uniquely exhibited a rebound at the morula stage (Fig. 1 D). Overall gene expression levels decreased from the oocyte to 2-cell stage, then rose sharply to 4-cell stage, surpassing oocyte levels, before declining again and showing a modest recovery between the 8-cell and morula stages (Fig. 1 E). In human embryos, the proportion of mitochondrial-encoded (MT) gene expression was initially lower than in pigs but increased significantly from the oocyte stage onward, ultimately exceeding porcine levels and displaying greater intragroup variability. In contrast, pigs exhibited more moderate fluctuations across developmental stages and reduced variability within groups (Fig. 1 F). Together, these findings indicate that as embryos progress through the ZGA process, the transcriptional programs of human and porcine embryos become increasingly coordinated, suggesting a convergence of developmental trajectories during this critical window of early embryogenesis. To delineate developmental stages through key TFs, and conversely to assess how species-specific and stage-specific programs shape TF dynamics, we examined the expression patterns of POU5F1 , SOX2 , KLF4 , NANOG , and MYC . In both species, POU5F1 exhibited relatively strong maternal deposition, with KLF4 showing notable maternal contribution only in pigs. The remaining four of the five TFs displayed marked upregulation during major ZGA (Fig. 1 G). Among these factors, MYC demonstrated the most pronounced induction and achieved the highest expression levels during ZGA. Both MYC and NANOG were present at low maternal mRNA levels but showed sharp increases during ZGA, peaking at the 4-cell stage in pigs and at the 8-cell stage in humans, followed by a decline toward the morula stage (Fig. 1 G). KLF4 and SOX2 were also significantly upregulated during major ZGA; notably, their expression remained elevated in pigs with little to no subsequent downregulation, whereas only modest decreases were observed overall (Fig. 1 G). In contrast, POU5F1 displayed distinct species-specific dynamics: in humans, it originated from a relatively low maternal reserve and increased steadily from the 4-cell stage onward, while in pigs, it was maintained at consistently high levels throughout early development (Fig. 1 G). We investigated cross-species correspondence between developmental stages using stage-specific marker genes, with the top five markers for each stage selected independently. Clustering based on the expression matrices of these markers revealed interspecies similarities at two paired developmental intervals: the 1-cell and 2-cell stages, and the 8-cell and morula stages. In contrast, the 1-cell stage in both pigs and humans functioned as a transitional state bridging these paired stages (Fig. 1 H). Conversely, clustering analyses based on downregulated marker genes primarily separated samples by species rather than by shared developmental stage. Within each species, developmental stages formed similar hierarchical relationships, with the 1-cell, 2-cell, 4-cell stages clustering together and the 8-cell and morula stages forming a second cluster. Notably, the 1-cell and 2-cell stages exhibited the highest degree of similarity within their respective clusters (Supplementary Fig. S1 C). Downregulated marker genes exhibited strong stage-wise similarity within each species. In humans, PRAMEF15 was specifically enriched at the 1-cell and 2-cell stages, whereas SULT1E1 characterized later developmental stages, including the 4-cell, 8-cell, and morula stages (Supplementary Fig. S1 B). A comparable pattern was observed in pigs, where RPS13 and TPT1 served as markers for the 1-cell and 2-cell stages, while RPL13A , RPL19 , and H3-3B were associated with later stages of development (Supplementary Fig. S1 B). Notably, several ribosomal protein genes ( RPS13/23 and RPL12/13A/19/26/36 ), along with H3-3B , which is involved in histone modification, were consistently and highly expressed across all developmental stages in human embryos (Supplementary Fig. S1 C). In contrast, PLIN2 and PLA2G15 , genes associated with lipid metabolism, showed elevated expression specifically at the 8-cell and morula stages in pig embryos (Supplementary Fig. S1 C). These divergent expression patterns likely reflect species-specific evolutionary adaptations. In humans, sustained ribosomal activity and histone modification may support complex organogenesis and enhanced epigenetic plasticity, features characteristic of primate development. Conversely, in pigs, increased lipid metabolic activity may correspond to dietary flexibility and the heightened energy storage requirements necessary to sustain rapid post-embryonic growth. Although no overlap was observed among upregulated marker genes across developmental stages or between species, clear functional parallels were evident at corresponding stages. During early development, TSPAN19 , associated with signal transduction, marked the human 1-cell stage, whereas CERK , involved in lipid-mediated signaling, characterized the porcine 2-cell stage (Fig. 1 H). Similarly, at the 1-cell stage, IFNE in humans and IL1R2 in pigs were both implicated in immune regulation and inflammatory responses (Fig. 1 H). As development progressed, FBP2 in humans and PCK1 in pigs, marking the 4-cell and 8-cell stages, respectively, were associated with carbohydrate metabolism (Fig. 1 H). At later stages, PLXNB3 expression at the human 8-cell stage and FLRT2 expression at the porcine 4-cell stage highlighted conserved roles in neural development and cell adhesion (Fig. 1 H). At the morula stage, CCL15 in humans and PTGS2 and IL6 in pigs were enriched and are known to play key roles in immune modulation and inflammatory signaling (Fig. 1 H). The preservation of functional themes despite divergence at the individual gene level suggests that humans and pigs have evolved distinct gene regulatory programs to accomplish similar developmental objectives. These differences may reflect species-specific life history traits, such as the greater emphasis on neural complexity in humans versus metabolic efficiency and rapid growth in pigs. Exploration of stage relationships between species To assess interspecies similarities in developmental progression, we performed differential expression analyses across consecutive stages within each species. Differentially expressed genes (DEGs) were identified using a threshold of |log2 fold change (FC)| > 1 and an adjusted p value < 0.05, whereas genes with |log2FC| < 0.25 were classified as developmentally conserved. Both species exhibited parallel trends in dynamic transcriptional reprogramming, defined here as similar directional changes (e.g., increases, peaks, and declines) in the number of DEGs at corresponding developmental milestones, despite differences in absolute developmental timing. Specifically, humans and pigs displayed comparable patterns in the numbers of stage-specific upregulated and downregulated DEGs, reflecting shared regulatory transitions between consecutive developmental stages. In addition, a substantial proportion of developmentally conserved genes was observed in both species, indicating the maintenance of stable gene expression programs across successive stages within each developmental trajectory. For upregulated DEGs, the number in human embryos declined sharply from approximately 2179 genes at the 1-cell versus oocyte comparison to 253 genes at the 2-cell versus 1-cell stage. This was followed by a moderate increase to 1005 genes at the 4-cell versus 2-cell stage and a pronounced surge to 3115 genes at the 8-cell versus 4-cell stage, as identified using DESeq2 with false discovery rate (FDR) correction (adjusted p value < 0.05). Porcine embryos exhibited a comparable trend, with substantial upregulation observed at both major ZGA-associated stages, reaching levels similar to those in humans, followed by a decline at the morula stage in both species (Fig. 2 A). In contrast, downregulated DEGs in both humans and pigs increased from the 2-cell versus 1-cell comparison and peaked at the 8-cell versus 4-cell stage, with a greater magnitude of downregulation observed in humans, before decreasing at the morula stage (Fig. 2 A). Notably, maternally derived transcripts (TPM ≥ 10) constituted the majority of downregulated genes at each developmental stage (Supplementary Fig. S3 A). The parallel fluctuations in DEG numbers and expression patterns across developmental transitions underscore a conserved framework of transcriptional reprogramming between human and porcine embryos, particularly during the ZGA phase, when coordinated activation of the embryonic genome occurs in both species. Fig. 2. Open in a new tab Differential analysis and enrichment of DEGs across consecutive developmental stages from oocyte to morula for human and pig. ( A ) a line chart shows conserved, upregulated, and downregulated gene counts from differential analysis between consecutive stages. ( B ) paired violin plots with overlaid simplified boxplots display the distribution of log₂ FC for DEGs between consecutive developmental stages (e.g., 1-cell vs. oocyte). The central dot and vertical line within each half-violin represent the median and IQR, respectively. Individual points beyond a log₂ FC of 15 are plotted as scatter points to highlight extreme outliers. ( C ) paired violin plots with overlaid simplified boxplots (as in B) show the expression level distribution (log₂[TPM+1]) of DEGs across adjacent stages. For each stage, the paired distributions are presented side-by-side with a gap. ( D ) an alluvial plot shows intersection sizes of DEGs for consecutive stages between species, distinguishing between upregulated (left) and downregulated (right) DEGs. ( E ) chord diagrams depict the intraspecies relationships of DEGs across adjacent stages, with pigs on the left and humans on the right. ( F ) Heatmap and word clouds present the clustering results of GO terms of U1. ( G ) bubble charts detail the participation of DEGs in KEGG pathways across adjacent developmental stages Notably, during the 1-cell versus oocyte transition, human embryos exhibited a substantially greater number of both upregulated and downregulated genes than porcine embryos (Fig. 2 A). However, the corresponding FCs and absolute expression levels were generally lower in humans than in pigs and did not approach the magnitude observed during human major ZGA (Fig. 2 B, C). In contrast, the number of developmentally conserved genes tended to show an inverse relationship with the extent of upregulation, particularly during the 2-cell versus 1-cell transition, where conserved gene counts were higher in humans (Fig. 2 A). Overall, although the human 4-cell versus 2-cell comparison displayed relatively limited transcriptional changes, it represented a marked reversal relative to the preceding 2-cell versus 1-cell transition, indicating a shift in the underlying regulatory dynamics. Differential expression FCs and overall expression patterns showed closer alignment between humans and pigs during the eight-cell versus four-cell transition, which corresponds to the peak transcriptional activation in humans. In contrast, pigs reached their peak during the four-cell versus two-cell transition (Supplementary Fig. S2 B, C). Consistent with this pattern, pronounced FC outliers emerged during this interval in humans, whereas in pigs such outliers appeared earlier, beginning at the 4-cell versus 2-cell stage, and were both more frequent and more extreme (Supplementary Fig. S2 B). In both species, genes with exceptionally high expression levels (loge[TPM + 1] > 9) were predominantly MT genes, with a higher prevalence observed in humans than in pigs. This pattern was particularly evident in humans at the 1-cell versus oocyte, 4-cell versus 2-cell, and morula versus 8-cell comparisons. Notably, 11 of the 13 MT genes were upregulated in humans at the 1-cell versus oocyte stage, whereas minimal overlap in MT DEGs was observed at the 4-cell and morula stages (Supplementary Fig. S3 B). Among non- MT outliers, FTL , identified as a DEG outlier at the 1-cell versus oocyte stage in pigs, exhibited high maternal transcript abundance in both species but followed opposite expression trajectories. Specifically, FTL was generally upregulated in pigs, whereas in humans it was rapidly downregulated after the 2-cell stage (Supplementary Fig. S3 C). Similarly, CKB , a DEG outlier at the morula versus 8-cell stage in pigs, showed a species-specific expression spike during the transition from the 8-cell stage to the morula. We further examined FTH1 , a closely related gene, which displayed substantial intrastage heterogeneity and consistently higher expression levels in pigs compared with humans. FTH1 followed a continuous upregulation trend in pigs, similar to FTL , whereas its upregulation in humans was comparatively modest. Importantly, both FTL and FTH1 reached extremely high expression levels in pigs (approximately log[TPM + 1] = 10) during major ZGA, FTL peaking at the 4-cell versus 2-cell transition and FTH1 at the 8-cell versus 1-cell transition, thereby markedly amplifying interspecies transcriptional differences during this critical developmental window (Supplementary Fig. S3 D). We further examined the expression trajectories of additional genes associated with ferroptosis (iron-dependent cell death) and observed marked differences between minor and major ZGA phases. LTF , although expressed at relatively low levels, is involved in iron sequestration and was detected exclusively in porcine embryos, where it underwent pronounced downregulation during major ZGA. Iron transport–related genes, including TF and TFR2 , exhibited higher expression levels in human embryos, whereas TFRC was more strongly expressed in pigs. Genes involved in redox homeostasis and ferroptosis regulation, such as GSS and GSR , showed more pronounced expression in pigs, suggesting enhanced redox activity during early development. In contrast, the core ferroptosis-regulating enzyme GPX4 displayed highly similar expression trajectories between humans and pigs, indicating conservation of this critical regulatory component across species (Supplementary Fig. S4 ). The pronounced differences between pigs and humans at the 1-cell versus oocyte stage warranted closer examination. In porcine embryos, the 1-cell stage was uniquely enriched for pathways related to lipid metabolism and ferroptosis, whereas human embryos exhibited enrichment primarily in core MT respiratory chain processes (Supplementary Fig. S5 A). Genes involved in lipid metabolic pathways generally displayed higher expression levels in pigs than in humans, with particularly notable differences observed for HACD4, HADH, HADHA, HADHB, ACSL1 , and HTD2 . Among these, HACD4 showed especially low expression in human embryos compared with pigs (Supplementary Fig. S5 B). These divergent expression patterns likely reflect species-specific evolutionary adaptations. In humans, elevated expression of iron-handling genes (e.g., TF) and enrichment of MT respiratory activities may correspond to the high metabolic demands associated with neurodevelopment and cognitive evolution. In contrast, in pigs, the upregulation of iron-storage genes (e.g., FTL ) and lipid metabolic pathways may facilitate rapid growth and dietary flexibility, traits that favor energy efficiency in porcine evolution (Supplementary Fig. S4 ). Beyond iron metabolism, broader transcriptional differences further underscore fundamental evolutionary divergence between the two species. Comparative analysis revealed significant correspondence of DEGs between species during major ZGA-associated transitions, specifically the human 8-cell versus 4-cell stage and the porcine 4-cell versus 2-cell stage (Fig. 2 D). The intersecting set of upregulated DEGs (designated U1) was more strongly enriched during the 4-cell versus 2-cell transition than during the 8-cell versus 4-cell transition, comprising 809 and 1242 genes, respectively, with a shared core of 300 genes (Supplementary Fig. S3 E). In pigs, transcriptomic correlations were strongest between the 8-cell versus 4-cell stage and its adjacent stages, with particularly tight coupling observed between major ZGA transitions. In contrast, in humans, comparable transcriptomic similarity was primarily evident between the 8-cell versus 4-cell stage and the morula stage (Fig. 2 E). Notably, the correlation between the human 1-cell versus oocyte transition and the human 8-cell versus 4-cell transition was comparable in magnitude to that observed among porcine major ZGA stages (Fig. 2 E). Furthermore, this early human transition also exhibited relatively strong correlations with major ZGA stages in pigs (Fig. 2 D), suggesting partial conservation of transcriptional programs across species despite differences in developmental timing. Ridge plot analysis revealed that U1 genes exhibited higher expression levels compared with other ZGA-associated genes, with more than 40% of U1 genes exceeding a loge (TPM + 1) value of 4 (Supplementary Fig. S3 F). Gene Ontology (GO) term clustering further demonstrated that U1 genes were significantly enriched in biological processes related to RNA splicing and ribosome biogenesis (Fig. 2 F and Supplementary Fig. S3 G). Enrichment significance was assessed using Benjamini–Hochberg FDR correction, with adjusted p values < 0.05 considered statistically significant. Consistent with these findings, Kyoto encyclopedia of genes and genomes (KEGG) pathway analysis indicated that U1 genes were actively involved in autophagy, mRNA surveillance, protein turnover, and cellular energy metabolism. Among these pathways, ribosome biogenesis emerged as the most prominently enriched pathway (Supplementary Fig. S3 H), with all reported pathways meeting the FDR significance threshold. In addition to analyzing shared ZGA-associated genes, we performed GO enrichment and clustering analyses on genes specific to the P4, P8, and H8 stages. The results showed that P4 and H8 were enriched in multiple categories of fundamental biological processes and exhibited substantial functional overlap, including pathways related to histone acetylation and chromatin modification. In contrast, P8 was associated with a more specialized functional profile, distinct from the other two stages, and was predominantly enriched in metabolic processes (Supplementary Fig. S6 A–C). By stratifying GO enrichment results using the keywords “metabolic” and “catabolic”, we identified 131, 257, and 499 metabolism-related genes and 17, 23, and 81 corresponding GO terms for H8, P4, and P8, respectively (Supplementary Fig. S6 D). Consistent with this finding, combined GO and KEGG analyses revealed that P8 exhibited a large number of metabolism-enriched pathways, encompassing lipid, amino acid, carbohydrate, glycan, vitamin, and coenzyme metabolism (Supplementary Fig. S6 F and Supplementary Fig. S7 B). In contrast, metabolism-related terms for P4 and H8 were largely similar and primarily associated with ATP synthesis and energy production; most of these pathways were also included among the enriched terms identified for P8 (Supplementary Fig. S7 A). When GO enrichment results were categorized using the keyword “histone”, H8, P4, and P8 were associated with 43, 69, and 7 genes, and 10, 18, and 1 corresponding GO terms, respectively (Supplementary Fig. S6 E). Chromatin remodeling processes were predominantly related to histone acetylation. While the enrichment profiles of P4 closely resembled those of H8, P4 uniquely exhibited enrichment in pathways associated with histone ubiquitylation (Supplementary Fig. S6 G). Network analysis revealed that chromatin remodeling in P4 was centered around KAT5 , EP300 , and BRCA1 , whereas in H8, the most highly connected genes included RNF168, SMARCA5 , and ATM . Analysis of transcriptional regulators further indicated that GATA3, THAP7, MIDEAS , and ARID5A were prominent TFs associated with P4-specific genes. In contrast, the P8-specific TF TP53 showed strong connectivity with chromatin remodeling–associated genes and appeared to regulate genes involved in metabolic processes or genes located in non-overlapping genomic regions (Supplementary Fig. S7 C–G). We defined chromatin segments containing ZGA-associated genes as transcriptionally active regions and examined their genomic distribution and coverage across developmental stages. From minor ZGA at the 1-cell versus oocyte stage to major ZGA at the porcine 4-cell (P4) and 8-cell (P8) stages, both the distribution density and chromosomal coverage of ZGA genes increased markedly (Supplementary Fig. S8 A–D). Although P4 and P8 exhibited comparable numbers of ZGA genes per chromosome, their chromosomal coverage patterns differed substantially. For instance, P4 showed markedly higher enrichment on chromosomes 1, 2, 3, 5, 7, 9, 11, 12, and 13, whereas P8 displayed stronger enrichment on chromosomes 8 and 10 (Supplementary Fig. S8 E). To further quantify chromosomal similarity between P4 and P8, each chromosome was partitioned into equal-sized bins, and the number of ZGA genes within each bin was used as a measure of segmental activity. Correlation analyses based on these activity profiles revealed generally low concordance between P4 and P8. Even for chromosomes exhibiting relatively high ZGA gene coverage, such as chromosomes 8, 9, 11, 12, 15, 16, and 18, correlation coefficients did not exceed 0.5, with chromosome 12 reaching only approximately 0.6 (Supplementary Fig. S8 F). Consistently, analysis using 1 − Jaccard coefficients (250-bin resolution) demonstrated frequent and pronounced fluctuations across chromosomes 1 through X, with sustained high values in specific genomic regions, indicating substantial segment-level divergence (Supplementary Fig. S8 G). We next examined the genomic distribution of chromatin remodeling–related genes at the P4 stage and observed that highly connected regulators, including EP300/EP400 and KAT5 , exhibited porcine-specific, contiguous enrichment across two major ZGA-associated chromatin segments accompanied by elevated expression levels. Using GTRD annotations, we further found that the target genes bound by these chromatin remodeling factors were significantly and non-randomly distributed near transcription start sites (TSSs), with closely aligned patterns observed at both the P4 and P8 stages (Supplementary Fig. S8 H, I). Upset plot analysis identified a total of 1040 intersecting DEGs shared between the human 1-cell versus oocyte transition and the porcine major ZGA stages. These genes were evenly distributed across both porcine major ZGA stages, among which 304 genes, designated U2, were uniquely upregulated during human major ZGA as well (Supplementary Fig. S9 A, B). Functional enrichment analysis revealed that U2 genes were predominantly associated with processes involving mRNA metabolism, MT function, and ribonucleoprotein complexes, with strong enrichment in RNA processing and MT protein transport pathways (Supplementary Fig. S9C–E). Within this group, key RNA-processing genes such as DDX56 and SRSF1 exhibited a distinctive two-stage upregulation pattern that was unique to human embryos (Supplementary Fig. S10 A). Among MT-related genes, TOMM40 displayed low baseline expression in both species but showed minimal expression in pigs until the onset of major ZGA (Supplementary Fig. S10 B). Genes directly involved in protein translocation into mitochondria, including TIMM50 and TOMM40 , exhibited broadly similar activation profiles; however, TIMM50 followed a conserved expression trajectory between species, whereas TOMM40 was not substantially expressed in pigs until the major ZGA period (Supplementary Fig. S10 B). Notably, two additional translocation-associated genes, TIMM8A and TIMM10 , displayed expression patterns in pigs similar to that of TOMM40 (Supplementary Fig. S10 C). In contrast to the functional emphasis of U2 genes, the set of genes commonly upregulated during both porcine major ZGA stages was predominantly enriched in pathways related to ATP synthesis, highlighting a species-specific focus on energy production during porcine embryonic genome activation (Supplementary Fig. S10 ). To further investigate functional correspondence between developmental stages across species, we analyzed KEGG pathway regulation patterns driven by DEGs. During major ZGA stages in both humans and pigs, we observed convergent regulation of pathways involved in fundamental cellular development and maintenance, including autophagy, ribosome biogenesis, and RNA surveillance (Fig. 2 G; KEGG database accessed on 6 August 2023; FDR < 0.05). Energy metabolism–related KEGG pathways became prominently regulated beginning at the 4-cell versus 2-cell transition and showed increasing numbers of upregulated terms, with enhanced cross-species similarity observed at the morula stage (Fig. 2 G). In contrast, pathways associated with cell–cell interaction and signal transduction, such as focal adhesion and MAPK signaling, exhibited divergent regulation patterns: these pathways were generally downregulated during human major ZGA but were upregulated in pigs at the adjacent morula stage (Fig. 2 G). Notably, the p53 and HIF-1 signaling pathways were significantly upregulated in pigs during this period. Despite these differences, both species demonstrated significant downregulation of MAPK and Ras signaling pathways at the eight-cell versus four-cell transition (Fig. 2 G), indicating a shared regulatory shift during this critical phase of embryonic genome activation. To identify evolutionarily divergent signaling pathways, we prioritized genes exhibiting statistically significant reversal expression patterns between species (FDR < 0.05) within highly interconnected biological networks. Several core signaling components, including MAPK1, MAPK9, MAPK12, AKT1, AKT2, PIK3R2, PRKCB, PLCG1, and PLCG2 , showed strong associations with these reversal patterns, characterized by downregulation in human embryos and upregulation in pig embryos around the morula stage (Supplementary Fig. S11 A, B). Among these genes, AKT1 displayed the highest network degree, indicating a central regulatory role (Supplementary Fig. S11 C). In addition, PIK3CB, HRAS, and RAC3 exhibited markedly higher upregulation in pigs than in humans, further reinforcing species-specific regulatory trends (Supplementary Fig. S11 A). Collectively, MAPK1, AKT1/2, PIK3R2, and PRKCB were found to influence a broad spectrum of downstream pathways. Conversely, pathway-level analysis revealed that the most affected signaling cascades included Ras signaling, neurotrophin signaling, Rap signaling, and growth hormone–related pathways (Supplementary Fig. S11 D), highlighting key axes of evolutionary divergence in post-ZGA developmental regulation. Collectively, the enrichment of RNA processing and MT functions in humans, together with the pronounced activation of metabolic pathways and growth-related signaling in pigs, likely reflects divergent life-history strategies: humans favoring neurological complexity and energy-intensive cognitive development, whereas pigs prioritize rapid growth and metabolic adaptability in resource-variable environments. Exploring the correlation in maternal mRNA for chromatin remodeling across species During the early MZT, embryonic development relies predominantly on maternally deposited mRNAs. Comparative transcriptomic analysis revealed both shared features and distinct differences between human and pig oocytes (Fig. 1 C). A TPM threshold of 2.5 was applied to distinguish abundantly expressed maternal transcripts from those with minimal expression (Supplementary Fig. S12 A), revealing that pigs possess approximately 10% fewer maternal mRNAs than humans (Supplementary Fig. S12 A). Differential transcriptome analysis between the two oocyte types enabled the classification of maternal mRNAs based on shared and species-specific expression patterns, identifying 5379 genes commonly expressed in both species (Fig. 3 A). Among these shared genes, 1539 exhibited higher expression in humans, whereas 1218 were more highly expressed in pigs (Fig. 3 A). In addition, 410 common genes showed no significant interspecies difference but displayed generally lower expression levels. Beyond the shared maternal transcript pool, 2524 and 1358 genes were predominantly expressed in human and pig oocytes, respectively (Fig. 3 A), with a subset exhibiting high absolute expression levels (Fig. 3 Band Supplementary Fig. S12 B). Fig. 3. Open in a new tab Differential and enrichment analysis of maternal mRNA in both human and pig. ( A ) a venn diagram illustrates the interrelationships among maternal mRNAs in humans and pigs, alongside DEGs that are upregulated or downregulated with pigs compared to humans. ( B ) maternal mRNAs are categorized into four sets based on their commonality and significance of difference between species. Fully standardized heatmap displays the expression profiles of these four sets as well as non-maternal mRNAs in oocytes of pigs and humans. ( C ) a heatmap and word clouds present the clustering results of GO terms for maternal mRNAs shared by both humans and pigs. Adjacent to it, a bar chart illustrates the number of GO terms and the proportion of involved genes for each corresponding cluster. ( D ) Bar charts highlight the top five significant items related to histone modification, mRNA, and DNA methylation within cluster d in panel C. ( E ) Bar charts show enrichment results for specific histone modifications within cluster d in panel C. ( F ) Bar charts detail the representative KEGG enrichment results for maternal mRNAs common to both humans and pigs. ( G ) protein-protein interaction (PPI) network analysis of the ATP-dependent chromatin remodeling pathway. The network displays experimentally validated interactions of core genes in the ATP-dependent chromatin remodeling pathway from the KEGG database. Node size reflects protein connectivity (degree centrality), with star-shaped nodes representing subunits of distinct chromatin remodeling complexes (color-coded by complex type), and yellow nodes indicating DNA-specific recognition TFs. Edge width indicates interaction confidence, based on experimental evidence from the STRING database. ( H ) PPI network of the PRC pathway. The network shows experimentally validated interactions among core members of the Polycomb protein family. Edge color depth denotes interaction types (blue: physical interaction; red: functional association) Despite containing numerous DEGs, the common gene subset exhibited the highest average expression levels compared with the other subgroups (Supplementary Fig. S12 B). Functional enrichment and clustering analyses indicated that these genes are primarily involved in the cell cycle, regulation of gene expression, organelle localization, stress responses, signal transduction, protein folding, and maintenance of cellular stemness (Fig. 3 C). Among the functional categories, Category D, encompassing GO terms related to mRNA processing and histone biology, contained the largest number of enriched GO terms and involved the highest proportion of genes (Fig. 3 C). Notably, 42 histone-associated GO terms were identified, with the top five predominantly related to histone H4 acetylation and other specific post-translational modifications, highlighting coordinated regulation through both acetylation and methylation mechanisms (Fig. 3 D). In parallel, 40 GO terms associated with mRNA regulation were enriched, with the top five terms emphasizing RNA splicing and mRNA metabolic processes (Fig. 3 D). Histone modifications represent pivotal events in chromatin remodeling during early embryogenesis and play a fundamental role in regulating subsequent developmental processes. Further extraction of histone modification–related GO terms from Category D revealed modifications involving histones H2A, H2B, H3, and H4, with the highest gene participation observed in H3 acetylation, H3K4 methylation, and H4 acetylation (Fig. 3 E). Among these, H3-associated modifications displayed the greatest diversity of enriched GO terms, including H3K9 modification, H3K14 acetylation, and methylation at H3K4 and H3K36 (Fig. 3 E). Modifications of H2 histones were predominantly characterized by ubiquitination of H2A and H2B, along with acetylation of H2A (Fig. 3 E). In contrast, H4-related modifications were mainly defined by acetylation at H4K16, H4K8, and H4K5, as well as methylation at H4K20 (Fig. 3 E). Notably, H3 acetylation exhibited gene involvement counts and enrichment significance comparable to those observed for H4 acetylation, underscoring their coordinated roles in early embryonic chromatin reorganization (Fig. 3 E). We assessed the conservation of histone modification–related regulatory mechanisms between humans and pigs by identifying core genes based on multiple criteria, including GO term involvement, expression abundance, and protein–protein interaction (PPI) network connectivity (Supplementary Fig. S12 C). Applying a GO term participation threshold of ≥15, we identified BRCA1, WDR5, KMT2A, KAT8, and SMARCB1 as key candidate genes. With the exception of BRCA1, which exhibited slightly higher expression in humans, the remaining four genes showed no statistically significant interspecies differences. Additional DEGs that were more highly expressed in humans also displayed moderate to high expression levels in pigs, indicating overall conservation of expression magnitude (Supplementary Fig. S12 B, F). Analysis of expression rankings revealed that the top 10 most highly expressed genes in both species included DNMT1, KPNA7 (the highest expressed gene in both), and RBBP7 (ranked among the top four in each species), none of which demonstrated significant interspecies differences (Supplementary Fig. S12 G). Likewise, genes exhibiting relatively lower expression in one species generally maintained moderate expression levels in the other (Supplementary Fig. S12 B, G). Using PPI network degree as an additional criterion, with a threshold of ≥70, we identified WDR5, HDAC1, EP300, KAT2A, and HDAC2 as highly connected hub genes. Although the latter four showed statistically significant interspecies differences, their connectivity values remained close to the threshold, suggesting functional conservation (Supplementary Fig. S12 D, E, H). Notably, the essential chromatin regulators HDAC1 and HDAC2 were specifically expressed at higher levels in pig oocytes (Supplementary Fig. S12 G), highlighting subtle but potentially meaningful species-specific regulatory biases. KEGG pathway enrichment analysis of common maternal mRNAs, performed using a hypergeometric test with Benjamini–Hochberg FDR correction (adjusted p -value < 0.05), revealed significant enrichment in pathways associated with RNA processing (e.g., RNA splicing), chromatin remodeling (e.g., Polycomb repressive complex (PRC) and ATP-dependent chromatin remodeling), and energy metabolism (e.g., oxidative phosphorylation, amino acid degradation, and fatty acid metabolism) (Fig. 3 F). Among these, pathways related to RNA regulation and chromatin remodeling were the most prominently enriched, with thermogenesis and the spliceosome encompassing the largest number of involved genes (Fig. 3 F). The PRC and ATP-dependent chromatin remodeling play central roles in chromatin reorganization during early embryonic development. Core genes associated with these processes were generally conserved between humans and pigs (Supplementary Fig. S13 A, B). Within ATP-dependent chromatin remodeling, both species recruited complexes from the ISWI, CHD, SWI/SNF, and INO80/SWR families (Supplementary Fig. S14 A). Notably, BRG1/BRM complexes were broadly represented in maternal mRNAs from both species, although SMARCA2 was detected only in humans (Supplementary Fig. S14 A). Key components such as SMARCA5 (ISWI family) and ACTL6A (BRG1/BRM complex) were present in both species and exhibited largely comparable expression patterns, with SMARCA5 showing relatively higher expression in pigs (Fig. 3 G; Supplementary Fig. S13 A). In contrast, human embryos displayed significantly higher maternal expression of BRG1/BRM core ATPases, including SMARCA4 and SMARCA2 , compared with pigs. SMARCA2 expression subsequently declined in humans, and expression levels of these genes converged between species during the major ZGA stages (Supplementary Fig. S13 C). In addition, ARID1B , a potential regulatory subunit of the SWI/SNF complex, was expressed below the maternal threshold specifically in pigs (Supplementary Fig. S13 D). Core components of the PRC, encompassing both PRC1 and PRC2 subtypes, appear to be functionally conserved between humans and pigs but display pronounced species-specific expression dynamics. For instance, the core PRC2 catalytic subunits EZH1 and EZH2 are conserved in both species, whereas RNF2 shows divergent regulation: it is markedly upregulated in pigs during early development but declines after the major ZGA stages, in contrast to humans, where RNF2 expression rises sharply at the 8-cell stage (Supplementary Fig. S14 B; Fig. 3 H; Supplementary Fig. S13 B, C). Another PRC1 core component, RING1, exhibited higher enrichment in humans and showed an expression pattern opposite to that of RNF2 (Supplementary Fig. S13 C). Further analysis of CBX family members within the PRC1 complex, which recognize the H3K27me3 mark, also revealed reciprocal expression patterns between species (Supplementary Fig. S14 B). In humans, CBX2 initially exhibited the highest expression, whereas CBX8 predominated in pigs. In human embryos, CBX8 expression increased after zygote formation but declined by the 8-cell stage, while in pigs it was strongly upregulated during major ZGA. Conversely, CBX2 expression in pigs did not increase until the major ZGA stages, whereas in humans it decreased at both the 1-cell and 8-cell stages. CBX4 displayed a broadly conserved trajectory in both species, with expression levels near the maternal mRNA threshold initially and pronounced upregulation during major ZGA, particularly in pigs (Supplementary Fig. S13 E). In contrast, CBX6 showed a species-specific expression pattern, being selectively induced during human major ZGA. Gene set variation analysis (GSVA) in both species confirmed that histone acetylation is generally more active than methylation during early embryonic development, remaining highly activated until the morula stage, with a notable decline in H4 acetylation observed in pig morulae. Overall, pigs exhibited higher levels of H3 acetylation than humans, whereas humans showed higher levels of H4 and H2A acetylation as well as H2A ubiquitination (Supplementary Fig. S13 F). In contrast, H3K9 methylation was markedly downregulated to low levels during the major ZGA period. mRNA decay drives the inter-species transcriptomes more similarly We observed that interspecies transcriptional divergence increased from the oocyte to the 4-cell stage but decreased markedly from the 4-cell to the morula stage. Based on this pattern, we performed interspecies difference analyses using a receiver operating characteristic (ROC)–based method for the 4-cell and 8-cell stages separately. Genes whose interspecies difference scores were reduced to less than half of their 4-cell values during the transition to the 8-cell stage were defined as difference-reducing genes, whereas genes whose difference scores increased more than twofold were classified as difference-expanding genes. The relationships between these gene categories and distinct subsets of maternal mRNAs, as well as major ZGA-associated genes at the 8-cell stage, are illustrated using flow diagrams (Fig. 4 A, B). Difference-reducing genes were predominantly enriched among non-shared maternal mRNAs exhibiting significant interspecies differences and were far more abundant than genes in other categories. In contrast, difference-expanding genes were mainly concentrated in the pig 8-cell (P8) stage. Consistent with the overall trend toward transcriptional convergence, the number of difference-reducing genes substantially exceeded that of difference-expanding genes (Fig. 4 A, B). Notably, the majority of difference-reducing genes corresponded to genes that were significantly downregulated at the 8-cell stage (Fig. 4 C). Functional enrichment analyses of these downregulated genes using GO and KEGG revealed predominant involvement in cell cycle regulation, protein degradation, autophagy, cell division, RNA processing, and related biological processes (Fig. 4 D; Supplementary Fig. S15 A–D). Fig. 4. Open in a new tab Analysis of genes responsible for the reduction and expansion of interspecific variation in 4- to 8-cells. ( A ) Sankey diagram showing distribution of genes with reduced (DOWN) versus expanded (UP) interspecific differences. Using a simplified layout with color gradients (blue: reduction, red: expansion) to visualize gene allocation across subsets of maternal mRNAs and major zygotic genome activation ZGA. ( B ) heatmap showing the expression of genes with reduced (DOWN) versus expanded (UP) interspecific differences in 4-cell and 8-cell embryos. ( C ) distribution of genes with reduced interspecies differences among significantly down-regulated genes at each time period. ( D ) bubble plot of GO enrichment for down-regulated genes co-enriched with pig maternal environment. Gene counts per GO term are added (bubble size), and enrichment significance is indicated by color intensity. ( E ) GO enrichment results for interspecific difference-expanding genes (enriched with the major ZGA genes of pigs) Vital TFs for major ZGA event for both human and pig TFs play a central role in regulating ZGA, likely orchestrating the conserved epigenetic activation programs, such as chromatin remodeling and histone modification, identified in the preceding analyses. Examination of the top 20 most highly expressed TFs during the major ZGA stages in both species revealed substantial overlap (Supplementary Fig. S16 A). Specifically, nine TFs were shared between the human 8-cell stage and the pig 4-cell stage: MYC, FOXR1, ATF4, LEUTX, KLF17, ZNF207, LYAR, GPBP1, and TFAM (Supplementary Fig. S16 A, B). Moreover, five of these TFs, ATF4, KLF17, ZNF207, LYAR, and TFAM , were common to the human major ZGA stage and both major ZGA stages in pigs (Supplementary Fig. S16 A, B). Notably, MYC ranked as the most highly expressed TF in humans, with FOXR1 and ATF4 among the top three. In contrast, LYAR consistently ranked second during both pig major ZGA stages, while KLF17 ranked first at the pig 4-cell stage and fifth at the pig 8-cell stage (Supplementary Fig. S16 B). We assessed TF activity across developmental stages in both species using SCENIC (v1.3.1), with the AUCell algorithm applied to score regulon activity, followed by soft clustering and pattern classification. In pigs, 61 and 45 regulons exhibited stage-specific activation at the 4-cell and 8-cell stages, respectively, before being subsequently downregulated. In humans, 33 regulons showed pronounced activation at the 8-cell stage (Fig. 5 A). In addition, 44 regulons in pigs were upregulated and maintained high activity from the 4-cell stage through the morula stage, while humans displayed 74 regulons with a comparable sustained activation pattern (Fig. 5 A). The intersecting set of regulons between humans and pigs was subsequently used for downstream analyses to identify conserved key TFs. Fig. 5. Open in a new tab Analysis of key TFs in the major ZGA event for both human and pig. ( A ) soft clustering results of regulons recognized by SCENIC for both species. ( B ) heatmaps display the activity of regulons highly active during major ZGA across various stages, with high correlations between species. ( C ) Multi-bar charts show the relative relationship between the regulatory elements in B and the main ZGA genes in humans and pigs after homologous screening. U1.P4/p8/c represents the independent pig 4-cell portion, 8-cell portion, and common portion between 4-cell and 8-cell in U1. ( D ) PPI networks display the relationships between various factors and highly expressed TFs in C. ( E ) heatmaps illustrate the correspondence of 7 significant regulons across subsets of major ZGA genes. ( F ) trajectory plots display the expression of TFs corresponding to 6 significant regulons at various developmental stages in both humans and pigs. ( G ) Bar charts show the enrichment participation of the MYC regulon based on GO enrichment results from the U1 set, highlighting only the most highly involved parts. The first two sets of numbers in the annotations on each bar represent the number of genes in the MYC regulon that were enriched and the number of genes in the background U1 set, respectively. The same rule applies to all other similar types of charts By further selecting regulons with highly correlated interspecies activity patterns within the intersecting set, we identified 13 TFs that exhibited closely aligned activity during major ZGA, among which MYC, ATF3 , and ATF1 possessed the largest regulon sizes (Fig. 5 B). A second group of 11 TFs displayed moderately high interspecies correlation and included prominent contributors such as ETV3, SP1, USF1, GATA6 , and MYCN . In contrast, 11 TFs exhibited species-specific activity patterns during major ZGA. Among these, KLF8 was pig-specific, whereas YY1, RORA, REST, ETV4 , and EGR2 , each with regulon sizes exceeding 100 target genes, were specific to humans (Fig. 5 B). In addition to general transcriptional cofactors such as POLR2A and BRF2 , several enriched cofactors, including RAD21 , EZH2, KDM5B, HDAC2 , and BCLAF1 , were identified with comparatively large regulons and showed similar major ZGA-associated activity patterns across both species (Fig. 5 B). These TFs and cofactors were subsequently screened for sequence homology to support downstream comparative analyses. We observed that the regulons of MYC , ETV3 , SP1 , ATF3 , and USF1 exhibited the most extensive overlap with major ZGA-associated genes, markedly exceeding that of other TFs. The cofactor BCLAF1 also showed a relatively strong interaction, with a bias toward the pig 4-cell stage (Fig. 5 C). Among these, the MYC regulon demonstrated the largest intersection, primarily encompassing the U1 gene set and showing a pronounced enrichment in the pig 4-cell (P4) subset, while also significantly intersecting with stage-specific ZGA genes unique to pig 4-cell and human 8-cell stages (P4 and H8, respectively) (Fig. 5 C). ETV3 and SP1 preferentially overlapped with pig-specific ZGA genes, particularly those enriched at the 4-cell stage, and also showed substantial interactions with genes independently differentially expressed during the pig 8-cell versus 4-cell transition (P8) (Fig. 5 C). PPI analysis of the intersecting gene set revealed a network centered on MYC , suggesting potential direct interactions with highly expressed TFs such as FOXR1 , as well as with other TFs possessing large regulons, including ATF3 , SP1 , USF1 , ETV3 , and BCLAF1 . Notably, the network analysis indicated putative exclusive interactions between MYC and ETV3 , and between MYC and BCLAF1 (Fig. 5 D). The intersections among the six key regulons within the major ZGA-associated DEGs were predominantly concentrated in the U1 gene set and were characterized by largely independent interaction patterns (Fig. 5 E). The MYC regulon showed prominent individual interactions with ETV3 , SP1 , ATF3 , and ATF1 . In contrast, in pigs, intersections identified during the 8-cell versus 4-cell DEG comparison were mainly characterized by ETV3 interacting independently with ATF3 and USF1 (Fig. 5 E). Notably, the majority of the BCLAF1 regulon within the U1 set overlapped with the MYC regulon, forming a more extensive interaction network than that observed between MYC and other regulatory factors (Fig. 5 E). TF expression trajectories generally mirrored their corresponding regulon activities and their interactions with species-specific major ZGA-associated DEGs. USF1 and ETV3 exhibited weaker upregulation during major ZGA in humans compared with pigs, whereas ATF1 , ATF3 , and MYC showed highly similar expression dynamics across species. ETV3 displayed higher maternal expression in pigs than in humans. ATF1 showed comparable maternal deposition in both species, while ATF3 , USF1 , and MYC exhibited little to no maternal mRNA reserve (Fig. 5 F). Notably, SP1 in humans showed substantial maternal deposition and maintained relatively stable expression levels with only modest upregulation at the 8-cell stage, in contrast to the pronounced upregulation observed in pigs (Fig. 5 F). Similarly, the cofactor BCLAF1 was extremely abundant in human maternal mRNA, at levels far exceeding those of SP1 , and exhibited less upregulation during major ZGA compared with pigs (Fig. 5 F). We observed a significant correlation between BCLAF1 and MYC in ChIP-seq datasets derived from two human cancer cell lines: GM12878 (a lymphoblastoid cell line immortalized with Epstein–Barr virus, derived from a chronic myelogenous leukemia patient) and K562 (an erythroleukemia cell line from a chronic myelogenous leukemia patient in blast crisis). These cell lines represent distinct biological contexts. GM12878 models transcriptionally active, immortalized lymphoid cells with relatively stable proliferative behavior, whereas K562 represents an aggressive, poorly differentiated myeloid leukemia characterized by strong oncogenic signaling. Despite their cancer origin, both datasets revealed co-localization of BCLAF1 and MYC at shared genomic regions (e.g., the U4/U1 region), along with overlapping target genes, including BCLAF1 , MYC , and genes involved in transcription ( POLR2A ), RNA processing ( PNN , DDX5 ), and oncogenic regulation ( ZNF217 ) (Supplementary Fig. S17 A, B). These findings suggest a potentially conserved regulatory interplay between BCLAF1 and MYC that may operate in rapidly proliferating or transcriptionally dynamic cellular systems. Nevertheless, extrapolation of these observations to early embryonic processes such as ZGA should be interpreted with caution, given the substantial differences in chromatin organization, transcriptional programs, and regulatory network architecture between cancer cells and early embryos. In pigs, where experimental ChIP-seq data are currently unavailable, motif scanning analyses also suggested potential co-localization of BCLAF1 and MYC on a subset of U4 genes (Supplementary Fig. S17 C). Although AlphaFold3 yielded low simulation scores for interactions between BCLAF1 and MYC/MAX , this interaction has been reported in the BioGRID database. Consistently, protein–protein docking predictions using HDOCK indicated a high-confidence interaction between BCLAF1 and MYC (Supplementary Fig. S18 A, B). Beyond its transcriptional regulatory role, BCLAF1 is an important post-transcriptional regulator. Previous studies have shown that BCLAF1 interacts with SRSF1 to regulate MYC mRNA maturation [ 24 – 26 ]. BioGRID also documents direct physical interactions between BCLAF1 and several other serine/arginine-rich splicing factor (SRSF) family members, including SRSF6 , SRSF7 , SRSF9 , and SRSF11 . Although AlphaFold3 predicted a low binding score (ipTM) for the BCLAF1–SRSF1 interaction, these two factors exhibited highly concordant expression trajectories during the major ZGA period (Supplementary Fig. S17 D). In parallel, HDOCK docking analysis predicted a high-confidence interaction between BCLAF1 and SRSF1 (Supplementary Fig. S18 C, D). In addition, we extracted putative BCLAF1 -binding RNAs from the POSTAR database and observed a strong enrichment within the U1 and U4 gene sets. Among these, BCLAF1 , EIF4G1 , MATR3 , and DDX5 were identified as the most likely RNA targets (Supplementary Fig. S17 E, F). Notably, several of these targets, including PNN and MATR3 , also exhibit PPIs with BCLAF1 . Collectively, these findings suggest that DDX5 may represent a stable BCLAF1 target conserved across species, potentially subject to BCLAF1 -mediated regulation at both transcriptional and post-transcriptional levels. Network analysis revealed that the MYC regulon, frequently correlated with other U1 genes within high-degree network regions, was associated with several highly interconnected nodes, including POLR2B , DHX15 , and IMP3 , which exhibited slightly higher degree distributions as illustrated by ridge plots (Supplementary Fig. S16 C, D). Consistent with U1 functional enrichment, genes within the MYC regulon were significantly associated with GO processes related to nucleosome complex assembly and RNA processing, both of which represent core functional themes of the U1 gene set. Notably, the BCLAF1 regulon also showed substantial enrichment in these same processes (Fig. 5 G). Genes regulated by MYC were involved across all U1-associated KEGG pathways, including autophagy and energy metabolism (Supplementary Fig. S16 E). In particular, the MYC regulon showed strong enrichment for spliceosome-related functions, encompassing genes encoding essential spliceosomal components. These included key U1 members such as DDX5 (p68), DHX15 (Prp43), DHX8 (Prp22), SNRPB (Sm), and SRSF1/3/6/7 (SR proteins), all predominantly represented in the U1.p4 subset. Among these, DHX8 was jointly targeted by MYC and BCLAF1 , whereas DDX5 was targeted by five regulators (excluding ETV3 and USF2 ). In contrast, DHX15 was predicted to be uniquely regulated by MYC (Supplementary Fig. S19 A, B). Furthermore, within the P4 and P8 gene subsets, the MYC regulon was significantly enriched for genes involved in ubiquitin-mediated proteolysis (Supplementary Fig. S20 A). In addition to RNA splicing, both SP1 and ETV3 regulons were significantly enriched in a broad range of metabolic processes, including RNA metabolism, glycoprotein metabolism, amino acid and steroid metabolism, as well as phospholipid and glycolipid metabolism (Supplementary Fig. S20 B, C). A total of 24 major ZGA genes were predicted to be co-regulated by at least four TFs from the set of significant regulons. Among these, KAT5 and FADS2 in the P8 subset were predicted to be regulated independently of MYC . The most highly influenced genes within the common U1 subset and the pig 8-cell-specific group, such as PLA2G15 and CALR (independent of ATF1 ), exhibited corresponding significant upregulation at their respective developmental stages (Supplementary Fig. S20 D, E). Additionally, 32 TFs exhibited species-specific upregulation during the major ZGA, including 17 pig-specific TFs (such as SP2 , SP8 , and MAX ). Among these, GTF2F1 and ZNF384 , within their highly homologous regions, showed significant interactions with the P4 gene set (Supplementary Fig. S21 A, B). In contrast, 15 TFs were human-specific, including SP6 and KLF5 , with TBL1XR1 displaying the most prominent interactions with the H8 gene set within its highly homologous region (Supplementary Fig. S21 A, B). Functional enrichment analysis revealed that, in addition to shared processes related to RNA and protein metabolism, ZNF384 at P4 was specifically associated with histone modification, whereas GTF2F1 was additionally involved in glycoprotein metabolism. Notably, enrichment of TBL1XR1 at H8 was predominantly linked to DNA repair–related processes (Supplementary Fig. S21 C). Conservation of significant maternal TFs during early phases Given the crucial role of maternal TFs in driving the minor ZGA phase and shaping the subsequent major ZGA, we systematically examined TFs present in maternal mRNAs. Using a threshold of log_e (TPM + 1) = 2 to distinguish actively expressed TFs from those with minimal expression, we found that active TFs accounted for approximately 34% of all TFs in humans and 27% in pigs (Supplementary Fig. S22 A), with 231 TFs shared between the two species (Fig. 6 A). Overall, the expression trajectories of these TFs were largely conserved between humans and pigs (Supplementary Fig. S22 B), displaying modest fluctuations during the minor ZGA phase. These TFs remained prevalent during early development, constituting 0.66% of total expression in humans and 0.87% in pig oocytes, and persisted until the onset of the major ZGA (at the 4-cell stage in pigs). Although their expression gradually declined thereafter, they continued to represent a substantial proportion of the transcriptome through the morula stage (Fig. 6 B). Fig. 6. Open in a new tab Analysis of high-expression TFs in maternal mRNA across species. ( A ) Euler diagrams display the intersections between significantly expressed genes in humans and pigs (left) and the top 20 expressed genes in both species, highlighting those with a log2 FC (log2FC) ≤ 1.5 (right). ( B ) column charts illustrate the proportions of significantly expressed TFs in humans and pigs within all TFs. ( C ) soft clustering results show the expression patterns of significantly expressed TFs from the oocyte to morula stages for both species. ( D ) venn diagrams illustrate the interrelationships among sets for similar patterns obtained from soft clustering in Fig. C and highly homologous TFs with no significant difference between humans and pigs. ( E ) expression trajectories of Set1 across developmental stages in both species. ( F ) expression trajectories of CTNNB1 and LEF1 across stages in both species Refined expression analysis highlighted TFs prominently active during early developmental stages in both humans and pigs. We observed a strong cross-species concordance in expression trends, characterized by early stability followed by downregulation approaching the major ZGA. Specifically, pig clusters 4 and 5 aligned with human clusters 5 and 6 (Fig. 6 C–D). Further filtering based on high homology and the absence of significant interspecies differences yielded 23 genes (Set1) (Supplementary Fig. S22 D), including CLOCK , MXD1 , SOX30 , E2F8 , LEF1 , SALL2 , TCF7L1 , and TADA2A . These genes exhibited highly consistent expression patterns across both species from the oocyte to morula stages (Fig. 6 E). Despite this overall consistency, PPI analysis revealed limited connectivity within Set1, with interactions restricted to two small groups: LEF1 , EBF1 , and TCF7L1 in one group, and TADA2A and CDC5L in the other (Supplementary Fig. S22 E). Notably, LEF1 emerged as the most interconnected factor within Set1 when considering interactions with other significant TFs, particularly those involved in chromatin remodeling and Polycomb complex regulation, including HDAC1 , EP300 , EZH2 , and ACTB (Supplementary Fig. S22 E). The top portion of Set1 exhibited extremely high expression levels, with the top 10 expressed factors being consistent between humans and pigs, notably including THYN1 and LEF1 among the top five in both species (Supplementary Fig. S22 F). As key components of the Wnt signaling pathway, LEF1 and TCF7L1 frequently interact with CTNNB1 (Supplementary Fig. S22 G). CTNNB1 showed high initial expression, particularly in pigs, followed by downregulation during the major ZGA phase, with a slight rebound observed in pigs thereafter (Fig. 6 F). In contrast to LEF1 , which continues to decline from the major ZGA to the morula stage, TCF7L1 in pigs displayed a pronounced rebound in expression at the morula stage, more closely paralleling the expression trajectory of CTNNB1 (Fig. 6 F). All three genes were commonly enriched in Wnt-related GO processes, with LEF1 additionally involved in histone modification, particularly acetylation, and DNA binding (Fig. 7 A). Fig. 7. Open in a new tab Functional and expression analysis of conserved maternal mRNAs across species. ( A ) Bar charts based on GO enrichment results from common maternal mRNAs showcase the top 6 significant items in categories involved by LEF1 . ( B ) Bar charts display the expression levels of the top 20 genes and their log2FC in the differential analysis of maternal mRNAs, in humans(left) and pigs(right). ( C ) expression trajectories of the intersection set of 12 genes with no significant difference from the top 20 genes for both species. ( D ) Bar charts based on GO enrichment results from common maternal mRNAs present the top 6 significant GO items related to acetylation and mRNA involved by genes of Set2 We explored TCF7 and TCF7L2 alongside TCF7L1 , noting higher expression levels in humans than in pigs. TCF7L2 exhibited the lowest initial expression in both species, while TCF7 showed intermediate expression in pigs (Supplementary Fig. S22 H). The expression trajectories of TCF7 and TCF7L1 were closely aligned, particularly in humans. All three genes, TCF7, TCF7L1, and TCF7L2 , underwent marked downregulation during the 8-cell vs. 4-cell transition in both species (Supplementary Fig. S22 H). We scrutinized the top 20 TFs present in maternal mRNAs, which markedly exceeded other TFs in abundance, accounting for more than 25% of all TFs in early human embryos and over 50% in pig embryos. These TFs were notably downregulated by the major ZGA stage (4-cell stage in pigs) (Supplementary Fig. S23 A–B; Supplementary Table S3 ). The intersection between species included TAX1BP1, LYAR, CSDE1, and OTX2 , along with additional factors, ARID4A, TADA2A, ZNF706, PINX1, HMGB2, THYN1, HHEX, and NCOR1 ,that showed no significant interspecies differences (defined as Set2). Among these, HHEX, THYN1, and TADA2A overlapped with Set1 (Figs. 6 A and 7 B; Supplementary Fig. S22 F). Notably, PAX6, MIXL1, and ID3 exhibited high interspecies expression variance. In particular, ID3 showed exceptionally high expression in pigs, clearly surpassing all other TFs in both humans and pigs, while being significantly downregulated during the major ZGA. In contrast to PAX6 and MIXL1 , ID3 maintained moderate expression levels in humans (Fig. 7 B; Supplementary Fig. S23 B–C). Set2 TFs, including TAX1BP1, LYAR, CSDE1, NCOR1, PINX1, and OTX2 , displayed strong expression correlations during early development in both species (Fig. 7 C). Functional enrichment analysis revealed that Set2 TFs were primarily involved in organelle localization, stemness maintenance, and cell-cycle regulation, with notable enrichment in pathways related to histone acetylation, mRNA processing, and metabolism (Fig. 7 D; Supplementary Fig. S23 D). We conducted a more detailed analysis of ID3 , whose homologs ID1, ID2, and ID4 showed no significant expression during early embryonic development (Supplementary Fig. S23 E). ID3 generally functions as a transcriptional inhibitor by interacting with other regulatory factors. In the PPI network of active TFs in pigs, ID3 directly interacts only with FOXO1, HIF1A, KLF2, TCF3/7, and E2F4 , but it indirectly influences numerous TFs with high network connectivity. Among these, STAT3 and YY1 , which exhibit the highest degrees of connectivity, are also closely associated with FOXO1, HIF1A, TCF3, and E2F4 (Supplementary Fig. S23 F). Although AlphaFold3 predictions yielded low confidence scores for most interactions except for TCF3 , the HDOCK predictions indicated high confidence scores for these PPIs (Supplementary Fig. S25 ). The TFs directly and indirectly influenced by ID3 represent common active regulators in both species and are involved in key biological processes, including chromosome segregation and nuclear division (the most significantly enriched terms), as well as histone modification, which involved the largest number of genes. These processes were prominently enriched in GO analyses (Supplementary Fig. S23 G and H). By comparing the expression patterns of ID3 ’s direct interaction factors across species (Supplementary Fig. S24 A), we selected downstream target genes for further analysis, with FOXO1, HIF1A, TCF3, and E2F4 defined as core regulators. Target genes were required to be regulated by at least two of these factors, to show positive expression correlations with at least two upstream regulators, and to be negatively correlated with ID3 . The resulting target gene set comprised nearly 1000 genes that were upregulated during the major ZGA period in pigs, with a strong bias toward the P8 stage. In comparison with humans, a greater number of these target genes were significantly upregulated at earlier developmental stages (1-cell and 2-cell) (Supplementary Fig. S24 B and C). Enrichment analysis of these ZGA-associated target genes revealed that, during early ZGA in humans and major ZGA in pigs, they were predominantly involved in metabolic processes and cell growth control, including oxidative phosphorylation, lipid metabolism, ferroptosis, apoptosis, and autophagy (Supplementary Fig. S24 D and Supplementary Fig. S26 A and B). Metabolic pathways were particularly prominent in the enrichment results for the pig major ZGA stage. Heatmap analysis showed that target genes involved in ferroptosis, lipid metabolism, and PRC–related pathways were generally positively correlated with YY1, TCF3, E2F4, and FOXO1 , while being negatively correlated or uncorrelated with HIF3A, STAT3, and ID3 . Representative genes exhibiting this pattern included JARID2, GCLM, and FASN . Additional target genes such as SMARCA4, FTH1, FTL, and HADHA displayed similar expression trends (Supplementary Fig. S26 D). These genes are predicted to be jointly regulated by STAT3, YY1, FOXO1, HIF1A, TCF3, and E2F4 (Supplementary Fig. S26 C). Furthermore, Wnt signaling–related genes, including LEF1, CTNNB1, and TCF7L1 , showed strong positive correlations with all seven factors, particularly with HIF1A, STAT3, and ID3 (Supplementary Fig. S26 C). By integrating regulon activity analyses, we identified TFs associated with early active regulons. Twelve factors exhibited high interspecies correlation in expression patterns and comparable activity during the minor ZGA period, including SOX30, TCF7L1, and BNC1 (Set 1), as well as notable factors such as FOXC1 and ELK1/3 . Following homology filtering, GRHL2 showed the largest overlap with early DEGs (Supplementary Fig. S27 A). Among species-specific early active factors, 18 TFs corresponded to pig-specific high-activity regulons, including FOXP1/N3, GATA4, NFKB1, YY2, and MLX . Of these, the regulons of HMGA2, ZBTB33, and MLX exhibited the greatest overlap with early DEGs, with HMGA2 and ZBTB33 showing a strong bias toward pig-specific DEGs (Supplementary Fig. S27 A and B). Similarly, humans displayed 17 unique early active regulons corresponding to TFs such as SOX1, FOXO3, REL, and PLAG1/2 . Among these, LHX4, FOXO3, and PLAG1 demonstrated the largest overlap with early DEGs and were preferentially associated with the human gene set (Supplementary Fig. S27 A and B). Despite these species-specific differences in regulon activity, PPI analysis revealed frequent interactions among these TFs. Notably, HMGA2 showed interactions with KAT2A and the highly connected factor NFKB1 , suggesting potential cross-regulatory or cooperative roles during early embryonic development (Supplementary Fig. S27 C). Pattern of MT genes between species Prior to implantation, embryos rely exclusively on MT energy production, and our results reveal clear interspecies differences in MT gene expression. Heatmap analyses show that the majority of the 13 MT genes exhibit higher expression levels in humans across nearly all developmental stages, with the exception of the oocyte stage. In contrast, ATP6 displays consistently higher expression in pigs, while ATP8 shows pronounced stage-specific heterogeneity uniquely in pig embryos (Fig. 8 Aand B). Fig. 8. Open in a new tab Expression patterns of MT genes across species. ( A ) heatmap showing the variation of 13 MT gene expression at various stages between species. ( B ) Characterization plot demonstrating the high heterogeneity of ATP8 in pigs. ( C ) PPI demonstrating the interactions of MT genes. ( D ) clustering of MT genes in pigs. ( E ) MT genes were compared with human ones in terms of expression patterns according to each of the four pig categories obtained from Fig. D These MT genes can be broadly classified into three functional categories, collectively contributing to the assembly of all respiratory chain complexes except complex II (Fig. 8 Cand Supplementary Fig. S28 A). Comparative analysis of MT gene expression revealed a more diverse and dynamic expression pattern in pigs, whereas human MT gene expression was relatively uniform across developmental stages (Fig. 8 D and Supplementary Fig. S22 B). In pigs, MT genes generally exhibited downregulation at the 1-cell versus oocyte stage. In contrast, human MT genes showed a distinct and pronounced upregulation at this stage, with 7 of the 13 MT genes not increasing in pigs until the 2-cell stage. Notably, almost all MT genes in both species were upregulated during the 4-cell versus 2-cell transition. This upregulation was particularly prominent for CO3 , ND2/ND6 , and ATP8 , which displayed major ZGA-specific activation (Fig. 8 D and Supplementary Fig. S28 B and C). Furthermore, we observed substantial maternal mRNA reserves of the upstream MT TF TFB1M in both species, while humans exhibited significantly higher maternal levels of TFAM , TFB2M , and MTERF4 compared with pigs. Importantly, all five MT regulatory factors were markedly upregulated during the major ZGA, beginning at the 4-cell stage in pigs (Supplementary Fig. S28 D). Additionally, despite differences in overall expression dynamics, specific MT gene pairs, such as ND1 and ND4 , as well as COX1/COX2 and ATP6 , exhibited strong intraspecies homogeneity in both expression levels and temporal expression patterns (Fig. 8 E). Discussion In this study, we integrated single-cell RNA sequencing data from 139 samples spanning developmental stages from the oocyte to the morula in both humans and pigs, and systematically compared transcriptomic dynamics during the early embryonic MZT. Our analysis revealed a high degree of gene expression similarity between human and pig embryos throughout the MZT, alongside pronounced species-specific differences. Notably, during the major ZGA stage, both species exhibited robust upregulation of key TFs and coordinated gene expression programs. Furthermore, we identified both conserved and species-specific genes associated with MT function, chromatin remodeling, and metabolic pathways, highlighting shared regulatory principles as well as divergent developmental strategies between the two species. The sets of ZGA genes identified at the two major ZGA stages in pigs (4-cell and 8-cell) show a high degree of overlap, whereas a comparable pattern in humans is observed between the 1-cell and 8-cell stages. Notably, three overlapping human genes also correspond to major ZGA genes in pigs, revealing a distinctive “two-stage rise” pattern that spans a broader developmental window in humans. Among these, SRSF1 , EFTUD2 , and IMP4 are of particular interest. SRSF1 has been identified as a maternally deposited gene with high expression in critical regions such as the head and central nervous system during early development [ 27 ]. Loss of SRSF1 results in severe embryonic abnormalities, including defects in cardiac and cartilage formation, underscoring its essential role in embryonic patterning [ 27 ]. Mutations in EFTUD2 are associated with congenital disorders such as microphthalmia and craniofacial abnormalities, and studies in zebrafish have highlighted its indispensable function during embryogenesis [ 28 ]. IMPα2 is required for the nuclear localization of Oct3/4, a core pluripotency factor essential for maintaining mouse embryonic stem cell identity [ 29 , 30 ]. In contrast, overexpression of IMPα4 leads to reduced Oct3/4 levels, suggesting a role in promoting cellular differentiation [ 29 ]. Although these genes are clearly involved in ZGA, their precise regulatory roles and species-specific mechanisms remain to be fully elucidated. Furthermore, it is noteworthy that the initiation of human ZGA exhibits a pronounced paternal genome bias, a defining characteristic that distinguishes human embryos from those of other mammals, including pigs and mice [ 31 ]. Central to our comparative analysis is the alignment of specific embryonic stages between pigs and humans, namely the pig 4-cell (P4) and 8-cell (P8) stages with the human 8-cell (H8) stage, as representative phases of major ZGA. This alignment is primarily based on conserved transcriptomic dynamics and the peak of de novo gene activation observed at these stages in each species, consistent with prior reports. However, we emphasize that this correspondence represents a functional alignment grounded in transcriptional similarity rather than a strict chronological equivalence. Important biological distinctions underlie this comparison. Notably, the initiation of ZGA in human embryos exhibits a pronounced paternal genome bias, a defining characteristic that distinguishes humans from other mammals, including pigs. Moreover, although these stages correspond to periods of extensive transcriptional activation, the underlying molecular triggers, epigenetic reprogramming kinetics, and developmental competencies may differ substantially between species. Consequently, while this alignment facilitates the identification of conserved and species-specific features of ZGA, interpretations regarding the homology of regulatory mechanisms or the precise timing of developmental events should be made with caution. During the major ZGA stages, pig and human embryos exhibit distinct species-specific characteristics. In pigs, the 4-cell and 8-cell stages display 1566 and 1883 uniquely expressed genes, respectively, whereas the human 8-cell stage contains 1364 unique gene expressions. Several chromatin remodeling–related genes, including EP300 , EP400 , MTA2 , and KAT5 , show pronounced species specificity during pig ZGA and maintain high expression levels across both major ZGA stages. EP300 plays a pivotal role in early embryonic development through its function in chromatin remodeling. It cooperates with the TF KLF4 to activate the pluripotency-associated gene POU5F1 [ 32 ] and enhances cellular reprogramming efficiency by inhibiting its bromodomain, thereby contributing to the suppression of somatic gene expression [ 33 ]. Additionally, interactions between EP300 and SMYD3 connect EP300 to epigenetic reprogramming pathways [ 34 ]. EP400 is essential for maintaining oocyte quality, and its deficiency results in abnormal embryonic development. It forms a complex with the TF NFYA , which regulates ZGA-associated gene expression [ 35 ]. Moreover, the Tip60/EP400 complex is critical for metabolic regulation and proliferative capacity in neural progenitor cells [ 36 ]. MTA2 is a core component of the NuRD complex, where it enhances chromatin remodeling capacity and contributes to transcriptional regulation during early embryogenesis [ 37 ]. Beyond development, MTA2 is also implicated in the invasion and metastasis of gastric cancer cells, with its expression regulated by the TF SP1 . SP1 binds directly to the MTA2 promoter, promoting its transcriptional activity [ 38 ]. KAT5 (also known as TIP60 ) is another key regulator of chromatin remodeling and self-renewal during early embryonic development [ 39 ]. Although KAT5 -deficient embryonic stem cells exhibit minimal global transcriptional changes, their differentiation into mesodermal and endodermal lineages is severely impaired [ 40 ]. Furthermore, the TIP55 isoform of KAT5 is indispensable for organogenesis, as mutant embryos display cardiac and neural tube defects at embryonic day 11.5 [ 41 ]. Analysis of GTRD data revealed that these genes are associated with chromatin-active regions; however, their expression patterns frequently exhibit uncoupling between the two ZGA stages. Based on these observations, we hypothesize that these chromatin remodeling genes may continue to function throughout ZGA by forming stage-specific epigenetic modification sites through interactions with stage-specific TFs, such as TP53 [ 42 ]. In human and pig embryos, shared maternal mRNA plays an important role in RNA processing and In both human and pig embryos, shared maternal mRNAs play critical roles in RNA processing and histone modification, particularly in histone H3 regulation. Among these modifications, methylation of H3K4 and H3K9 is repeatedly implicated during early embryonic development. Most genes associated with core histone modifications are conserved between the two species. Notably, the expression levels of HDAC1 and HDAC2 are significantly higher in pig embryos than in humans, whereas EP300 is more abundant in human maternal reserves but becomes activated during the major ZGA stage in pigs, thereby reducing interspecies differences. The elevated expression of HDAC1 and HDAC2 in pig embryos suggests that stronger deacetylation control may be required during porcine embryogenesis. Histone deacetylation regulates chromatin structure by removing acetyl groups from histones, resulting in chromatin compaction and transcriptional repression. Previous studies have demonstrated that HDAC1/2 are essential for ZGA in mouse and bovine embryos. During ZGA, HDAC1/2 help maintain an appropriate balance between chromatin activation and repression by ensuring controlled deacetylation [ 43 ]. Loss of HDAC1/2 activity impairs the removal of acetylation marks such as H3K27ac, disrupts gene activation, and ultimately leads to embryonic developmental arrest [ 44 ]. Thus, HDAC1/2 -mediated deacetylation is indispensable for regulating early embryonic gene expression and maintaining a balanced chromatin state conducive to normal development. Furthermore, the similar reprogramming patterns of key histone modifications, including H3K4me3 and H3K27me3, observed in both pigs and humans underscore the epigenetic relevance of pig embryos as a model for studying early human development [ 45 ]. KEGG pathway analysis revealed enrichment of ATP-dependent chromatin remodeling complexes and PRCs in both species. Several co-regulated genes, including CBX2/8 and SMARCA2/4 , exhibit interspecies activity-switching patterns. Specifically, CBX2 and SMARCA2 display high expression during minor ZGA in humans followed by degradation during major ZGA, whereas the opposite pattern is observed in pig embryos. We hypothesize that these genes may adapt to align with the species-specific timing of ZGA gene activation. CBX2/8 and SMARCA2/4 show opposing regulatory dynamics across species, with important implications for ZGA-associated gene regulation. CBX2/8 participate in chromatin remodeling, while SMARCA2/4 function as core components of ATP-dependent chromatin remodeling complexes, each exhibiting distinct regulatory strategies during ZGA in different species. For instance, SMARCA5 has been shown to initiate ZGA in mouse embryos by modulating chromatin accessibility [ 46 ], whereas CBX2/8 may regulate gene expression through alternative mechanisms in other species. These interspecies differences highlight the species-specific timing and regulatory mechanisms underlying ZGA. Overall, the initiation of embryonic gene expression is tightly linked to the temporal activity switching of regulatory proteins, emphasizing the multilayered and evolutionarily diversified nature of ZGA regulation [ 47 ]. TFs play a pivotal role in ZGA during early embryonic development. During the minor ZGA stage in both humans and pigs, WNT pathway–related genes, including LEF1 , TCF7L1 , and CTNNB1 , are co-upregulated and implicated in chromatin remodeling, suggesting that they may perform conserved regulatory functions during early embryogenesis. LEF1 and TCF7L1 (formerly TCF3 ) belong to the TCF/LEF TF family and typically bind β-catenin to mediate WNT signaling–dependent transcriptional regulation. Previous studies have shown that LEF1 activity in osteosarcoma is governed by complex regulatory mechanisms; for example, circMYO10 promotes chromatin remodeling and enhances the transcriptional activity of the β-catenin/ LEF1 complex by regulating multiple molecular intermediates, providing insight into the potential regulatory mechanisms of LEF1 [ 48 ]. In contrast, binding of TCF7L1 to β-catenin does not activate transcription; instead, β-catenin inhibits TCF7L1 activity by facilitating its dissociation from DNA and promoting subsequent protein degradation [ 49 , 50 ]. TCF7 functions similarly to LEF1 and TCF7L1 within the WNT pathway, generally acting as a β-catenin coactivator. However, in pig embryos, TCF7 expression is markedly reduced, suggesting that its role may be functionally compensated by other regulatory factors. This shift may alter chromatin remodeling dynamics and, consequently, influence the process of genome activation [ 51 ]. Gain-of-function mutations in CTNNB1 (e.g., S33Y, S37A, T41A, and S45F) are known to induce sustained β-catenin activation and promote oncogenesis [ 52 ]. Nevertheless, despite similar expression patterns of these WNT pathway components during early embryogenesis, their precise regulatory mechanisms may differ across species due to variations in regulatory networks and cellular contexts [ 53 ]. During the minor ZGA stage, the expression of the inhibitory protein ID3 is significantly upregulated in pig embryos. ID3 likely exerts its effects by binding to other TFs and preventing their interaction with DNA, thereby maintaining embryonic cells in an undifferentiated state, delaying premature differentiation, and contributing to chromatin remodeling regulation [ 54 ]. In contrast, ID3 expression remains relatively low in human embryos, whereas ID2 displays higher activity, potentially compensating for the reduced ID3 expression and sustaining the transcriptional regulatory balance required for early development [ 55 ]. Thus, elevated ID3 expression in pig embryos appears to regulate cell differentiation primarily through inhibition of TF activity, while ID2 activity in human embryos may serve a compensatory role. Together, the differential expression and functional interplay of ID family members contribute to the coordinated regulation of embryonic development across species. During the major ZGA period, MYC , ETV3 , and SP1 were identified as highly active and key TFs in both human and pig embryos, exhibiting broadly similar expression patterns. Among these, MYC functions as a central regulator of ZGA. It not only plays a pivotal role in transcriptional activation but also participates in the regulation of mRNA processing and ribosome biogenesis, particularly targeting genes closely associated with the U1 small nuclear ribonucleoprotein complex [ 56 ]. Notably, MYC shows relatively limited overlap in its regulatory targets with other TFs, such as BCLAF1 and ETV3 , suggesting that MYC governs gene expression through distinct and specialized regulatory pathways. BCLAF1 physically interacts with MYC and serves an important role in post-transcriptional regulation, especially in the maturation and processing of MYC mRNA [ 57 ]. Although computational predictions indicate a relatively low binding confidence for this interaction, it is likely that effective regulation requires cooperative interactions with additional cofactors. Functional distinctions were also observed between ETV3 and SP1 . SP1 appears to contribute to the initiation of ZGA by activating the transcription of specific early-response genes, such as HSP70.1 [ 58 ], whereas ETV3 is more strongly associated with the fine-tuning of gene expression, particularly during defined developmental windows [ 59 ]. As the core regulatory factor during ZGA, MYC enhances the complexity and robustness of the gene regulatory network through its interactions with cofactors such as BCLAF1 . In this context, BCLAF1 not only modulates transcriptional activity but also supports MYC -mediated regulation via post-transcriptional mechanisms, while SP1 and ETV3 exert distinct yet complementary roles in orchestrating precise gene expression programs during early embryonic development. Pig embryos primarily rely on lipid metabolism, whereas human embryos tend to depend more on carbohydrate metabolism during early development [ 60 ]. Enrichment analysis of minor ZGA genes revealed that pig embryos exhibit significant upregulation of lipid metabolism – and ferroptosis-related genes, including HADH , HADHA , HADHB , and HACD4 . The expression levels of these genes are markedly higher in pigs than in humans, particularly for the ferritin genes FTL and FTH1 , which show substantially elevated expression in pigs and opposite developmental trends between the two species. The pronounced upregulation of HADH , HADHA , HADHB , and HACD4 in pig embryos indicates a stronger reliance on fatty acid β-oxidation for energy production [ 61 ]. These genes encode key enzymes involved in lipid catabolism, facilitating the breakdown of long-chain fatty acids to generate ATP [ 62 ]. This metabolic preference is consistent with the abundant lipid droplets present in pig oocytes and embryos, which serve as endogenous energy reserves supporting early embryonic development [ 63 ]. In contrast, human embryos rely more heavily on carbohydrate metabolism, primarily through glycolysis and the tricarboxylic acid cycle, and depend to a greater extent on exogenous glucose to meet the energetic demands of rapid cell division and growth [ 64 ]. In addition, the significantly higher expression of FTL (ferritin light chain) and FTH1 (ferritin heavy chain) in pig embryos suggests a species-specific strategy for iron storage and ferroptosis regulation. FTH1 , which possesses ferroxidase activity, facilitates the conversion of Fe 2 + to Fe 3 + , enabling safe iron sequestration, while FTL primarily contributes to iron storage capacity [ 65 ]. The elevated expression of these genes implies a greater requirement for iron buffering in pig embryos, potentially mitigating oxidative stress and ferroptosis during early development [ 66 ]. This iron-centered metabolic strategy appears distinct from that in human (and mouse) embryos and may enhance developmental robustness by coupling iron storage with antioxidant defense mechanisms [ 67 ]. Supporting this notion, additional iron metabolism–related genes, such as TFRC and STEAP3 , also exhibit divergent expression patterns between pigs and humans. Regarding MT regulation, ZGA-associated genes in human embryos are predominantly linked to MT function, with significant upregulation of MT genes ( MT genes ) such as TFAM , TFB1M , and TFB2M . These genes regulate MT biogenesis, DNA replication, and protein transport, underscoring the central role of mitochondria in supporting early human embryonic development. In contrast, MT gene expression in pig embryos displays greater diversity and temporal variation. The coordinated upregulation of MT genes in human embryos is accompanied by increased MT mass and activity, ensuring adequate energy supply during critical stages of cell proliferation and differentiation [ 68 ]. Previous studies have also shown that increased MT DNA ( mtDNA ) content and elevated expression of MT TFs, such as TFAM and NRF1 , enhance MT protein import and ATP production [ 69 ], reinforcing the importance of MT regulation in early embryogenesis. Collectively, these findings highlight fundamental metabolic differences between species during early embryonic development, particularly the reliance of pig embryos on lipid and iron metabolism versus the carbohydrate- and mitochondria-centered metabolic strategy in humans. These distinct metabolic programs likely reflect species-specific evolutionary adaptations that support early developmental requirements [ 70 , 71 ]. In summary, this study integrated single-cell RNA sequencing data from human and pig embryos spanning the oocyte to morula stages. Our results demonstrate that transcriptomic differences between the two species are most pronounced at the 4-cell stage and progressively diminish as development proceeds. Although stage-specific marker genes differ between species, corresponding developmental stages retain substantial functional similarities. Core chromatin regulators are functionally conserved but display species-specific expression dynamics; for example, RNF2 and several CBX family genes exhibit opposite regulatory trends during the major ZGA stage in humans and pigs. Active maternal mRNA decay during major ZGA (8-cell in humans and 4-cell in pigs), particularly involving pathways related to the cell cycle and RNA processing, emerges as a key mechanism underlying the increasing transcriptomic convergence between the two species. In addition, we identified critical TFs shared during major ZGA, such as MYC, alongside species-specific regulators, including KLF8 in pigs and YY1 and RORA in humans. MT gene expression also exhibits species-specific patterns: pigs show more dynamic MT transcriptional changes, whereas humans display relatively uniform expression with a distinctive upregulation at the 1-cell stage. Early embryonic development in pigs and humans is governed by largely conserved biological events, reflecting evolutionary preservation of core regulatory programs and underscoring universal principles of mammalian embryogenesis. However, the specific genes involved and the temporal dynamics of these events differ between species. These divergences highlight lineage-specific adaptive modifications that underlie species-specific developmental timing and regulatory fine-tuning, such as differences in the onset and progression of ZGA. Collectively, these findings strengthen the value of porcine embryos as a translational model for human embryology. Nevertheless, there are also several limitations of this study. First, the present work is based entirely on integrative computational analyses of existing single-cell transcriptomic datasets and therefore remains predominantly descriptive and correlative in nature. While the observed expression dynamics and cross-species similarities provide valuable insights, they do not establish causal or mechanistic relationships. Experimental validation will be required in future studies to directly assess the functional roles of the identified genes and regulatory pathways. Second, the conservation discussed in this study is confined to the transcriptomic level. Similar expression trajectories or shared pathway enrichments between human and pig embryos do not necessarily imply conserved regulatory mechanisms or equivalent biological functions. Instead, our analyses define a framework of transcriptional correspondence that serves as a foundational resource for hypothesis generation and for guiding future functional and mechanistic investigations across species. Materials and Methods Ethics statement All animal-related procedures followed the guidelines set forth by the Institutional Animal Care and Use Committees (IACUCs) of Northeast Agriculture University (NO. NEAUEC20190118). In vivo embryo collection and single cell separation Pig embryos used in this study were collected from the same batch of Bama sows with comparable and appropriate mating ages (Supplementary Table S1). After two rounds of mating at a 12 h interval, uteri were harvested, and embryos were recovered by uterine flushing at defined developmental time points. Oocytes were collected 24 h after the onset of estrus. 1-cell–stage embryos were collected 19 h after mating and were identified by the presence of sperm attached to the zona pellucida and a visible second polar body within the perivitelline space. 2-cell–stage embryos were collected 40 h after mating. 4-cell–stage embryos were collected 60–70 h after mating. Eight-cell–stage embryos were collected 3.5–4 days after mating, and morula-stage embryos were collected 4.5–5 days after mating. Single-blastomere isolation was subsequently performed. Embryos at each developmental stage were treated with an acidic Tyrode’s solution to remove the zona pellucida. After zona removal, embryos were washed twice with calcium- and magnesium-free DPBS and then incubated in Accutase solution for approximately 30 min. Individual blastomeres were gently dissociated by repeated aspiration and expulsion using a mouth-controlled glass pipette with a diameter appropriate for single-cell isolation. Isolated blastomeres were washed three times in DPBS and transferred into microcentrifuge tubes containing lysis buffer. Samples were immediately frozen and stored at − 80 °C until further processing. Single-cell transcriptome amplification and RNA sequencing Single-cell samples were amplified using the Smart-seq2 protocol. First-strand cDNA synthesis was performed in a reaction mixture containing reverse transcriptase, reaction buffer, template-switching oligonucleotide primers, and oligo(dT) primers with a common sequence. Second-strand synthesis and amplification were achieved by adding PCR reagents and ISPCR primers with a common sequence to the first-strand cDNA products. The resulting cDNA was fragmented to an average length of approximately 300 bp using a Bioruptor® sonication system (Diagenode Inc., Denville, NJ, USA). Fragmented cDNA then underwent end repair, A-tailing, sequencing adapter ligation, and PCR amplification. PCR products with fragment sizes ranging from 350 to 450 bp were size-selected and purified. Library concentrations were quantified using a Qubit 2.0 fluorometer (Life Technologies, Carlsbad, CA, USA), and library integrity was assessed with an Agilent 2100 Bioanalyzer using the High-Sensitivity DNA Assay Kit (Agilent Technologies, Santa Clara, CA, USA). Final libraries were sequenced on an Illumina HiSeq X Ten platform using 150 bp paired-end reads. Single-cell RNA sequencing data acquisition and processing Raw sequencing data from a total of 60 individual human embryonic cells spanning the oocyte to morula stages (derived from in vitro fertilization) were retrieved from the NCBI gene expression omnibus (GEO; http://www.ncbi.nlm.nih.gov/geo/ ; accessed on 19 July 2023) under accession number GSE36552 . Pig embryo sequencing data were generated in our laboratory and stored on a local server. Detailed metadata for all samples are provided in Supplementary Table S1. Raw reads were first processed using fastp (version 0.23.2) [ 72 ] to remove adapter sequences and filter low-quality reads. The resulting high-quality reads were then aligned to the human and pig reference genomes (GRCh38.107 and Sscrofa11.1, respectively) using the STAR aligner (version 2.7.10a) [ 73 ]. Gene expression levels, including raw counts, TPM, and fragments per kilobase of transcript per million mapped reads, were quantified using RSEM (version 1.3.3) [ 74 ] with default parameters. Data integration, sample selection and preliminary analysis Orthologous genes between humans and pigs were retrieved from Ensembl ( http://www.ensembldb.ensembl.org/ (accessed on 19 July 2023). To minimize potential gene omissions, orthology assignments were further verified using BLASTp, with canonical human and pig protein FASTA sequences obtained from UniProt ( https://www.uniprot.org/ ). Homology was defined using a minimum sequence coverage of 75% and an E-value threshold < 0.05. All gene identifiers were standardized to human gene symbols based on annotations from Ensembl and NCBI ( https://www.ncbi.nlm.nih.gov/ ). TPM values were log-transformed (log-e (TPM + 1)) and analyzed using Seurat (version 5.0.1) [ 75 ] for dimensionality reduction. A mild application of Seurat’s RPCA integration strategy was performed with caution, primarily to provide a reference for cross-species developmental stage alignment rather than full data integration. Pseudotime trajectories for each species were inferred independently using Monocle3 (version 1.3.1) [ 76 ] based on the reduced dimensionality space, and three-dimensional visualizations were generated using the Plotly package (version 4.10.4). Inter-sample correlations were assessed using Pearson correlation coefficients calculated from the top 4000 highly variable genes identified by Seurat for each species. To further refine the analysis, samples exhibiting low intragroup homogeneity or ambiguous positioning between adjacent developmental stages were excluded from all subsequent analyses. Stage-specific differential expression markers were identified using Seurat’s FindAllMarkers function. Stable upregulated markers were defined by the following criteria: log2FC > 1, adjusted p -value < 0.05, and consistent expression across all cells within the corresponding developmental stage. In contrast, downregulated markers were identified using log2FC < −1 with a more permissive significance threshold ( p -value < 0.05), a strategy adopted to accommodate limited sample sizes at certain embryonic stages while maintaining consistency of identified markers across stages. Differential expression analysis in adjacent stages within species and in oocytes across species For intra-species analyses, DESeq2 (version 1.38.3) [ 77 ] was used to identify DEGs between adjacent developmental stages based on orthologous gene count data. DESeq2 applies a negative binomial generalized linear model and performs significance testing using the Wald test, with p -values adjusted by the Benjamini–Hochberg FDR method. Genes were defined as DEGs using thresholds of |log2 (FC)| > 1 and adjusted p -value < 0.05. In addition, genes with minimal expression change (|log2FC| < 0.25) were considered conserved, serving as a reference for expression stability across stages. For interspecies oocyte comparisons involving limited sample sizes, we employed the ROC-based differential expression method implemented in Seurat on logₑ-transformed TPM data. Significant DEGs were identified using stringent criteria (|AUC| = 1 and |log2FC| > 1.5), which accommodates batch effects and cross-species variability by prioritizing robust, directionally consistent expression differences. Analysis of intra-species differences variation between 4-cell and 8-cell stages To identify genes exhibiting reduced or increased interspecies divergence between the 4-cell and 8-cell stages, we applied the ROC method implemented in Seurat to log-transformed TPM expression data from human and pig embryos. Genes showing significant interspecies differences at the 4-cell stage were first selected using stringent criteria (|AUC| = 1 and |log2FC| > 1.5). For each selected gene, we calculated the ratio of the avg_diff value (derived from the ROC analysis) at the 4-cell stage to that at the 8-cell stage. Genes with a ratio ≤ 0.5 were defined as interspecies difference-reducing genes, indicating that the expression divergence decreased to less than half of its original magnitude during the 4-cell-to-8-cell transition. In contrast, genes with a ratio ≥ 2 were classified as interspecies difference-expanding genes, reflecting an increase in divergence to more than twice the initial level. Enrichment analysis Gene enrichment analyses, including GO and KEGG pathway analyses, were conducted using the clusterProfiler package (version 4.6.2) [ 78 ]. KEGG pathway enrichment was assessed using a hypergeometric test, with gene annotation supported by the org.Hs.eg.db package (version 3.16.0) to ensure annotation accuracy. Enriched terms were considered significant at an adjusted p -value < 0.05 following Benjamini–Hochberg FDR correction. To reduce redundancy and improve interpretability of GO enrichment results, the simplifyEnrichment package (version 1.8.0) [ 79 ] was applied. Visualization of KEGG pathways (accessed on 6 August 2023) was performed using the pathview package (version 1.42.0) [ 80 ], with pathway nodes colored according to normalized log₂ FC (log2FC) values. Analysis of TFs and factors’ binding sites The list of human TFs and co-TFs was obtained from AnimalTFDB ( http://bioinfo.life.hust.edu.cn/AnimalTFDB4/#/version , accessed on 10 September 2023) and used as the reference for subsequent analyses. Soft clustering of TF expression patterns was performed using ClusterGVis (version 1.42.0). TF activity was inferred using SCENIC (version 1.3.1) [ 81 ] with human annotation files. To evaluate correlations in TF expression patterns between humans and pigs, the logₑ(TPM + 1) expression matrix was transformed as follows: the 4-cell stage was excluded from the human dataset, whereas both the 4-cell and 8-cell stages were retained for pigs, resulting in two pig expression matrices. Pearson correlation coefficients were calculated between the human matrix and each pig matrix, and the higher coefficient was selected as the final correlation value for each TF. Orthology-based conservation scoring was further conducted using the MetaPhOrs repository ( http://orthology.phylomedb.org/version 2.5; accessed on 10 September 2023) [ 82 ] For conservation analysis of highly activated TFs, a recommended consistency score threshold of > 0.5 was applied. Putative binding sites of TFs and other regulatory factors were retrieved from the GTRD database ( http://gtrd20-06.biouml.org/ ) using the following criteria: TF binding site location within promoter regions (−1000 to + 100 bp relative to the TSS) and organism restricted to human ( Homo sapiens ). Protein–protein interaction analysis General PPIs were investigated using the STRING database ( https://string-db.org/ ; accessed on 10 September 2023) and visualized with Cytoscape (version 3.10.1) [ 83 ]. Protein sequences were obtained from UniProt ( https://www.uniprot.org/ ), and PPIs were further predicted using HDOCK ( http://hdock.phys.hust.edu.cn/ ) [ 84 ]. For each protein pair, only the results with the highest confidence score were retained. The confidence score is an HDOCK-defined metric representing the likelihood of interaction, ranging from 0 to 1, with values closer to 1 indicating higher interaction probability. The UniProt accession numbers for all protein sequences used in this study are provided in Supplementary Table S2. In addition, the potential probability of direct PPIs was evaluated using AlphaFold 3 ( https://golgi.sandbox.google.com/ ). Experimentally supported PPI information from BioGRID ( https://thebiogrid.org/ ) was incorporated as supplementary evidence. Analysis of MYC and BCLAF1 binding probabilities near TS sin specific gene sets Using human MYC ChIP-seq data from ENCODE ( https://www.encodeproject.org/ ), including narrowPeak files ENCSR000DLZ (K562) and ENCSR000DKU (GM12878), and BCLAF1 ChIP-seq data ENCSR000BKH (K562) and ENCSR000BJZ (GM12878), we calculated the distribution of binding probabilities for these two factors within ±1,000 bp of TSSs. The analyzed gene set comprised 107 genes, defined as the intersection of BCLAF1 and MYC regulons among conserved major ZGA genes across species. Binding probabilities were computed by weighting ChIP-seq signal intensities (narrowPeak score column divided by 1000) and averaging across the number of TSSs per gene. Analyses were performed using the R packages Biostrings (v2.66.0) and ChIPseeker (v1.38.0) [ 85 ]. TSS coordinates were obtained from BioMart (v2.58.2) [ 86 ]. For genes with multiple TSSs located within 100 bp of each other, these sites were merged and represented by their midpoint as the final TSS. Due to the lack of experimental ChIP-seq data in pigs, motif-based analyses were conducted using SCENIC-inferred best motifs: “dbcorrdb__POLR2A__ENCSR000EAU_1__m1” for MYC and “dbcorrdb__BCLAF1__ENCSR000BKH_1__m1” for BCLAF1. The distribution of predicted binding sites within ±500 bp of TSSs was calculated using the R package universalmotif (v1.20.0), applying the following thresholds: score.pct > 0.75, p value < 0.05, and FDR < 0.2. Position weight matrix files were obtained from Factorbook ( https://www.factorbook.org/ ). Porcine TSS coordinates were retrieved from BioMart, and ±500 bp genomic sequences surrounding TSSs were extracted using Biostrings in conjunction with the BSgenome.Sscrofa.UCSC.susScr11 reference genome (v1.4.2). Definition of chromatin active regions and relationship between two major ZGA phases in pigs Chromatin active regions during the two major ZGA phases in pigs were defined based on the chromosomal locations of independently identified major ZGA gene subsets. Genomic positional information was retrieved from BioMart, and gene annotations were based on the BSgenome.Sscrofa.UCSC.susScr11 reference genome. Chromosomal distributions were visualized using karyoploteR (v1.28.0). Each chromosome was partitioned into equal-sized bins, and the number of active genes within each bin was quantified for both major ZGA phases. The relationship between chromatin activity patterns across the two phases was evaluated using Pearson correlation coefficients and the 1 − Jaccard index based on gene count distributions. Genomic coverage was calculated as the proportion of the chromosome length occupied by active regions. Bins with a Jaccard index ≤ 0.5 were classified as low-similarity regions. RNA binding information The RNA binding information for BCLAF1 was retrieved from the POSTAR3 database ( http://111.198.139.65/index.html ) [ 87 ]. A threshold of 200 binding site records was set to select the genes for analysis. Electronic supplementary material Below is the link to the electronic supplementary material. Supplementary Material 1 (227.2MB, docx) Acknowledgements This work was supported by National Natural Science Foundation of China (32272885), the Key Project of Natural Science Foundation of Heilongjiang Province of China (ZD2022C004). Author contributions Jianlin Fan, Qinjian Li, Zhongyu Yuan, Conceptualization, Data curation, Investigation, Methodology, Writing – original draft; Chaoqian Jiang, Investigation, Methodology; Xiaokang Xu, Xueqing Liu, Methodology; Yanshuang Mu, Writing – review and editing; Zhonghua Liu, Conceptualization, Supervision, Funding acquisition, Writing – review and editing. Data availability This study utilized single-cell transcriptome data from human and pig oocytes through to morula stages. The in vivo pig embryo data were generated in-house and are available in the NCBI GEO database (Accession: GSE139512 ). Human embryo data was also obtained from the NCBI GEO database (Accession: GSE36552 ). Declarations Competing interest The authors declare that no competing interests exist. Footnotes Publisher’s Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Jianlin Fan, Qinjian Li and Zhongyu Yuan contributed equally to this work. Contributor Information Zhonghua Liu, Email: [email protected]. Yanshuang Mu, Email: [email protected]. References 1. Schulz KN, Harrison MM. Mechanisms regulating zygotic genome activation. Nat Rev Genet. 2019;20(4):221–34. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 2. Lee MT, Bonneau AR, Giraldez AJ. Zygotic genome activation during the maternal-to-zygotic transition. Annu Rev Cell Dev Biol. 2014;30(1):581–613. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 3. Tadros W, Lipshitz HD. The maternal-to-zygotic transition: a play in two acts. Development. 2009;136(18):3033–42. [ DOI ] [ PubMed ] [ Google Scholar ] 4. Sha QQ, et al. Dynamics and clinical relevance of maternal mRNA clearance during the oocyte-to-embryo transition in humans. Nat Commun. 2020;11(1):4917. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 5. Sha QQ, et al. Characterization of zygotic genome activation-dependent maternal mRNA clearance in mouse. Nucleic Acids Res. 2020;48(2):879–94. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 6. Boroviak T, et al. Single cell transcriptome analysis of human, marmoset and mouse embryos reveals common and divergent features of preimplantation development. Development. 2018;145(21). [ DOI ] [ PMC free article ] [ PubMed ] 7. Liu TB, et al. Cross-species single-cell transcriptomic analysis reveals pre-gastrulation developmental differences among pigs, monkeys, and humans. Cell Discov. 2021;7(1). [ DOI ] [ PMC free article ] [ PubMed ] 8. !!! INVALID CITATION! [8]. 9. Vastenhouw NL, Cao WX, Lipshitz HD. The maternal-to-zygotic transition revisited. Development. 2019;146(11):dev161471. [ DOI ] [ PubMed ] 10. Magnani L, Johnson CM, Cabot RA. Expression of eukaryotic elongation initiation factor 1A differentially marks zygotic genome activation in biparental and parthenogenetic porcine embryos and correlates with in vitro developmental potential. Reprod Fertil Dev. 2008;20(7):818–25. [ DOI ] [ PubMed ] [ Google Scholar ] 11. Zhou W, et al. Nuclear accumulation of pyruvate dehydrogenase alpha 1 promotes histone acetylation and is essential for zygotic genome activation in porcine embryos. Biochim Et Biophys Acta (BBA)-Mol Cell Res. 2020;1867(4):118648. [ DOI ] [ PubMed ] [ Google Scholar ] 12. Du ZQ, et al. Single cell RNA-seq reveals genes vital to in vitro fertilized embryos and parthenotes in pigs. Sci Rep. 2021;11(1):14393. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 13. Wang XY, et al. Single-cell transcriptomic and cross-species comparison analyses reveal distinct molecular changes of porcine testes during puberty. Commun Biol. 2024;7(1). [ DOI ] [ PMC free article ] [ PubMed ] 14. Tian Y, et al. The transcriptomic signature of donkey ovarian tissue revealed by cross-species comparative analysis at single-cell resolution. Animals. 2025;15(12). [ DOI ] [ PMC free article ] [ PubMed ] 15. Meurens F, et al. The pig: a model for human infectious diseases. Trends Microbiol. 2012;20(1):50–57. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 16. Lunney JK. Advances in swine biomedical model genomics. Int J Biol Sci. 2007;3(3):179. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 17. Virant-Klun I, Krijgsveld J. Proteomes of animal oocytes: what can we learn for human oocytes in the in vitro fertilization programme? Biomed Res Int. 2014;2014(1):856907. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 18. Lucas CG, et al. Applications of omics and nanotechnology to improve pig embryo production in vitro. Mol Reprod Devel. 2019;86(11):1531–47. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 19. Redel BK, Spate LD, Prather RS. In vitro maturation, fertilization, and culture of pig oocytes and embryos. Comparative embryo culture: methods and protocols, 2019. 93–103. [ DOI ] [ PubMed ] 20. Lee JW, Tian XC, Yang X. Failure of male pronucleus formation is the major cause of lack of fertilization and embryo development in pig oocytes subjected to intracytoplasmic sperm injection. Biol Reproduction. 2003;68(4):1341–47. [ DOI ] [ PubMed ] [ Google Scholar ] 21. Lunney JK, et al. Importance of the pig as a human biomedical model. Sci Transl Med. 2021;13(621):eabd 5758. [ DOI ] [ PubMed ] 22. Yan L, et al. Single-cell RNA-Seq profiling of human preimplantation embryos and embryonic stem cells. Nat Struct Mol Biol. 2013;20(9):1131–39. [ DOI ] [ PubMed ] [ Google Scholar ] 23. Zheng CB, Ballard EB, Wu J. The road to generating transplantable organs: from blastocyst complementation to interspecies chimeras. Development. 2021;148(12). [ DOI ] [ PMC free article ] [ PubMed ] 24. McPherson JP, et al. Essential role for Bclaf1 in lung development and immune system function. Cell Death Differ. 2009;16(2):331–39. [ DOI ] [ PubMed ] [ Google Scholar ] 25. Sarras H, Azami SA, McPherson JP. In search of a function for BCLAF1. ScientificWorldjournal. 2010;10(1):1450–61. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 26. Zhou X, et al. Heat shock protein 90α-dependent B-cell-2-associated transcription factor 1 promotes hepatocellular carcinoma proliferation by regulating MYC proto-oncogene c-MYC mRNA stability. Hepatology. 2019;69(4):1564–81. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 27. Lee SH, et al. The splicing factor SRSF1 modulates pattern formation by inhibiting transcription of tissue specific genes during embryogenesis. Biochem Bioph Res Co. 2016;477(4):1011–16. [ DOI ] [ PubMed ] [ Google Scholar ] 28. Deml B, et al. EFTUD2 deficiency in vertebrates: identification of a novel human mutation and generation of a zebrafish model. Birth Defects Res Clin Mol Teratol. 2015;103(7):630–40. https://onlinelibrary.wiley.com/doi/10.1002/bdra.23397 . [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 29. Young JC, et al. Distinct effects of importin α2 and α4 on Oct3/4 localization and expression in mouse embryonic stem cells. FASEB J. 2011;25(11):3958–65. [ DOI ] [ PubMed ] [ Google Scholar ] 30. Zhou C, et al. The transcription factor GABPA is a master regulator of naive pluripotency. Nat Cell Biol. 2025;27(1):48–58. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 31. Yuan S, et al. Human zygotic genome activation is initiated from paternal genome. Cell Discov. 2023;9(1):13. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 32. Wang WP, et al. The EP300, KDM5A, KDM6A and KDM6B chromatin regulators cooperate with KLF4 in the transcriptional activation of POU5F1. PLoS One. 2012;7(12):e52556. [ DOI ] [ PMC free article ] [ PubMed ] 33. Ebrahimi A, et al. Bromodomain inhibition of the coactivators CBP/EP300 facilitate cellular reprogramming. Nat Chem Biol. 2019;15(5):519–28. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 34. Fasano C, et al. The chromatin remodeling factors EP300 and TRRAP are novel SMYD3 interactors involved in the emerging ‘nonmutational epigenetic reprogramming’ cancer hallmark. Comput Struct Biotechnol J. 2023;21:5240–48. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 35. Tian Q, et al. Chromatin modifier EP400 regulates oocyte quality and zygotic genome activation in mice. Adv Sci. 2024;11(20):2308018. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 36. Gehlen-Breitbach S, et al. The Tip60/Ep400 chromatin remodeling complex impacts basic cellular functions in cranial neural crest-derived tissue during early orofacial development. Int J Oral Sci. 2023;15(1):16. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 37. Wu M, et al. The MTA family proteins as novel histone H3 binding proteins. Cell & Biosci. 2013;3(1):1. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 38. Zhou C, et al. MTA2 promotes gastric cancer cells invasion and is transcriptionally regulated by Sp1. Mol Cancer. 2013;12(1):102. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 39. Li X, et al. LINE-1 transcription activates long-range gene expression. Nat Genet. 2024;56(7):1494–502. [ DOI ] [ PubMed ] [ Google Scholar ] 40. Acharya D, et al. KAT-Independent gene regulation by Tip60 promotes ESC self-renewal but not pluripotency. Cell Rep. 2017;19(4):671–79. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 41. Acharya D, et al. TIP55, a splice isoform of the KAT5 acetyltransferase, is essential for developmental gene regulation and organogenesis. Sci Rep. 2018;8(1):14908. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 42. Zhang J, et al. Super RNA pol II domains enhance minor ZGA through 3D interaction to ensure the integrity of major transcriptional waves in late-ZGA mammals. Cell Genom. 2025;5(6):100856. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 43. Dang Y, et al. The lysine deacetylase activity of histone deacetylases 1 and 2 is required to safeguard zygotic genome activation in mice and cattle. Development. 2022;149(11):dev200854. [ DOI ] [ PubMed ] 44. Shao R, et al. The impact of selective HDAC inhibitors on the transcriptome of early mouse embryos. BMC Genomics. 2024;25(1):143. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 45. Bu G, et al. Coordination of zygotic genome activation entry and exit by H3K4me3 and H3K27me3 in porcine early embryos. Genome Res. 2022;32(8):1487–501. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 46. Kubinyecz ON, et al. Maternal SMARCA5 is required for major ZGA in mouse embryos bioRxiv, 2023. 2023.12.05.570276. 47. Wei J, et al. Temporospatial hierarchy and allele-specific expression of zygotic genome activation revealed by distant interspecific urochordate hybrids. Nat Commun. 2024;15(1):2395. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 48. Chen J, et al. CircMYO10 promotes osteosarcoma progression by regulating miR-370-3p/RUVBL1 axis to enhance the transcriptional activity of β-catenin/LEF1 complex via effects on chromatin remodeling. Mol Cancer. 2019;18(1):150. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 49. Shy BR, et al. Regulation of Tcf7l1 DNA binding and protein stability as principal mechanisms of wnt/β-catenin signaling. Cell Rep. 2013;4(1):1–9. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 50. Jiang H, et al. ID proteins promote the survival and primed-to-naive transition of human embryonic stem cells through TCF3-mediated transcription. Cell Death Disease. 2022;13(6):549. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 51. Weise A, et al. Alternative splicing of Tcf7l2 transcripts generates protein variants with differential promoter-binding and transcriptional activation properties at Wnt/β-catenin targets. Nucleic Acids Res. 2010;38(6):1964–81. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 52. Katoh M. Multi‑layered prevention and treatment of chronic inflammation, organ fibrosis and cancer associated with canonical WNT/β‑catenin signaling activation. Int J Mol Med. 2018;42(2):713–25. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 53. Hrckulak D, et al. TCF/LEF transcription factors: an update from the internet Resources. Cancers. 2016;8(7):70. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 54. Han J, et al. Expression and regulation of inhibitor of DNA binding proteins ID1, ID2, ID3, and ID4 at the maternal-conceptus interface in pigs. Theriogenology. 2018;108:46–55. [ DOI ] [ PubMed ] [ Google Scholar ] 55. Hong SH, et al. ID1 and ID3 represent conserved negative regulators of human embryonic and induced pluripotent stem cell hematopoiesis. J Cell Sci. 2011;124(9):1445–52. [ DOI ] [ PubMed ] [ Google Scholar ] 56. Liu Z, et al. Mutations in the RNA splicing factor SF3B1 promote Tumorigenesis through MYC stabilization. Cancer Discov. 2020;10(6):806–21. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 57. Odero MD, et al. A novel gene, MDS2, is fused to ETV6/TEL in a t(1;12)(p36.1;p13) in a patient with myelodysplastic syndrome. Genes Chromosomes Cancer. 2002;35(1):11–19. [ DOI ] [ PubMed ] [ Google Scholar ] 58. Bevilacqua A, Fiorenza MT, Mangia F. A developmentally regulated GAGA box-binding factor and Sp1 are required for transcription of the hsp70.1 gene at the onset of mouse zygotic genome activation. Development. 2000;127(7):1541–51. [ DOI ] [ PubMed ] [ Google Scholar ] 59. Abe A, et al. ETV6-LPXN fusion transcript generated by t(11;12)(q12.1;p13) in a patient with relapsing acute myeloid leukemia with NUP98-HOXA9. Genes Chromosomes Cancer. 2016;55(3):242–50. [ DOI ] [ PubMed ] [ Google Scholar ] 60. Sturmey RG, et al. Role of fatty acids in energy provision during oocyte maturation and early embryo development. Reprod Domest Anim. 2009;44(3):50–58. [ DOI ] [ PubMed ] [ Google Scholar ] 61. Wang X, et al. Abnormal expression of HADH, an enzyme of fatty acid oxidation, affects tumor development and prognosis (review). Mol Med Rep. 2022;26(6):355. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 62. Arêdes DS, et al. Silencing of mitochondrial trifunctional protein a subunit (HADHA) increases lipid stores, and reduces oviposition and flight capacity in the vector Insect rhodnius prolixus. Front Insect Sci. 2022;2:885172. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 63. Kajdasz A, et al. Lipid Stores and lipid metabolism associated gene expression in porcine and bovine parthenogenetic embryos revealed by fluorescent staining and RNA-seq. Int J Mol Sci. 2020;21(18):6488. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 64. Mota-Rojas D, et al. Foetal and neonatal energy metabolism in pigs and humans: a review. Veterinární medicína. 2011;56(5):215–25. [ Google Scholar ] 65. Umathum V, et al. Distribution of ferritin complex in the adult brain and altered composition in neuroferritinopathy due to a novel variant in the ferritin heavy chain gene FTH1 (c.409_410del; p.H137Lfs*4). Brain Pathol. 2024;34(1):e13176. [ DOI ] [ PMC free article ] [ PubMed ] 66. Shieh JT, et al. Heterozygous nonsense variants in the ferritin heavy-chain gene FTH1 cause a neuroferritinopathy. HGG Adv. 2023;4(4):100236. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 67. Tian Y, et al. FTH1 inhibits ferroptosis through ferritinophagy in the 6-OHDA model of Parkinson’s Disease. Neurotherapeutics. 2020;17(4):1796–812. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 68. Ren L, et al. High-resolution profiles of gene expression and DNA methylation highlight mitochondrial modifications during early embryonic development. J Reprod Dev. 2017;63(3):247–61. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 69. May-Panloup P, et al. Increase of mitochondrial DNA content and transcripts in early bovine embryogenesis associated with upregulation of mtTfA and Nrf1 transcription factors. Reprod Biol Endocrinol. 2005;3(1):65. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 70. Festuccia N, et al. Nr5a2 is dispensable for zygotic genome activation but essential for morula development. Science. 2024;386(6717):eadg 7325. [ DOI ] [ PubMed ] 71. Lu T, et al. Decoding transcriptional identity in developing human sensory neurons and organoid modeling. Cell. 2024;187(26):7374–93.e28. [ DOI ] [ PubMed ] [ Google Scholar ] 72. Chen S, et al. Fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics. 2018;34(17):i884–90. [ DOI ] [ PMC free article ] [ PubMed ] 73. Dobin A, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 74. Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinf. 2011;12:1–16. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 75. Hao Y, et al. Dictionary learning for integrative, multimodal and scalable single-cell analysis. Nat Biotechnol. 2024;42(2):293–304. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 76. Cao J, et al. The single-cell transcriptional landscape of mammalian organogenesis. Nature. 2019;566(7745):496–502. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 77. Love MI, Huber W, Anders S. Moderated estimation of Fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:1–21. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 78. Yu G, et al. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics. 2012;16(5):284–87. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 79. Gu Z, Hübschmann D. simplifyEnrichment: a Bioconductor package for clustering and visualizing functional enrichment results. Genomics, Proteomics Bioinf. 2023;21(1):190–202. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 80. Luo W, Brouwer C. Pathview: an R/Bioconductor package for pathway-based data integration and visualization. Bioinformatics. 2013;29(14):1830–31. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 81. Aibar S, et al. SCENIC: single-cell regulatory network inference and clustering. Nat Methods. 2017;14(11):1083–86. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 82. Chorostecki U, et al. MetaPhOrs 2.0: integrative, phylogeny-based inference of orthology and paralogy across the tree of life. Nucleic Acids Res. 2020;48(W1):W553–57. [ DOI ] [ PMC free article ] [ PubMed ] 83. Shannon P, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13(11):2498–504. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 84. Yan Y, et al. The HDOCK server for integrated protein-protein docking. Nat Protoc. 2020;15(5):1829–52. [ DOI ] [ PubMed ] [ Google Scholar ] 85. Wang Q, et al. Exploring epigenomic datasets by ChIPseeker. Curr Protocol. 2022;2(10):e585. [ DOI ] [ PubMed ] 86. Durinck S, et al. Mapping identifiers for the integration of genomic datasets with the R/Bioconductor package biomaRt. Nat Protoc. 2009;4(8):1184–91. [ DOI ] [ PMC free article ] [ PubMed ] [ Google Scholar ] 87. Zhao W, et al. POSTAR3: an updated platform for exploring post-transcriptional regulation coordinated by RNA-binding proteins. Nucleic Acids Res. 2022;50(D1):D287–94. [ DOI ] [ PMC free article ] [ PubMed ] Associated Data This section collects any data citations, data availability statements, or supplementary materials included in this article. Supplementary Materials Supplementary Material 1 (227.2MB, docx) Data Availability Statement This study utilized single-cell transcriptome data from human and pig oocytes through to morula stages. The in vivo pig embryo data were generated in-house and are available in the NCBI GEO database (Accession: GSE139512 ). Human embryo data was also obtained from the NCBI GEO database (Accession: GSE36552 ). Articles from Journal of Translational Medicine are provided here courtesy of BMC ACTIONS View on publisher site PDF (9.3 MB) Cite Collections Permalink PERMALINK Copy RESOURCES Similar articles Cited by other articles Links to NCBI Databases Cite Copy Download .nbib .nbib Format: AMA APA MLA NLM Add to Collections Create a new collection Add to an existing collection Name your collection * Choose a collection Unable to load your collection due to an error Please try again Add Cancel Follow NCBI NCBI on X (formerly known as Twitter) NCBI on Facebook NCBI on LinkedIn NCBI on GitHub NCBI RSS feed Connect with NLM NLM on X (formerly known as Twitter) NLM on Facebook NLM on YouTube National Library of Medicine 8600 Rockville Pike Bethesda, MD 20894 Web Policies FOIA HHS Vulnerability Disclosure Help Accessibility Careers NLM NIH HHS USA.gov Back to Top