- Research
- Open access
- Published:
Integrated m6A RNA methylation and transcriptomic analysis of Apostichopus japonicus under combined high-temperature and hypoxia stress
BMC Genomics volume 26, Article number: 363 (2025)
Abstract
Background
Global climate change has significantly increased environmental stress in marine ecosystems, with rising sea surface temperatures and declining dissolved oxygen (DO) levels. These stressors pose critical challenges to aquaculture, particularly for Apostichopus japonicus, an economically significant species in China. A. japonicus is highly sensitive to combined high-temperature and hypoxia stress, which disrupts physiological processes, suppresses immune responses, and increases mortality. While epigenetic mechanisms such as N6-methyladenosine (m6A) RNA modifications are known to regulate stress adaptation, their role under dual stressors in A. japonicus remains poorly understood.
Results
This study integrates m6A methylation sequencing (MeRIP-seq) and transcriptomic analysis (RNA-seq) to investigate molecular responses in A. japonicus under combined high-temperature (32 °C) and hypoxia (DO = 2 mg/L). Results show that approximately 90% of genes had 1–3 m6A peaks, with single peaks being the most frequent (∼ 60%). Genes with m6A modifications exhibited varying expression levels, with some showing significantly higher expression, suggesting a complex relationship between m6A methylation and stress-responsive gene expression. GO and KEGG enrichment analyses revealed that m6A-modified genes regulate pathways associated with oxidative stress, protein homeostasis, and energy metabolism, such as the PI3K-Akt and MAPK signaling pathways. Key stress-responsive genes, including HSP70, NOX5, and SLC7A11, exhibited dynamic m6A methylation changes, highlighting their roles in redox homeostasis and cellular resilience. Comparative analysis across experimental groups revealed distinct molecular responses to hypoxia, high-temperature stress, and their combination, with combined stress inducing more pronounced changes in m6A methylation and gene expression.
Conclusion
In this study, we explored the central regulatory role of m6A RNA methylation in the response of A. japonicus to the dual environmental stress of high-temperature and hypoxia. The findings show that m6A modification regulates the expression of key genes, allowing A. japonicus to effectively adapt to harsh environmental conditions. This study not only provides an important new perspective on the molecular stress recovery mechanism of marine invertebrates in the face of complex environmental stress, but it also provides theoretical support for aquaculture practice, assisting in the development of more stress-resistant aquaculture systems to deal with the severe challenges posed by global climate change.
Introduction
Due to the ongoing effects of global climate change, environmental pressures on Marine ecosystems have increased significantly. Since the 1970s, sea surface temperatures have shown a clear upward trend, rising by an average of about 0.13 °C per decade, significantly impacting marine environments and biomes [1]. The Yellow and Bohai Seas, critical places for Chinese aquaculture, have warmest at rates faster than the world-wide average, with extreme high-temperature events becoming more frequent during the summer [2]. Simultaneously, dissolved oxygen levels in these regions have consistently decreased due to stratification and eutrophication, leading to critically low oxygen levels below 2 mg/L, threatening aquatic life and ecosystem stability [3, 4]. These environmental changes pose significant threats to marine organisms by impairing physiological processes, suppressing immune responses, and reducing survival rates [5].
According to the Food and Agriculture Organization (FAO), extreme climate events have caused annual global economic losses of billions of dollars [6]. In China, the extreme hot weather season has become directly linked to increased mortality incidences in animal species, including the sea cucumber (Apostichopus japonicus) [7]. With an annual production exceeding 200,000 tons and a market value of over $4 billion, A. japonicus is one of the most important aquaculture species in China, particularly in Shandong and Liaoning provinces [8]. The optimum growth temperature of A. japonicus is 16 ∼ 21℃, when temperatures exceed 26 °C, causing oxygen consumption to increase. However, as the temperature rises further, metabolism is disrupted, leading to a decrease in oxygen consumption [7]. Furthermore, A. japonicus has shown a marked decrease in survival rates when exposed to temperatures exceeding 30 °C. Studies have demonstrated that at 32 °C, sea cucumbers experience significant physiological stress, including increased metabolic demand and reduced immune function, which may contribute to mortality in aquaculture systems [9, 10]. When the temperature of the water body is high and the oxygen content is low, the aquaculture system will be seriously affected, and this adverse environment will cause adverse effects on the health and survival of aquatic animals, such as reduced immune function, increased physiological burden and reduced survival rate [6].
N6-methyladenosine (m6A) is a common RNA modification essential for organisms to adapt to environmental stress. It is widely present in organismal mRNA and plays a key role in regulating protein synthesis, RNA stability, translation efficiency, and stress response, enabling cells to adjust to external changes without disrupting normal processes [11, 12]. Recent studies in marine invertebrates have shown that m6A modifications regulate immune responses under stress in species such as Mytilus coruscus [13], where m6A influences immune function during stress conditions. In Marsupenaeus japonicus, m6A modifications are involved in metabolic reprogramming during viral infections and regulate molting processes [14]. Additionally, exposure to environmental stressors like ammonia nitrogen has been shown to alter m6A modification patterns in Litopenaeus vannamei [11]. While m6A’s role in stress adaptation has been explored in other marine species, research in A. japonicus is limited, particularly regarding combined high-temperature and hypoxia stress. The response mechanisms of m6A modification to combined high-temperature and hypoxia stress may involve complex cellular signaling pathways related to energy metabolism, immune regulation, and oxidative stress, all of which are crucial for the organism’s survival under harsh environmental conditions [12, 15].
This study utilizes m6A methylation sequencing (MeRIP-seq) and RNA-seq to investigate the role of m6A modification in A. japonicus under combined high-temperature and hypoxia stress. We reveal that m6A modification regulates immune and oxidative stress-related genes in A. japonicus, particularly in oxygen-sensitive respiratory tree tissues, enhancing the species’ adaptability and survival in hypoxia and high-temperature environments. These findings provide new insights into how marine organisms adapt to complex environmental stresses, highlighting the epigenetic role of m6A modification.
Results
RNA sequencing data quality and mapping statistics
Each sequencing library generated around 40 million high-quality reads, with Q30 scores consistently surpassing 92.80%. The alignment rate to the A. japonicus reference genome ranged between 68.42% and 76.46%, with more than 45% of the reads being uniquely mapped (Supplementary Table S1). These outcomes confirm that the sequencing data met the necessary quality criteria for subsequent analyses.
A total of 12,768, 8,679, and 10,983 m6A peaks were identified in the T18DO7, T18DO2, and T32DO2 groups, respectively (Fig. 1B). Shared and unique peaks between groups are presented in Fig. 1C, with 4,188 peaks overlapping to all conditions. Approximately 90% of the genes had 1–3 m6A peaks, with single peaks being the most frequent, accounting for roughly 60% of the genes across all groups (Fig. 1D). Specifically, the T18DO2 group exhibited the highest proportion of genes with a single m6A peak (60.39%), followed by the T32DO2 group (56.62%) and the T18DO7 group (54.82%). Genes with more than three m6A peaks were rare, making up less than 10% in all groups (Fig. 1D).
Illustrates the experimental design and the general distribution patterns of m6A peaks under different stress conditions in A. japonicus. (A) Represents the experimental treatment groups, detailing the environmental stress conditions applied, including variations in temperature and dissolved oxygen levels. (B) Summarizes the total number of m6A peaks identified in each group, reflecting the overall methylation landscape under each condition. (C) Shows a Venn diagram displaying shared and unique m6A peaks among the different treatment groups, while (D) Categorizes the proportion of genes containing different numbers of m6A peaks. (E) Visualizes the distribution of m6A peaks across transcript regions, including 3’UTR, 5’UTR, CDS, start codon, and stop codon, emphasizing functional localization patterns. (F) Depicts the density distribution of m6A peaks along transcripts, highlighting positional preferences and regions with higher methylation density. (G) Identifies consensus motifs associated with m6A peaks, revealing sequence-specific regulatory features linked to stress adaptation. (H) and (I) present Venn diagrams comparing shared and unique m6A peaks and m6A-associated genes across group pairs, respectively. Finally, (J) displays the volcano plot illustrating the differential m6A peaks between groups and their associated genes
m6A peak distribution across gene functional elements
The distribution of m6A peaks across gene functional elements is illustrated in Fig. 1E. Most m6A peaks were localized within coding sequences (CDS), with the T18DO7, T18DO2, and T32DO2 groups showing 41.11%, 41.08%, and 42.44% of peaks in the CDS, respectively. In contrast, the 5’ untranslated region (UTR) had the fewest m6A peaks (3.41–3.66%). Peak counts were highest in T18DO7 (12,768 peaks), followed by T32DO2 (10,983 peaks), and lowest in T18DO2 (8,679 peaks).
The distribution of m6A peaks across different gene functional elements revealed specific trends (Fig. 1F). In the 5’ UTR, the proportion of m6A peaks increased gradually, reaching a slight plateau across all groups. A sharp rise in peak density was observed at the start of the CDS region, followed by a decline in the middle portion of the CDS. At the terminal region of the CDS, a secondary increase in peak density was noted, which then sharply declined. In the 3’ UTR, a moderate rise in m6A peak density was observed, followed by a final decrease towards the end. While similar trends were observed in all three groups, T18DO7 consistently exhibited higher m6A enrichment across all regions compared to T18DO2 and T32DO2, indicating that hypoxia and combined high-temperature stress may reduce overall m6A enrichment levels.
A motif refers to a short, recurring sequence of nucleotides in RNA, typically associated with the binding of specific proteins or modifications, such as m6A. In this study, motifs like ARAAACWKT (T18DO2) and RRAAACWTT (T32DO2) were significantly enriched in the stressed groups (Fig. 1G). The enrichment of these motifs was notably higher in the stressed groups compared to the control group, indicating a shift in gene regulation due to the stress. The motifs in the stressed conditions showed high similarity, likely reflecting the regulatory convergence of stress-related pathways activated by environmental stressors, such as heat and hypoxia.
Comparative m6A peak analysis between groups
A comprehensive analysis of m6A peaks revealed significant differences across the experimental groups. In the T18DO7_VS_T18DO2 comparison, 299 peaks were identified (Fig. 1H), with 71 peaks upregulated and 228 downregulated (P < 0.05,|log2FC| > 1) (Fig. 1J), including genes such as MMP14 (Matrix metallopeptidase 14), Sstr4 (Somatostatin receptor 4), NOX5 (NADPH oxidase 5), and others. The comparison between T18DO7 and T32DO2 identified 651 peaks (Fig. 1H), with 337 upregulated and 314 downregulated peaks (Fig. 1J), including TNR (Tenascin R), Slc6a5 (Solute carrier family 6 member 5), PRSS12 (Serine protease 12), among others. In the T18DO2_VS_T32DO2 comparison, 626 peaks were identified (Fig. 1H), with 466 upregulated and 160 downregulated (Fig. 1J), including genes such as NOX5, TNR, TMEM246 (Transmembrane protein 246), and others. Gene analysis related to m6A peaks showed that the T18DO7_VS_T18DO2 comparison revealed 178 distinct genes (Fig. 1I), while the T18DO7_VS_T32DO2 comparison showed 377 distinct genes (Fig. 1I), and T18DO2_VS_T32DO2 had 340 distinct genes (Fig. 1I). Notably, only 15 genes were common across all three comparisons (Fig. 1I), suggesting that hypoxia and high-temperature stress impact different gene sets, although some overlap exists. These findings highlight the significant impact of high-temperature and hypoxia stress on m6A methylation. High-temperature stress in T32DO2 caused more upregulated peaks, indicating a stronger influence on m6A modification compared to hypoxia. This suggests distinct stress-induced reprogramming of m6A modifications under different conditions.
Functional enrichment analysis of differential m6A peaks
To further investigate the role of m6A modification in A. japonicus under high-temperature (T32) and severe hypoxia (DO2) conditions, GO enrichment analysis was conducted on differentially m6A peak-associated genes across the three comparison groups: T18DO7_VS_T18DO2, T18DO7_VS_T32DO2, and T18DO2_VS_T32DO2. The results revealed substantial overlap in GO terms among the three groups (Fig. 2A), particularly in categories such as metabolic processes (e.g., primary metabolic processes, organic substance metabolic processes, and nitrogen compound metabolic processes), while also highlighting unique GO terms specific to each group, reflecting distinct biological responses under different stress conditions. In the biological process (BP) category, enriched terms were primarily associated with metabolic and cellular processes. The cellular component (CC) category showed significant enrichment in cellular cytoplasm, organelles, and membrane-bounded organelles. In the molecular function (MF) category, key terms included binding activities (e.g., protein binding and anion binding) and catalytic activities (e.g., transferase and hydrolase activities). These findings underscore the pivotal role that m6A modification plays a critical role in regulating both general fundamental biological functions and condition-specific responses in A. japonicus under high-temperature and hypoxia stress.
KEGG pathway enrichment analysis was performed for the comparison groups T18DO7_VS_T18DO2, T18DO7_VS_T32DO2, and T18DO2_VS_T32DO2 (Fig. 2B). The analysis identified a considerable overlap in enriched pathways across the groups, including metabolic pathways, the PI3K-Akt signaling pathway, and protein processing within the endoplasmic reticulum, which suggests shared mechanisms for responding to stress. Nonetheless, unique KEGG pathways were identified in each group, highlighting physiological responses tailored to specific conditions. In the comparison between T18DO7 and T18DO2, several unique pathways were enriched, including various types of N − glycan biosynthesis, proteoglycans in cancer, dopaminergic synapse, cAMP signaling pathway, PI3K − Akt signaling pathway, and focal adhesion. For T18DO7 vs. T32DO2, spliceosome was the unique pathway identified. In the comparison between T18DO2 and T32DO2, unique pathways such as the NOD − like receptor signaling pathway and shigellosis were observed.
Demonstrates the functional enrichment analysis of genes associated with differential m6A peaks across stress conditions. (A) GO annotation results for biological processes (BP), molecular functions (MF), and cellular components (CC) of genes associated with differential m6A peaks. (B) KEGG pathway enrichment analysis of m6A peak-associated genes
Shared and unique DEGs in A. japonicus under stress
RNA-seq analysis revealed distinct clustering of the three treatment groups in the PCA plot (Fig. 3A), with the T32DO2 group showing the greatest separation from the other two groups. This separation indicates a unique gene expression profile under extreme conditions of combined heat and hypoxia stress.
Visualizes transcriptional differences and functional characteristics of DEGs under varying stress conditions. (A) PCA plot shows distinct clustering of the three experimental groups. (B) Venn diagram summarizes shared and unique DEGs among the groups. (C) Bar plot displays the number of upregulated and downregulated DEGs in each pairwise comparison. (D) Volcano plots illustrate the magnitude and significance of transcriptional changes, identifying key stress-responsive genes. (E) GO annotation categorizes DEGs by biological functions. (F) KEGG pathway analysis highlights pathways involved in stress adaptation
The shared and unique differentially expressed genes (DEGs) among the groups are summarized in Fig. 3B. A core set of five DEGs (DMBT1(Deleted in malignant brain tumors 1), Cyp2u1 (Cytochrome P450 family 2 subfamily u member 1), ITIH3 (Inter-alpha-trypsin inhibitor heavy chain 3), hsp90ab1 (Heat shock protein 90 alpha family class B member 1) and slc5a8 (Solute carrier family 5 member 8)) was identified across all comparisons, highlighting their central roles in general stress responses. In addition, 26 DEGs were shared between the T18DO7_VS_T18DO2 and T18DO2_VS_T32DO2 comparisons, 29 DEGs were shared between T18DO7_VS_T18DO2 and T18DO7_VS_T32DO2, and 923 DEGs were shared between T18DO7_VS_T32DO2 and T18DO2_VS_T32DO2, reflecting transitional and temperature-specific molecular adaptations.
The DEGs identified in each pairwise comparison was detailed in Fig. 3C and D. From the comparison of T18DO7_VS_T18DO2, a total of 117 DEGs including 23 upregulated and 94 downregulated genes were detected. T18DO7_VS_T32DO2 revealed 1,257 DEGs with 361 upregulated and 896 downregulated genes. Through T18DO_VS_T32DO2 1,921 DEGs were found, with 1,070 being the upregulated and 851 the downregulated genes. For instance, in the T18DO7_VS_T18DO2 comparison, among the identified genes were OLFM4 (Olfactomedin 4), DMBT1 and SFTPD (Surfactant pulmonary defensin, whereas those like SLC2A5 (Solute carrier family 2 member 5) and HSP70IV (Heat shock protein 70 family member 4) were conspicuously expressed in T18DO7_VS_T32DO2. Also, there were such genes as AQP9 (Aquaporin 9), SLC2A5, and Hsp67B (Heat shock protein 67 kDa, family B) that stood out in the T18DO2_VS_T32DO2 comparison.
The Gene Ontology (GO) enrichment analysis (Fig. 3E) revealed commonalities among the three treatment groups, with differentially expressed genes (DEGs) predominantly enriched in biological processes like oxidation-reduction, protein folding, and transcriptional regulation, indicating a broad molecular response to environmental stress. In contrast, KEGG pathway analysis (Fig. 3F) highlighted group-specific variations. For instance, the T18DO7_VS_T18DO2 comparison showed enrichment in pathways related to steroid biosynthesis and DNA replication, suggesting shifts in metabolic and cellular activities under normothermic hypoxia. In the T18DO7_VS_T32DO2 comparison, pathways associated with protein digestion and lipid metabolism were enriched, reflecting adaptations to the combined effects of heat and hypoxia. Additionally, the T18DO2_VS_T32DO2 comparison revealed significant enrichment in pathways such as glutathione metabolism, emphasizing oxidative stress responses under combined high-temperature and hypoxia.
Role of m6A methylation in gene expression regulation
In the T18DO2, T18DO7, and T32DO2 groups, genes with m6A modifications consistently showed significantly higher expression compared to those without m6A modification (Fig. 4A). Integrated analysis of m6A-seq and RNA-seq data revealed that, under combined high-temperature and hypoxia stress conditions, genes with m6A modifications showed varying expression changes, either higher or lower. Across the treatment groups, m6A-modified genes were more likely to exhibit synchronized changes in both m6A methylation and mRNA expression levels. A higher proportion of genes fell into the “hyper-up” (increased m6A and mRNA levels) or “hypo-down” (decreased m6A and mRNA levels) categories. Specifically, in the T18DO7_VS_T32DO2 comparison, 56 genes showed significant changes, including 1 hyper-up gene and 53 hypo-down genes. In the T18DO7_VS_T18DO2 comparison, 160 genes were differentially affected, with 35 hyper-up and 80 hypo-down genes. Meanwhile, in the T18DO2_VS_T32DO2 comparison, 186 genes exhibited significant changes, including 104 hyper-up and 34 hypo-down genes. The distribution of these genes across the four categories, as illustrated in the four-quadrant plots and Venn diagrams (Fig. 4B and C), highlights the critical role of m6A methylation in modulating gene expression in response to environmental stress in A. japonicus.
Integrated Analysis of m6A Methylation and mRNA Expression. (A) Displays the expression levels of genes with and without m6A modifications. (B) Venn diagram shows the overlap between differentially m6A-methylated genes and DEGs. (C) Categorizes genes into four groups based on changes in m6A methylation and expression levels, reflecting diverse regulatory responses. (D) Functional enrichment analysis identifies key biological processes and pathways associated with genes showing coordinated changes in methylation and expression. (E) Presents representative genes involved in stress adaptation, emphasizing their roles in cellular responses
GO and KEGG pathway enrichment analyses were performed using the key genes identified from the integrated analysis of m6A-seq and RNA-seq data (Fig. 4D). In the comparison of T18DO7_VS_T18DO2, the abundance of ATP binding, RNA helicase activity, and nucleic acid binding was observed in the genes of the biological process, which is linked to RNA processing, and energy metabolism which are the underlying mechanisms working. During the comparison of T18DO7_VS_T32DO2, the genes were dense in the processes responsible for protein folding and the mitochondria, signifying the cellular stress responses. Similarly, in the T18DO2_VS_T32DO2 comparison, processes such as Hsp70 protein binding and SRP-dependent protein targeting were enriched in the RNA processing and energy metabolism, and protein synthesis, etc., under stress conditions. The fact that m6A modifications in gene expression regulation become extremely threatened particularly in environmental stress could be highly supported by shared GO terms across all comparisons, such as methyltransferase activity, hydrolase activity, and ribosome biogenesis. Furthermore, with regard to KEGG pathway analysis, such pathways as glycosaminoglycan degradation and ribosome biogenesis in eukaryotes emerged as significantly enriched particularly in the T18DO2_VS_T32DO2 comparison implying m6A’s involvement in cellular metabolism and protein synthesis under high-temperature and hypoxic conditions.
The expression levels of representative genes in T18DO7, T18DO2, and T32DO2 are shown in Fig. 4E. Genes such as Mettl16 (Methyltransferase like 16), RNF19A (Ring finger protein 19 A), and Slc6a5 are involved in various biological processes linked to stress adaptation. Specifically, these genes play critical roles in RNA processing (Mettl16, TRMT2A (tRNA methyltransferase 2 A)), protein folding and degradation (RNF19A, PRSS12), antioxidative responses (SLC7A11 (Solute carrier family 7 member 11), DHCR24 (24-Dehydrocholesterol reductase)), and mitochondrial function (AMOT (Angiomotin), HAO1(Hydroxyacid oxidase 1)). Additionally, genes such as Bambi (BMP and activin membrane-bound inhibitor) and Trib2 (Tribbles pseudokinase 2) are involved in signaling pathways related to cellular survival under stress. Collectively, these genes regulate key pathways associated with high-temperature and hypoxic stress responses, contributing to cellular homeostasis, metabolism, and stress adaptation.
Discussion
The combined effects of high temperature and hypoxia on marine aquaculture
Rising temperatures and decreasing oxygen levels are increasingly affecting marine ecosystems [16, 17], particularly in areas like the Yellow and Bohai Seas, where the combined effects of high temperature and low dissolved oxygen (DO) pose significant threats to marine life, including A. japonicus. These stressors impair physiological functions, such as metabolic processes and immune responses, and can drastically affect the survival of aquatic organisms [18, 19].
In this study, we simulate typical environmental conditions of high temperature (32 °C) and hypoxia (2 mg/L DO) to investigate the physiological and molecular response mechanisms of A. japonicus. Previous studies on other marine invertebrates, such as shrimp (Penaeus vannamei) [20, 21] and sea urchins (Strongylocentrotus purpuratus) [22,23,24], have shown similar stress responses, including oxidative damage, immunosuppression, and metabolic disruptions under high-temperature and hypoxic conditions. Our study aims to deepen the understanding of how these combined stresses impact the survival and adaptability of A. japonicus in these challenging environments.
In combined high-temperature and hypoxia stress, antioxidant-related genes such as NOX5 and DHCR24 were significantly upregulated, indicating that A. japonicus may have a special antioxidant mechanism that is relatively rare in other marine invertebrates. At the same time, the expression levels of metabolic regulatory genes ALPL (Alkaline phosphatase) and AMOT changed significantly, reflecting the metabolic reprogramming of A. japonicus in lipid metabolism and energy generation pathways. These findings are consistent with the research results of Litopenaeus vannamei, further emphasizing the key role of metabolic regulation in coping with environmental stress [25,26,27].
Gene regulation and molecular adaptation under high-temperature and hypoxia stress in A. japonicus
A. japonicus exhibits regulatory mechanisms in response to combined stress of high-temperature and hypoxia, particularly in terms of gene expression and m6A methylation. These results provide different perspectives on the adaptive strategies of organisms. Genes related to heat shock proteins (HSPs), antioxidant defense, and metabolic regulation showed significant changes in expression levels and m6A modification patterns. Under stress conditions, the products and m6A methylation patterns of heat shock proteins such as HSP70 and HSP90 (Heat shock protein 90) increased significantly. When the environment becomes harsh, heat shock proteins play a vital role in stress environments by maintaining protein stability and cellular homeostasis [28, 29]. Results show that m6A modification is involved in post-transcriptional regulation. When A. japonicus respond to environmental stress, they increase post-transcriptional epigenetic modifications, which promotes the expression regulation of related genes. Similar regulatory patterns have also been reported in shrimps, such as the significant increase in the expression of heat shock proteins under high-temperature stress and hypoxia [30].
This study found that the expression of antioxidant-related genes NOX5 and SLC7A11 changed significantly under the dynamic regulation of m6A methylation. These genes maintain redox balance by regulating the level of reactive oxygen species (ROS) in cells [31, 32]. Notably, the NADPH oxidase encoding gene NOX5 showed a significant upregulation in A. japonicus, suggesting a molecular adaptation mechanism under the combined effects of high temperature and hypoxia-induced oxidative stress. Additionally, changes in traditional antioxidant enzymes, such as SOD and CAT, were also observed under stress conditions in A. japonicus, and similar regulatory patterns have been reported in other marine organisms like shrimps and sea urchins [33, 34].
Metabolic regulation has been identified as an important component of the stress response in A. japonicus [35, 36]. The finding that genes such as ALPL and AMOT underwent significant changes in expression levels and m6A methylation in this study. ALPL, which plays a critical role in bone mineralization and energy metabolism, and AMOT, involved in cell migration and angiogenesis, are key regulators of metabolic pathways. Suggests that metabolic reprogramming is an important component of the adaptive mechanism of A. japonicus. This metabolic adjustment is mainly focused on oxidative phosphorylation and fatty acid metabolism pathways, meeting the higher energy demand during stress [37]. Similar metabolic changes have also been reported in other marine invertebrates, with studies in shrimp showing enhanced lipid mobilization under conditions of high temperature and hypoxia [38, 39].
Functional roles of signaling pathways in stress adaptation
Signaling pathways play an important role in the process of cell adaptation to environmental stress. This study revealed that A. japonicus uses a specific signaling network to cope with the combined stress of high temperature and hypoxia. GO and KEGG pathway analysis showed that oxidative stress response, protein homeostasis, and metabolism-related pathways are important in the adaptation process of A. japonicus, reflecting the complex adaptive strategies of organisms in response to adverse living conditions.
This study found that the PI3K-Akt and MAPK signaling pathways were significantly enriched under various stress conditions, indicating that they play a key role in cell survival, immune response and metabolic regulation. The PI3K-Akt pathway prevents cell apoptosis by alleviating oxidative stress and maintaining metabolic homeostasis, a function that has been widely confirmed in a variety of organisms [40]. In A. japonicus, the activation of this pathway may provide important support for cell adaptation to high temperature and hypoxia stress by reducing oxidative damage and supporting energy metabolism. As a key stress signal transduction regulator, the significant upregulation of the MAPK signaling pathway may be related to the activation of immune function and tissue repair processes. Similar adaptation patterns have also been observed in Litopenaeus vannamei and Strongylocentrotus purpuratus, where the MAPK pathway is considered a conserved mechanism for coping with heat or hypoxia stress [41,42,43,44].
Pathways related to protein processing were significantly enriched in the endoplasmic reticulum (ER), and high levels of expression of heat shock proteins (HSP70 and HSP90) highlighted the necessity of cells to maintain protein quality under stress conditions, suggesting that cells respond to extreme environments by preventing protein misfolding and aggregation. Consistent with studies in other marine organisms (such as A. japonicus and bivalves), these results further verify the importance of protein homeostasis, which is essential for A. japonicus to survive under stress conditions [45,46,47].
Adaptations of A. japonicus to environmental stress and future directions
This study revealed important adaptations in signaling pathways related to environmental stress. For example, the significant upregulation of calcium signaling-related genes emphasizes its crucial role in stress adaptation in A. japonicus. Calcium signaling regulates essential cellular processes, including apoptosis, immune responses, and oxidative stress, and has been shown to play a role in stress responses in other marine invertebrates, such as corals and shrimp. However, the specific regulatory mechanisms of calcium signaling in A. japonicus remain poorly understood [48, 49]. Additionally, the enrichment of neuroactive ligand-receptor interaction pathways points to potential cross-talk between metabolic and immune signaling, warranting further exploration.
Under the combined stress of high temperature and hypoxia, signaling pathways related to oxidative stress regulation (such as glutathione metabolic pathways) were significantly enriched. Genes such as SLC7A11 in this pathway play a central role in neutralizing ROS and maintaining redox homeostasis, and are of key significance in the response of A. japonicus to stress. Changes in m6A methylation of these genes suggest that cells may adapt to environmental stress by precisely regulating gene expression. These findings suggest that signaling pathways have dual functions in alleviating oxidative damage and promoting cell recovery.
Conclusions
This study explored the regulatory role of m6A RNA methylation in A. japonicus response to combined stress of high temperature and hypoxia. The results showed that key genes such as HSP70, NOX5, and SLC7A11, as well as PI3K-Akt and MAPK signaling pathways, played a key role in regulating oxidative stress, metabolic reprogramming, and enhancing immune response in A. japonicus. The complex relationship between m6A methylation and gene expression further suggests that this epigenetic mechanism helps A. japonicus adapt to environmental stress more efficiently by finely regulating gene expression. This mechanism enables A. japonicus to better cope with oxidative damage, maintain metabolic balance, and improve immunity, thereby increasing survival in extreme environments. Next, research will focus on exploring the long-term effects of stress and expanding to additional tissues to gain a more comprehensive understanding of the stress response mechanisms in A. japonicus.
Materials and methods
Experimental design and stress treatment protocols
A total of 90 A. japonicus individuals (mean weight 106 ± 15.3 g, mean ± SD) were obtained from Shandong Anyuan Aquatic Products Co., Ltd. All specimens were bred under identical growth conditions. Before the experiments, the A. japonicus were acclimatized under controlled conditions (18 ± 0.5 °C, salinity 30‰, and a 12 h light:12 h dark photoperiod) for 14 days. During this period, they were fed a commercial diet twice daily, and uneaten feed was removed to maintain water quality. The experiment included three treatment groups, each consisting of 3 biological replicates, each containing 10 individuals housed in separate 50 L tanks with a recirculating water supply, biofiltration system, and aeration. The stress conditions were: T18DO7 (control: 18 °C, DO = 7 mg/L), T18DO2 (hypoxia: 18 °C, DO = 2 mg/L), and T32DO2 (combined high temperature and hypoxia: 32 °C, DO = 2 mg/L). A gradual temperature increase of 1 °C every 12 h was achieved using temperature control rods. To establish low oxygen conditions, nitrogen gas was bubbled into the water to reduce oxygen levels. The DO was continuously monitored using a DO meter to ensure that the required levels were maintained during the exposure period. The stress conditions were maintained for 24 h, simulating environmental stress scenarios (Fig. 1A). At the end of the exposure period, three individuals from each treatment groups were randomly sampled to ensure unbiased selection. Respiratory tree tissues were rapidly collected, frozen in liquid nitrogen, and stored at − 80 °C for molecular analyses.
Methylated RNA immunoprecipitation sequencing (MeRIP-Seq) and RNA sequencing
Total RNA was extracted from respiratory tree tissues using TRIzol reagent (Invitrogen, USA) according to the manufacturer’s protocol. RNA integrity and concentration were assessed using agarose gel electrophoresis, NanoDrop ND-1000, and Bioanalyzer 2100. Poly(A) RNA was isolated from 50 µg of total RNA using Dynabeads™ Oligo(dT)25 (Thermo Fisher, USA). To enhance m6A peak resolution, RNA was fragmented to ∼ 100 nt. Immunoprecipitation was performed with 10 µg fragmented RNA using 2 µg m6A antibody (Synaptic Systems, Germany) in IP buffer (50 mM Tris-HCl, 750 mM NaCl, 0.5% Igepal CA-630) at 4 °C for 2 h. RNA-antibody complexes were captured using Dynabeads™ Protein A (Thermo Fisher, USA), eluted, and recovered via ethanol precipitation. Both input and immunoprecipitated RNA were reverse-transcribed using SuperScript™ II (Invitrogen, USA). cDNA libraries were prepared with dual-index adapters and amplified via PCR. The libraries were sequenced on an Illumina NovaSeq™ 6000 platform (LC-Bio Technology Co., Ltd., China).
Analysis of MeRIP-Seq and RNA-Seq data
Raw sequencing data were quality-filtered using fastp (v0.18.0) to remove low-quality bases and adapter sequences. Reads were mapped to the ribosomal RNA (rRNA) database with Bowtie2 (v2.2.8), and rRNA-mapped reads were excluded. The remaining clean reads were aligned to the A. japonicus reference genome (NCBI Assembly Version ASM275485v1) using HISAT2 (v2.1.0). For MeRIP-Seq, m6A-enriched regions (peaks) were identified using the R package exomePeak2 (v1.12.0) with a default p-value threshold of < 1e-5. Peaks were annotated using ANNOVAR to classify their genomic distribution (5’ UTR, CDS, or 3’ UTR). Significant peaks were defined as those meeting|log2(fold change)| > 1 and false discovery rate (FDR) < 0.05. Motif discovery was performed with HOMER (v4.11), and visualizations were generated using ggseqlogo and ggplot2. Differential peak analysis was conducted using DiffBind (v2.8.0), and overlapping peaks were identified as common peaks with the nearest gene considered peak-related.
For RNA-Seq, clean reads were aligned to the reference genome with HISAT2 (v2.1.0) and assembled into transcripts using StringTie (v1.3.4). Gene expression was quantified as FPKM using RSEM (v1.2.19). Differentially expressed genes (DEGs) were identified using DESeq2 (v1.20.0) with thresholds of|log2(fold change)| ≥ 2 and FDR < 0.05. Functional enrichment analysis of peak-associated genes and DEGs was performed using GO.db (v3.14.0) and the KEGG database (Release 101), with significance defined at adjusted P < 0.05.
Integrated analysis of differential m6A methylation and transcriptomic changes
Differentially expressed genes and m6A methylated peaks were analyzed using edgeR (v4.4.1) on biological replicates, applying thresholds of|log2(fold change)| > 1 and P < 0.05 to identify significant differences. Genes were categorized into four groups based on methylation and expression changes: “hyper-up” (increased methylation and expression), “hypo-down” (decreased methylation and expression), “hyper-down,” and “hypo-up.” A four-quadrant diagram was generated to visualize these relationships. GO and KEGG enrichment analyses were performed on genes in each category to identify pathways potentially regulated by m6A methylation.
Statistical analysis
Statistical analyses were applied to both MeRIP-Seq and RNA-Seq datasets. Differences between two groups were assessed using a two-tailed Student’s t-test, * and *** indicate significant difference at P < 0.05 and P < 0.001. Data are presented as mean ± SD (n = 3, biological replicates).
Data availability
The dataset supporting the conclusions of this article is available in the National Center for Biotechnology Information (NCBI) repository with the accession number of PRJNA1203921.
Abbreviations
- m6A:
-
N6-methyladenosine
- CDS:
-
Coding sequence
- 3’UTR:
-
3’ untranslated region
- 5’UTR:
-
5’ untranslated region
- MeRIP-seq:
-
Methylated RNA immunoprecipitation sequencing
- GO:
-
Gene Ontology
- KEGG:
-
Kyoto Encyclopedia of Genes and Genomes
- DEGs:
-
Differentially expressed genes
References
Li Y, et al. High-quality sea surface temperature measurements along Coast of the Bohai and yellow seas in China and their long‐term trends during 1960–2012. Int J Climatol. 2019;40:63–76.
Shi W, Wang M. Satellite views of the Bohai Sea, yellow Sea, and East China Sea - ScienceDirect. Prog Oceanogr. 2012;104(10):30–45.
Chen Y, et al. Seasonal variability in dissolved oxygen in the Bohai Sea, China. J Oceanol Limnol. 2022;40(1):78–92.
Yang F, et al. Nitrogen from agriculture and temperature as the major drivers of deoxygenation in the central Bohai sea. Sci Total Environ. 2023;893:164614.
Kang B, et al. Climate change impacts on China’s marine ecosystems. Rev Fish Biol Fish. 2021;31:599–629.
Chen Y, et al. Assessment of fish vulnerability to climate change in the yellow sea and Bohai sea. Mar Freshw Res. 2019;71(7):729–36.
Huo D, et al. Impact of hypoxia stress on the physiological responses of sea cucumber Apostichopus Japonicus: respiration, digestion, immunity and oxidative damage. PeerJ. 2018;6:e4651.
Xu M, et al. Coastal aquaculture farms for the sea cucumber Apostichopus japonicus provide spawning and first-year nursery grounds for wild black rockfish, Sebastes schlegelii: A case study from the Luanhe river estuary, Bohai Bay, the Bohai Sea, China. Front Mar Sci. 2022;9:911399.
Huo D, et al. Time course analysis of immunity-related gene expression in the sea cucumber Apostichopus japonicus during exposure to thermal and hypoxic stress. Fish Shellfish Immunol. 2019;95:383–90.
Wang F, et al. Effects of acute temperature or salinity stress on the immune response in sea cucumber, Apostichopus japonicus. Comp Biochem Physiol Part A: Mol Integr Physiol. 2008;151(4):491–8.
Lei Y, et al. Dynamic N6-methyladenosine RNA methylation landscapes reveal epi-transcriptomic modulation induced by ammonia nitrogen exposure in the Pacific whiteleg shrimp Litopenaeus vannamei. J Hazard Mater. 2023;458:131996.
Sun Y, et al. Transcriptome-wide methylated RNA Immunoprecipitation sequencing profiling reveals m6A modification involved in response to heat stress in Apostichopus japonicus. BMC Genomics. 2024;25(1):1071.
Li H, et al. Nanopore direct RNA sequencing Deciphering the m6A epitranscriptome and its role in the innate immune response of the mussel mytilus Coruscus: a prospective study of RNA modification in marine molluscs. Aquaculture. 2025;598:742040.
Lei Y, et al. A first glimpse into the M6A modification machinery of shrimp: genomic features, expression patterns and potential roles in molting regulation. Aquaculture Rep. 2023;29:101493.
Ji Q, et al. A characterization of the RNA modification response to starvation under low temperatures in large yellow croaker (Larimichthys crocea). Fishes. 2024;9(1):41.
Shoergashova S, et al. Assessment of Temporal and Spatial variations of water quality parameters in the Zarafshan river basin. CLEAN–Soil Air Water. 2024;52(8):2300454.
Qian C, et al. Hypoxia and warming take sides with small marine protists: an integrated laboratory and field study. Sci Total Environ. 2023;882:163568.
Rumondang A, Khobir ML, Mutiaragusti M. HIGH INFLUENCE OF WATER AND ENVIRONMENTAL QUALITY ON RESPONSE, BEHAVIOUR AND GROWTH BATAK FISH FRY (Neolissochilus thienemanni). Jurnal Perikanan Unram. 2024;14(4):1892–901.
Heriyati E, Ilmi MB, Suryono CA. Impact of water quality and phytoplankton on juvenile vannamei shrimp growth in Low-Salinity ponds. J Mar Biotechnol Immunol. 2024;2(3):6–11.
Jia X, et al. Immune responses of Litopenaeus vannamei to thermal stress: a comparative study of shrimp in freshwater and seawater conditions. Mar Freshw Behav Physiol. 2014;47(2):79–92.
Barajas-Sandoval DR, et al. Effect of Temporal thermal stress on Penaeus vannamei: growth performance and physiological plasticity. Comp Biochem Physiol A: Mol Integr Physiol. 2024;295:111653.
Johnstone J, et al. Effects of elevated temperature on gonadal functions, cellular apoptosis, and oxidative stress in Atlantic sea urchin arbacia punculata. Mar Environ Res. 2019;149:40–9.
Liu A, et al. High temperature influences DNA methylation and transcriptional profiles in sea urchins (Strongylocentrotus intermedius). BMC Genomics. 2023;24(1):491.
Pinsino A, Matranga V. Sea urchin immune cells as sentinels of environmental stress. Dev Comp Immunol. 2015;49(1):198–205.
Wang Z, et al. Taurine metabolism is modulated in Vibrio-infected Penaeus vannamei to shape shrimp antibacterial response and survival. Microbiome. 2022;10(1):213.
Zhang X, et al. The potential role of eyestalk in the immunity of Litopenaeus vannamei to vibrio infection. Fish Shellfish Immunol. 2022;121:62–73.
Cui C, et al. Single-cell RNA-seq revealed heterogeneous responses and functional differentiation of hemocytes against white spot syndrome virus infection in Litopenaeus vannamei. J Virol. 2024;98(3):e01805–23.
Kravats AN et al. Functional and physical interaction between yeast Hsp90 and Hsp70. Proceedings of the National Academy of Sciences. 2018;115(10):E2210-E2219.
Mitra S, et al. Bipartite role of heat shock protein 90 (Hsp90) keeps CRAF kinase poised for activation. J Biol Chem. 2016;291(47):24579–93.
Sun X, et al. Transcriptome-wide identification and analysis reveals m6A regulation of metabolic reprogramming in shrimp (Marsupenaeus japonicus) under virus infection. BMC Genomics. 2024;25(1):1103.
Xu Y, et al. METTL3 promotes lung adenocarcinoma tumor growth and inhibits ferroptosis by stabilizing SLC7A11 m6A modification. Cancer Cell Int. 2022;22(1):11.
Ji F-H, et al. FTO prevents thyroid cancer progression by SLC7A11 m6A methylation in a ferroptosis-dependent manner. Front Endocrinol. 2022;13:857765.
Parrilla-Taylor DP, Zenteno-Savín T. Antioxidant enzyme activities in Pacific white shrimp (Litopenaeus vannamei) in response to environmental hypoxia and reoxygenation. Aquaculture. 2011;318(3–4):379–83.
Muralisankar T, et al. Growth, biochemical, antioxidants, metabolic enzymes and hemocytes population of the shrimp Litopenaeus vannamei exposed to acidified seawater. Volume 239. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology; 2021. p. 108843.
Xia B, et al. Explore the gene regulation mechanism under high temperature and hypoxia stress in the sea cucumber Apostichopus japonicus by high-throughput RNA-seq technology. Aquaculture Rep. 2024;36:102151.
Yang Q, et al. Transcriptome analysis of sea cucumber (Apostichopus japonicus) in Southern China under heat stress. Aquaculture Rep. 2024;36:102036.
Jiang M, Fang H, Tian H. Metabolism of cancer cells and immune cells in the initiation, progression, and metastasis of cancer. Theranostics. 2025;15(1):155–88.
González-Ruiz R, et al. The combination of hypoxia and high temperature affects heat shock, anaerobic metabolism, and Pentose phosphate pathway key components responses in the white shrimp (Litopenaeus vannamei). Cell Stress Chaperones. 2023;28(5):493–509.
Olsen L, Thum E, Rohner N. Lipid metabolism in adaptation to extreme nutritional challenges. Dev Cell. 2021;56(10):1417–29.
Gao S, et al. Effects of WSSV and ammonia nitrogen stress on the enzyme activity and transcriptome of Litopenaeus vannamei. Aquaculture Rep. 2024;39:102412.
Lulu Y, et al. Transcriptomic analysis of Crassostrea Sikamea × Crassostrea angulata hybrids in response to low salinity stress. PLoS ONE. 2017;12(2):e0171483.
Zhang M, et al. Effects of flow velocity on the growth and survival of haliotis discus Hannai larvae in the recirculating upflow system from the point of energy metabolism. Front Mar Sci. 2021;8:763269.
Ruocco N, et al. New inter-correlated genes targeted by diatom-derived polyunsaturated aldehydes in the sea urchin paracentrotus lividus. Ecotoxicol Environ Saf. 2017;142:355–62.
Polinski JM, et al. Unique age-related transcriptional signature in the nervous system of the long-lived red sea urchin mesocentrotus franciscanus. Sci Rep. 2020;10(1):9182.
Kong X, et al. PDIA6 involves the thermal stress response of Razor clam, sinonovacula constricta. Fish Shellfish Immunol. 2022;131:766–74.
Hu Z, et al. Mechanisms of heat and hypoxia defense in hard Clam: insights from transcriptome analysis. Aquaculture. 2022;549:737792.
Chen J, et al. Integrated transcriptomic and metabolomic analysis sheds new light on adaptation of pinctada fucata martensii to short-term hypoxic stress. Mar Pollut Bull. 2023;187:114534.
Rathinam RB, et al. The immune system of marine invertebrates: earliest adaptation of animals. Comp Immunol Rep. 2024;7:200163.
Chu L, et al. TORC1 regulates thermotolerance via modulating metabolic rate and antioxidant capacity in scallop argopecten irradians irradians. Antioxidants. 2024;13(11):1359.
Funding
This work was supported by the Key Technology Research and Development Program of Shandong Province (2023CXGC010410, 2023LZGC019, 2022TZXD005), Agriculture Industry System of Shandong Province (SDAIT-22-03), Natural Science Foundation of Shandong Province (ZR2021MC023).
Author information
Authors and Affiliations
Contributions
WQ carried out the molecular genetic studies, participated in the sequence alignment and drafted the manuscript. SQZ, SYL, XHH and XHX participated in the sequence alignment. WQ and SQZ participated in the design of the study and performed the statistical analysis. WQ, YWF, JMY and GHS conceived of the study, and participated in its design and coordination. All authors read and approved the final manuscript.
Corresponding author
Ethics declarations
Ethics approval and consent to participate
The animal study was approved by the Institutional Animal Care and Use Committee of the Ludong University (protocol number LDU-IRB20210308NXY). The study was conducted in accordance with the local legislation and institutional requirements.
Consent for publication
Not applicable.
Competing interests
The authors declare no competing interests.
Additional information
Publisher’s note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Electronic supplementary material
Below is the link to the electronic supplementary material.
Rights and permissions
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/.
About this article
Cite this article
Wang, Q., Zhang, S., He, X. et al. Integrated m6A RNA methylation and transcriptomic analysis of Apostichopus japonicus under combined high-temperature and hypoxia stress. BMC Genomics 26, 363 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12864-025-11532-x
Received:
Accepted:
Published:
DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12864-025-11532-x