RESULTS

Transcriptome-wide mapping of m6A in bovine oocytes and preimplantation embryos by SAC-seq

SAC-seq enables single-nucleotide-resolution mapping of m6A sites,21 overcoming the limitations of antibody-based MeRIP-seq by detecting both canonical (DRACH) and non-canonical motifs without sequence bias. Here, we employed SAC-seq to simultaneously analyze the transcriptome and m6A epitranscriptome of bovine oocytes at the germinal vesicle (GV) and meta-phase II (MII) stages, as well as preimplantation embryos at the 2-cell, 4-cell, 8-cell, 16-cell, morula, and blastocyst stages (Figure 1A). For each sample, 200 oocytes or embryos were used, and experiments were performed twice. On average, we obtained approximately 39.3 ± 18.8 million uniquely mapped reads per SAC-seq sample (m6A epitranscriptome) and 32.7 ± 9.9 million reads per input control (transcriptome). Principal component analysis (PCA) of both the transcriptome and m6A epitranscriptome data showed high consistency between biological replicates and revealed stage-specific patterns of gene expression and m6A profiles (Figure 1A).

Motif analysis revealed strong enrichment of the GA motif across developmental stages (Figure 1B), consistent with previously reported SAC-seq performance,21 supporting its utility for identifying m6A sites in oocytes and early embryos. We next grouped m6A sites based on whether they were located within the canonical DRACH motif. On average, 33.5% ± 2.4% of identified m6A sites fell within the DRACH motif, ranging from 29.9% at the 4-cell stage to 36.6% at the MII stage (Data S1). Metagene analysis revealed that DRACH-associated m6A sites were strongly enriched near stop codons (Figure 1C), consistent with observations in mouse12 and humans.15 In contrast, m6A deposition in non-DRACH motifs was enriched within coding DNA sequences (CDSs) (Figure 1D), highlighting expanded sequence diversity and functional potential beyond canonical motif regions. To validate these observations, independent low-input eTAM-seq analysis confirmed m6A sites in both DRACH and non-DRACH sequence contexts for selected transcripts, supporting the authenticity of non-canonical m6A sites (Figures S1A and S1B). Motif frequency analysis across different read coverage thresholds showed consistent proportions of DRACH and RRACH motifs, indicating that enrichment of non-DRACH sites is not driven by sequencing depth or calling bias (Figures S1C and S1D).

The total number of identified m6A sites across preimplantation stages ranged from 60,678 (GV), 42,825 (MII), 8,598 (2-cell), 27,695 (4-cell), 10,683 (8-cell), 10,263 (16-cell), 8,593 (morula), to 18,320 (blastocyst) (Figure 1E). These sites correspond to 6,779 genes in GV, 5,879 in MII, 2,263 in 2-cell, 4,353 in 4-cell, 2,634 in 8-cell, 2,379 in 16-cell, 2,037 in morula, and 3,156 in blastocyst stages (Figure 1F). Notably, maternal transcripts in oocytes were highly m6A-modified compared to later stages (Figures 1E and 1F). These modifications sharply declined after fertilization at the 2-cell stage, indicating a rapid epitranscriptomic switch coincident with MZT. This was followed by two waves of m6A re-establishments, one at the 4-cell stage, just prior to bovine major ZGA,17 suggesting a potential role in minor ZGA, and another at the blastocyst stage, coinciding with lineage specification. Sankey analysis tracking transcripts with or without m6A across preimplantation development showed that most m6A-marked transcripts in GV oocytes were retained in MII oocytes, but were largely lost at the 2-cell stage (Figure 1G). A subset of maternal transcripts without m6A became methylated at later embryonic stages, particularly in blastocysts (Figure 1G). Overall, these data indicate that maternal transcripts are extensively modified by m6A, likely contributing to stability and timely clearance, whereas demethylation of the majority of genes at 8-cells and re-methylation at blastocysts support ZGA and lineage differentiation, respectively. Together, these findings highlight dynamic and tightly regulated m6A deposition as a key layer of post-transcriptional control during bovine early embryogenesis.

m6A epitranscriptomic regulation in differentially expressed gene clusters

To explore the relationship between m6A and gene expression, we analyzed differentially expressed genes (DEGs) across bovine oocytes and preimplantation embryos. A significantly higher proportion of DEGs were m6A-modified compared to non-DEGs, with strong enrichment in protein-coding genes (Figure 2A), suggesting preferential marking of dynamically regulated transcripts.

k-Means clustering of protein-coding DEGs revealed six distinct expression patterns during preimplantation development (Figure 2B). The fraction of m6A-modified genes varied across clusters, with clusters 1, 2, and 6 showing the highest proportions (Figure 2C), indicating functional diversity in stage-specific regulation. Similar proportions were observed for DRACH and non-DRACH sites (Figure 2C). Integration of Figures 2B and 2C revealed that oocyte- and early cleavage-enriched genes are highly m6A modified, followed by a sharp decrease at the 8-cell stage and a recovery at the morula and blastocyst stages.

Gene Ontology (GO) analysis of m6A-marked DEGs revealed a sequential functional enrichment during developmental progression (Figure 2D). Cluster 1 genes, highly expressed in oocytes, were enriched in dephosphorylation, carboxylic acid biosynthetic process, and protein localization. Cluster 2 genes, active during early cleavage stages, were associated with oogenesis, cell cycle progression, and gamete generation. Clusters 3–6 showed distinct enrichments: cluster 3 in cell-cell signaling and development, cluster 4 in transcriptional regulation and response to extracellular stimuli, cluster 5 in extracellular matrix organization and differentiation, and cluster 6 in translational initiation, ribonucleoprotein complex biogenesis, and mitochondrial metabolism. In contrast, non-m6A-marked genes were more frequently associated with signaling and transcriptional regulation pathways (Figure 2D). These differences suggest that m6A selectively marks transcripts associated with specific developmental programs.

Dynamics of m6A epitranscriptomic regulation

Transcripts with high m6A density are known to be subject to post-transcriptional regulation, including increased decay, altered translation efficiency, and localization changes, depending on the context and interaction with specific m6A reader proteins.6,31,32 We therefore analyzed m6A density dynamics within DEG clusters using subcluster-level analysis (Figures 2E and S2). Distinct temporal patterns of m6A enrichment were observed across subclusters (Figure 2F). Some subclusters showed higher m6A density during early cleavage stages (Figure 2F, Cluster 1 sub1, 3, and 5, Cluster 2 sub2 and 5, and Cluster 3 sub4), whereas others exhibited increased m6A toward morula and blastocyst stages (Figure 2F, Cluster 3 sub3–5 and Cluster 6 sub3–5), indicating stage-dependent regulation. Comparison of m6A density and DRACH motif enrichment (Figures 2F and S3A) revealed heterogeneous and stage-dependent relationships rather than a consistent global correlation. In some subclusters, higher m6A density was accompanied by lower DRACH motif enrichment; however, this pattern was not uniformly observed. These findings suggest that both canonical and non-canonical targeting mechanisms may contribute to m6A deposition in a context-dependent manner. Consistent with this, DRACH-enriched subclusters also showed variable m6A density, indicating that motif presence alone does not determine methylation extent.

To control for expression bias, m6A density was normalized by transcript length and expression level (sites per kb per reads per kilobase of transcript per million mapped reads (RPKM)) and significant differences across subclusters were retained (Figure S3B). These results indicate that the observed m6A dynamics are not driven by transcript abundance. Integration with transcriptomic and translatomic data (Figure S2) further supports coordinated but stage-specific regulation during development.

Stage-specific m6A regulation of non-coding RNAs

To extend our analysis beyond protein-coding genes, we examined m6A dynamics in non-coding RNAs (ncRNAs) (Figure 3A). m6A was most prevalent in small nucleolar RNAs (snoRNAs), followed by long non-coding RNAs (lncRNAs) and small nuclear RNAs (snRNAs) (Figure 3B). While DRACH and non-DRACH motifs contributed similarly to lncRNA methylation, non-canonical motifs were significantly enriched in snoRNAs and snRNAs (Figure 3B). Distinct sequence motif enrichments were observed across ncRNA classes, with subtle differences in the flanking nucleotides of the m6A sites (Figure 3C).

Analysis of expression and m6A abundance across developmental stages revealed distinct patterns for each type of ncRNA (Figures 3D-3F). lncRNAs were highly expressed and m6A modified in GV, MII (Clusters 1 and 2), and blastocyst stages (Cluster 5) (Figure 3D), suggesting stage-specific roles. In contrast, snRNAs showed variable m6A independent of expression changes (Figure 3E). snRNAs also exhibited asynchronous relationships between expression and m6A modification in the majority of subclusters (Figure 3F), with increased m6A from 8-cell to blastocyst stages. To further examine whether ncRNA m6A dynamics were driven by canonical or non-canonical sites, we stratified lncRNA and snoRNA m6A patterns by DRACH and non-DRACH motifs. Both classes contributed to stage-resolved ncRNA dynamics, indicating that the observed changes are not restricted to canonical motif sites (Figures S4A and S4B). We would like to point out that upregulation of snoRNAs reflects their established roles in rRNA modification and ribosome biogenesis, which are critical for expanding protein synthesis capacity in rapidly dividing embryonic cells.33,34 Overall, these findings indicate heterogeneous and class-specific regulation of ncRNAs by m6A during early embryogenesis.

Ribosomal protein gene-associated m6A regulation

Given the central role of translation control during MZT,35 we performed correlative analysis of transcriptome, m6A epitranscriptome, and translatome of bovine preimplantation development, focusing on RPGs (GO:0005840). Most RPGs increased expression after the 8-cell stage (Figure 3G). Based on m6A dynamics, they were divided into two sub-clusters (sub1 and sub2) (Figure 3G). sub1 showed increased m6A after ZGA, whereas sub2 displayed stable m6A across oocytes and preimplantation embryos.

Mitochondrial RPGs formed four subclusters (sub1–sub4) (Figure 3H). Sub1 and sub2 showed high expression and m6A levels during GV and MII stages, suggesting maternal origin and early marking; however, their transcript and translation levels declined after the 2-cell stage and remained relatively low at later stages. In contrast, sub3 and sub4 showed increased expression and translation from the 8-cell stage onward. Notably, these subclusters did not show prevalent m6A modifications, indicating uncoupling between m6A and expression or translation efficiency.

We further assessed whether specific m6A sequence contexts contribute to differential regulation. Motif analysis revealed a higher fraction of canonical DRACH-associated sites in sub2 compared with sub1 in cytoplasmic RPGs (Figure 3I), and similarly in sub3 and sub4 compared to sub1 and sub2 for mitochondrial RPGs (Figure 3K). This difference suggests that these transcripts may be preferentially targeted by the core m6A methyltransferase complex, which typically recognizes motif-dependent sites. Indeed, sequence logo analysis further confirmed these observations, revealing subtle but distinct motif differences between subclusters in both cytoplasmic and mitochondrial ribosomal protein groups (Figures 3J and 3L). Notably, sub2 motifs closely matched the canonical DRACH consensus, but sub1 motifs displayed a different motif pattern, particularly at the +1 and +2 positions flanking the methylated adenosine (Figure 3J). These findings suggest that functionally or temporally distinct subsets of ribosomal transcripts may be associated with different m6A writer or reader complexes, adding an additional layer of epitranscriptomic specificity during early embryonic development.

Together, these results indicate that m6A dynamics in RPGs may contribute to post-transcriptional regulation in a context-dependent manner, rather than directly determining transcript or translation levels. To further investigate functional relevance, we selected RPL12 for mechanistic study.

Site-specific m6A modification of RPL12 is required for ZGA and blastocyst formation

As global m6A profiling revealed distinct regulatory patterns among RPGs during preimplantation development (Figures 3G-3L), we next investigated whether a single site-specific m6A could have functional consequences for embryonic development. We selected RPL12 for detailed analysis (Figure 4A), as it contains a highly conserved, high-confidence m6A site at the 148th adenosine (148th A) within the coding sequence, with increased methylation from 8-cell to blastocyst stages (Figure 4B), coinciding with major ZGA and lineage specification. RPL12 encodes a component of the 60S ribosomal subunit involved in protein synthesis, and previous studies have linked ribosome function to embryonic development.36-38 However, its role in ZGA, and specifically the function of this m6A site (m6A-RPL12), has not been previously characterized.

To investigate this, we used CRISPR-based prime editing (PEMax)39-41 to precisely mutate the methylated adenosine to guanosine (A148G) in zygotes (see STAR Methods), thereby preventing m6A deposition at this site (Figure 4C). Amplicon sequencing confirmed efficient and precise editing with minimal off-target effects (Figures 4D and 4E). Embryos carrying the A148G mutation developed normally to the 2- to 8-cell stage but exhibited impaired development beyond 8-cell stage and a marked reduction in blastocyst formation compared with controls (Figures 4F-4I). These results indicate that loss of this single m6A site compromises ZGA and blastocyst formation.

To confirm that the observed phenotype was specifically due to loss of m6A rather than the nucleotide substitution itself, we performed site-specific m6A demethylation using dCas13Rx-ALKBH5 (m6A demethylase) delivered into zygotes (see STAR Methods) to selectively erase m6A at the same 148th A site in RPL12 mRNA (Figure 4J). This demethylation (ALKBH5-WT) treatment also impaired embryonic development, although the phenotype was less severe than that observed with A148G mutation, likely due to incomplete demethylation efficiency. In contrast, embryos injected with the catalytically inactive mutant (ALKBH5-mut) developed normally (Figures 4K and 4L). Site-specific m6A-qPCR confirmed a significant reduction in methylation at the target site following ALKBH5 treatment (Figure 4M), which was further supported by eTAM-seq analysis showing consistent partial demethylation (Figure S5).

To determine whether restoring RPL12 expression could rescue the developmental defect, we co-injected zygotes with the PEMax system targeting A148G and exogenous wild-type RPL12 mRNA lacking m6A (Figure 4N, see STAR Methods). Despite increased total RPL12 transcript levels, supplementation of unmodified mRNA did not rescue blastocyst formation or developmental progression compared with the edited controls (Figures 4O and 4P), indicating that RPL12 transcript abundance alone is insufficient and that the m6A modification itself is required for functional rescue.

Together, these results demonstrate that m6A modification of RPL12 is essential for normal embryonic development, and that disruption of this site either by base substitution or by site-specific demethylation impairs developmental progression.

Loss of m6A on RPL12 impairs protein synthesis in early embryos

To investigate the mechanism underlying the role of m6A-RPL12 in preimplantation development, we first assessed whether loss of this modification affects global protein synthesis, given the function of RPL12 in translation. We performed a global translation assay using the Click-iT L-homopropargylglycine (HPG) labeling (see STAR Methods) at the 8-cell stage (Figure 5A). Embryos carrying the A148G mutation showed markedly reduced HPG incorporation, indicating a strong decrease in nascent protein synthesis. Signal intensity was comparable to that observed in cycloheximide-treated embryos, which serve as a negative control for active translation (Figure 5B). Quantification confirmed a significant reduction in global translation following loss of m6A at RPL12 (Figure 5C).

To validate this result using an independent approach, we examined embryos subjected to dCas13Rx-ALKBH5-mediated demethylation. Consistently, ALKBH5-WT-treated embryos exhibited reduced global protein synthesis, whereas embryos treated with catalytically inactive ALKBH5-mut showed normal translation levels (Figures 5D and 5E).

Together, these findings demonstrate that loss of m6A at RPL12 impairs protein synthesis, underscoring the critical role of site-specific m6A in regulating ribosomal protein function during early embryogenesis.

m6A-RPL12 regulates post-transcriptional stability and translation of ribosomal transcripts

To examine the broader impact of m6A-RPL12 disruption on embryonic gene expression, we performed RNA sequencing (RNA-seq) analysis on control and PEMax-edited embryos at the 4- and 8-cell stages (Figure 6A, Data S1). PCA revealed clear separation between control and edited groups, as well as between developmental stages, indicating a global transcriptomic shift associated with loss of m6A-RPL12 (Figure 6B). Notably, m6A-RPL12 disruption reduced the number of DEGs between the 4- and 8-cell transition (Figure 6C), suggesting impaired activation of ZGA-associated transcriptional programs.

Given the observed reduction in protein synthesis (Figures 5B and S6A-S6D), we focused on genes involved in translation, including ribosomal proteins (GO:0005840) and translation initiation and elongation factors (GO:0006412). In control embryos, only a small number of RPGs (n = 5) (Figure S6E) and translational factors (n = 4) (Figure S6F) were modestly upregulated between the 4- and 8-cell stages, consistent with the onset of translational activation during this window.35 In contrast, m6A-RPL12-edited embryos showed aberrant upregulation of 22 RPGs at the 4-cell stage, followed by a failure to further increase expression at the 8-cell stage (Figure S6G), indicating disruption of normal temporal regulation. Translation initiation and elongation factors remained relatively stable (Figure S6H).

Direct comparison between control and edited embryos revealed that at the 4-cell stage, 64 RPGs and 13 translation-related factors were significantly downregulated in edited embryos (Figures 6D and 6E). This suppression persisted at the 8-cell stage, with 60 RPGs and 6 translation factors remaining significantly reduced (Figures 6F and 6G). Venn diagram analysis revealed minimal overlap in DEGs between stages, supporting the stage-specific effects of m6A-RPL12 disruption (Figures 6H and 6I). GO analysis further showed that mitochondrial organization, metabolic precursor generation, and translation were among the most affected pathways at both stages (Figures 6J and 6K). Pathways related to meiosis, chromosome organization, and cell cycle progression were also disrupted, suggesting broader developmental consequences beyond translation alone.

Together, these data demonstrate that m6A-RPL12 is required for proper activation of translation-related transcripts and maintenance of translational capacity during ZGA, thereby ensuring normal embryonic development.