1. Introduction
Ruditapes philippinarum is an important bivalve for large-scale intertidal aquaculture in China, possesses high economic and ecological value. Copper sulfate is a common fishery drug used in aquaculture to control parasites, bacteria, and cyanobacteria; however, its misuse, overuse, and abuse can easily lead to the accumulation of copper in water bodies and sediments, making it a persistent source of pollution.
To comprehensively evaluate the impact of environmental stress on Ruditapes philippinarum, this study employed an integrated multi-index analysis method. The Integrated Biomarker Response (IBR) represents a robust approach for integrating multi-biomarker responses in the evaluation of ecological risk and stress exposure. The star plot illustrates the biomarker transformation data in the hepatopancreas of Ruditapes philippinarum and Chlamys farreri under microplastic exposure in Laizhou Bay. This graphical method allows intuitive visualization of response differences among multiple biomarkers across various concentration groups, and the enclosed area corresponds to the IBR value 1.
Elucidating the organism’s response mechanisms at the molecular level is equally crucial. Transcriptomics, a branch of molecular biology, focuses on the entire RNA transcriptome of specific cells or tissues within a biological entity at a given time. It enables a comprehensive analysis of gene expression characteristics, revealing differentially expressed genes and their impact on cellular function 2–4. Transcriptomics techniques have been widely applied to studies on stress responses across various organisms. Transcriptomics provides a new perspective on the hypoxic stress response of aquatic organisms 5. Association analysis of transcriptomic and proteomic data from the hepatopancreas of Eriocheir sinensis identified the expression patterns of differentially expressed genes and proteins6.GO and KEGG functional enrichment analyses of differentially expressed genes predicted the key genes and molecular pathways underlying the response to high-temperature stress in Mizuhopecten yessoensis, laying a theoretical foundation for the breeding of thermo-tolerant varieties of this scallop7.
This study systematically elucidated the dynamic response mechanisms of Ruditapes philippinarum to copper stress, filling a gap in existing research by integrating physiological and molecular data.
2. Materials and Methods
2.1. Biological Materials and Experimental Reagents
The experimental Ruditapes philippinarum were purchased from Tianjin Hengqian Aquaculture Co., Ltd. and transported to the laboratory under ice temperature conditions for temporary acclimation. First, the clams were rinsed with seawater to remove surface debris, dead individuals were removed, and those with intact shells, uniform specifications, and vigorous activity were selected for temporary acclimation for 5-6 days. The temporary acclimation container was a 100 cm × 80 cm × 50 cm transparent glass tank. The temporary seawater was prepared using Haiding Technology aquaculture special synthetic lobster salt (salinity 25-28), and was continuously aerated for 24 hours before use. During the acclimation period, 2/3 of the seawater was replaced every 48 hours, and 300 mL of Chlorella was fed to each tank after water change. The environmental parameters for temporary acclimation were controlled as follows: pH 7.20 ± 0.16, dissolved oxygen 8.45 ± 0.56 mg/L, water temperature (22 ± 2)°C. Clam activity was observed daily, and dead individuals were promptly removed. Subsequent experiments were conducted when the mortality rate was below 5%. Feeding was stopped 24 hours before the experiment, and individuals with vigorous vitality and sensitive stress responses were selected as experimental subjects.
The heavy metal copper used in the experiment was analytical grade copper sulfate pentahydrate (CuSO₄·5H₂O) crystals from Tianjin Zhonglian Chemical Reagent Co., Ltd. A 1 g/L copper sulfate stock solution was prepared and diluted to the required concentration during the experiment. Artificial seawater was prepared using Haiding Technology brand synthetic lobster salt produced by Jiangxi Haizhi Technology Co., Ltd. All experimental water was deionized water.
2.2. Experimental Design
Previous acute toxicity experiments in our study determined that the 48h median lethal concentration (LC50) of copper sulfate to Ruditapes philippinarum was 5.68 mg/L. Based on this, two stress concentrations were selected: 0.0568 mg/L (low concentration, 1/100 of 48h LC50) and 0.568 mg/L (high concentration, 1/10 of 48h LC50). Each experimental tank (specification 1 m × 0.8 m × 0.5 m) was filled with 200 L of artificial seawater. A semi-static exposure method was used for the stress experiment. During the experiment, continuous aeration was maintained (aeration intensity sufficient to keep dissolved oxygen ≥ 7.5 mg/L), and 1/2 volume of fresh seawater with the same copper concentration was replaced every 24 hours to maintain a stable pollutant concentration. Three treatment groups were used: a low concentration group (0.0568 mg/L), a high concentration group (0.568 mg/L), and a blank control group. Each treatment group had 3 replicates, and 30 vigorous Ruditapes philippinarum were stocked in each replicate. The stress experiment lasted for 14 days. Water quality parameters were strictly controlled and kept stable during the experiment: water temperature (20 ± 0.5)°C, pH (7.8 ± 0.5), salinity 30, and dissolved oxygen (8.0 ± 0.2) mg/L, consistent with the parameters during the acclimation period to minimize interference from environmental factors.
2.3, Sample Collection and Enzyme/Lipid Metabolite Determination
After the start of the experiment, samples were taken on Day 0, 2, 7, and 14 (D0, D2, D7, D14). For each experimental group, 10 Ruditapes philippinarum were randomly selected; 5 were used for enzyme activity, tissue lipid substance, and protein content determination, and the other 5 were used for transcriptome analysis sample preparation. During sampling, the shell surface was wiped dry with gauze, placed on an ice tray, and rapidly dissected. The hepatopancreas tissue was stripped and placed in a cryotube, frozen in liquid nitrogen for 6 hours, and then transferred to a -80°C freezer for storage.
The processing procedure for tissue samples used for biochemical indicator analysis was as follows: the frozen sample was removed and placed in a centrifuge tube. Physiological saline (0.9%) was added at a 1:9 (m/V) ratio to the tissue sample. Glass beads were added simultaneously, and the system pH was adjusted to 7.2. Tissue homogenization was performed under ice bath conditions (12,000 r/min, 60 Hz, 10 s each time, continuous 4-5 times). After homogenization, the homogenate was placed in an ice bath centrifuge and centrifuged at 4°C, 4,000 r/min for 10 min. The supernatant was taken and centrifuged again under the same conditions. After the second centrifugation, the supernatant was taken for indicator determination. The indicators measured in this experiment included the activities of 7 enzymes: Superoxide Dismutase (SOD), Catalase (CAT), Acid Phosphatase (ACP), Alkaline Phosphatase (AKP), Hexokinase (HK), Pyruvate Kinase (PK), Lactate Dehydrogenase (LDH), and 2 lipid metabolism indicators: Total Cholesterol (TCHO) and Triglycerides (TG). All kits used for determination were purchased from Nanjing Jiancheng Bioengineering Institute, and experimental procedures were performed strictly according to the kit instructions.
2.4. RNA Extraction and Transcriptome Sequencing
Hepatopancreas tissue samples were collected, and total RNA was extracted using the TRIzol method, following the kit instructions. After extraction, multiple tests were performed on RNA quality: RNA purity was determined using a NanoPhotometer spectrophotometer, assessed by the OD260/OD280 ratio (1.8~2.0 as qualified) and OD260/OD230 ratio (≥ 2.0 as qualified); RNA concentration was precisely quantified using a Qubit 2.0 Fluorometer; RNA integrity was detected using an Agilent 2100 Bioanalyzer (RIN value ≥ 7.0 as qualified).
The transcriptome sequencing library construction process was as follows: qualified long-chain RNA was fragmented to the target length using Fragmentase® enzymatic digestion to complete RNA fragmentation; the first strand of cDNA was synthesized using reverse transcriptase, then using the first strand of cDNA as a template, the second strand was synthesized by replacement synthesis to obtain the double-stranded cDNA basic product; the double-stranded cDNA was subjected to end repair and Klenow fragment A-tailing; TrueSeq® UDI adapter ligation method was used for adapter ligation, while adding Unique Molecular Identifier (UMI), and the adapter concentration was gradient optimized to inhibit dimer formation; finally, the constructed library was quality tested and quantified, and qualified libraries were sent to the sequencing platform to complete sequencing. The sequencing work was commissioned by Guangzhou GeneDenovo Biotechnology Co., Ltd.
2.5. Gene Expression and Differential Expression Analysis
Gene expression levels were quantified using StringTie v2.2.1 and calculated as fragments per kilobase of transcript per million mapped reads (FPKM) based on transcript assembly and the estimation of mapped read abundance8. Differential expression analysis was performed using the DESeq2 v1.34.0 R package9. Specifically, library size differences were corrected by size factor normalization, read counts were fitted to a negative binomial generalized linear model, and significance was assessed using the Wald test. Differentially expressed genes (DEGs) were identified based on the criteria of |log2(Fold Change)| > 1 and an adjusted P-value < 0.05 (Benjamini-Hochberg method)10. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were conducted using the clusterProfiler v4.2.2 R package. The GO analysis used the hypergeometric test to identify significantly enriched GO terms among DEGs, while the KEGG pathway analysis used KEGG ortholog annotations to identify significantly over-represented biological pathways. A Benjamini-Hochberg corrected false discovery rate (FDR) < 0.05 was set as the significance threshold for both enrichment analyses.
2.6. Statistical analysis
Statistical analyses were performed using IBM SPSS Statistics 26. Prior to analysis, data were tested for normality and homogeneity of variance to ensure the assumptions were met. A two-way analysis of variance (ANOVA) was employed to evaluate the effects of copper concentration, stress duration, and their interaction on the physiological indices. When significant main effects or interactions were detected (p < 0.05), Tukey’s honestly significant difference (HSD) test was conducted for post-hoc multiple comparisons to identify specific differences between groups. All data are presented as mean ± standard deviation (SD), and the level of statistical significance was set at p < 0.05.
3. Results
3.1. Enzyme Activities and Lipid Metabolite Results
As shown in Figure 1A, the LDH activity in the hepatopancreas of clams in both the high concentration group (0.568 mg/L) and the low concentration group (0.0568 mg/L) reached a maximum value on Day 2, and the difference between the two treatment groups and the control group was significant (p < 0.05). During the period from Day 2 to Day 7, the LDH content in the hepatopancreas of the two concentration groups showed a downward trend and was significantly different from the control group (p < 0.05). By Day 14, it had returned to a level close to that of the control group, and the difference between the two concentration groups and the control group was not significant (p > 0.05).
As shown in Figure 1B, the HK activity in the hepatopancreas of clams in both the high concentration group (0.568 mg/L) and the low concentration group (0.0568 mg/L) reached a maximum value on Day 2, which was significantly different from the control group (p < 0.05). It showed a linear downward trend from Day 2 to Day 7, with a significant difference from the control group (p < 0.05). From Day 7 to Day 14, the difference between the low concentration group (0.0568 mg/L) and the high concentration group (0.568 mg/L) compared with the control group was significant (p < 0.05).
As shown in Figure 1C, the PK activity in the hepatopancreas of clams in both the low-concentration (0.0568 mg/L) and high-concentration (0.568 mg/L) groups peaked on day 2, showing significant differences compared to the blank control group (p < 0.05). On day 7, the low-concentration group (0.0568 mg/L) still exhibited a significant difference from the blank control group (p < 0.05). In contrast, the PK activity in the high-concentration group (0.568 mg/L) decreased to near control levels on day 7, with no significant difference observed (p > 0.05). By day 14 of stress, neither concentration group showed significant differences compared to the blank control group (p > 0.05).
As shown in Figure 2A, SOD activity in the hepatopancreas of clams in the high concentration group (0.568 mg/L) and the low concentration group (0.0568 mg/L) increased on Day 2 of stress, with significant differences between groups (p < 0.05). On Day 7 of exposure, the two concentration groups were significantly different from the blank control group (p < 0.05). From Day 7 to Day 14, the low concentration group (0.0568 mg/L) was significantly different from the blank control group (p < 0.05), while the high concentration group (0.568 mg/L) was not significantly different from the blank control group (p > 0.05).
As shown in Figure 2B, in both the high-concentration (0.568 mg/L) and low-concentration (0.0568 mg/L) groups, CAT activity in the hepatopancreas of clams increased significantly on day 2 of stress compared to the blank control group (p < 0.05). Both groups reached their peak values on day 7, after which a downward trend was observed; notably, CAT activity in both groups remained significantly higher than that of the control group (p < 0.05). On day 14 of stress, both concentration groups still exhibited significant differences relative to the control group (p < 0.05).
As shown in Figure 3A, ACP activity in the hepatopancreas of clams in both the high-concentration (0.568 mg/L) and low-concentration (0.0568 mg/L) groups reached its lowest value on day 2, which was significantly lower than that of the control group (p < 0.05). Subsequently, ACP activity in both groups increased significantly from day 2 to day 14, with differences from the control group throughout this period (p < 0.05).
As shown in Figure 3B, the AKP activity in the hepatopancreas of clams in the high concentration group (0.568 mg/L) reached a maximum value on Day 2, while that in the low concentration group (0.0568 mg/L) reached a minimum value on Day 2. Both concentration groups were significantly different from the control group (p < 0.05). Both concentration groups showed an overall downward trend from Day 2 to Day 7. Overall, the AKP activity in the high concentration group (0.568 mg/L) was higher than that in the low concentration group (0.0568 mg/L), and both concentration groups were significantly different from the control group (p < 0.05). On Day 14, the difference between the high concentration group and the control group was not significant (p > 0.05), whereas that between the low concentration group and the control group was significant (p < 0.05).
As shown in Figure 4A, the TG activity in the hepatopancreas of clams in both the high concentration group (0.568 mg/L) and the low concentration group (0.0568 mg/L) experienced a sharp drop on Day 2, with a significant difference (p < 0.05). It rebounded from Day 2 to Day 7, returning to normal levels on Day 7; the difference between the two concentration groups and the control group was not significant (p > 0.05). The change from Day 7 to Day 14 was relatively gentle, with no significant difference from the control group (p > 0.05).
As shown in Figure 4B, the TCHO content in the hepatopancreas of clams in both the high concentration group (0.568 mg/L) and the low concentration group (0.0568 mg/L) was similar to the TG content change, both dropping to the lowest point on Day 2, with a significant change compared to the control group (p < 0.05). Between Day 2 and Day 7, both concentration groups showed an upward trend, and the difference between the two concentration groups and the control group was significant (p < 0.05). Between Day 7 and Day 14, they tended to be stable and were significantly different from the control group (p < 0.05).
3.2. Integrated Biomarker Response
After the hepatopancreas of Ruditapes philippinarum was exposed to two copper concentrations (0.0568 mg/L, 0.568 mg/L) for 2, 7, and 14 days, the star chart (Figure 5) integrally displayed the transformation data of all biomarkers. The Integrated Biomarker Response (IBR) index was positively correlated with biological toxic stress, with higher IBR values indicating greater toxic stress pressure.The results showed that on Day 2 of exposure to the two copper concentrations, the biomarkers in the hepatopancreas of Ruditapes philippinarum generally indicated Lactate Dehydrogenase (LDH), Pyruvate Kinase (PK), and Hexokinase (HK). This suggests that under stress, the organism meets emergency energy supply needs by maintaining or enhancing glycolysis pathways. The upregulation of the activities of these enzymes is a key response for cells to meet energy and biosynthesis requirements. On Day 7 of exposure, the biomarkers generally indicated Superoxide Dismutase (SOD) and Catalase (CAT), suggesting that, under continuous copper stress, clams have initiated an oxidative damage defense system to cope with oxidative stress. On Day 14 of exposure, the biomarkers of the two concentration groups generally pointed to Acid Phosphatase (ACP), Triglycerides (TG), and Total Cholesterol (TCHO), suggesting that under long-term stress, clams’substance metabolism and cell damage-related pathways have undergone significant changes.
3.3. Sequencing Data Quality Assessment
The raw data of transcriptome sequencing were quality-controlled using fastp software to obtain high-quality clean reads for subsequent analysis. Specific quality control standards were as follows: discard reads containing adapter sequences; remove reads where unknown bases (N) accounted for more than 10%; discard reads consisting entirely of adenine (A) bases; filter out low-quality reads, defined as sequences where bases with a quality value Q ≤ 20 accounted for 50% or more of the entire read length (Table 1).
3.4. Analysis of Differentially Expressed Genes in Ruditapes philippinarum
Figure 6 shows the volcano plot comparing differentially expressed genes (DEGs) between the 0.0568 mg/L concentration group and the control group. Red dots represent significantly upregulated differentially expressed genes, blue dots represent significantly downregulated differentially expressed genes, and black dots represent genes with no significant difference. On Day 2 of exposure (D2), 432 DEGs were identified, including 408 significantly upregulated and 24 significantly downregulated. On Day 7 of exposure (D7), 439 DEGs were identified, including 407 significantly upregulated and 32 significantly downregulated. On Day 14 of exposure (D14), there were 5,337 DEGs, including 5,088 significantly upregulated and 249 significantly downregulated. The total number of differentially expressed genes increased with increasing exposure time, peaking on Day 14, indicating that clams’ gene expression was most sensitive to copper stress at this time.
Figure 7 shows the volcano plot comparing differentially expressed genes between the 0.568 mg/L concentration group and the control group. On Day 2 of exposure (D2), there were 84 DEGs, including 53 significantly upregulated and 31 significantly downregulated. On Day 7 of exposure (D7), there were 4,325 DEGs, including 4,108 significantly upregulated and 217 significantly downregulated. On Day 14 of exposure (D14), there were 4,145 DEGs, including 3,850 significantly upregulated and 295 significantly downregulated. The total number of differentially expressed genes showed an upward trend with the extension of exposure time, peaking on Day 14, indicating that clams’ gene expression was most sensitive to copper stress at this time.
3.5. GO Classification of Differentially Expressed Genes
As shown in Figures 8, in the GO functional enrichment analysis of the low concentration copper exposure group (0.0568 mg/L), the gene functions at different exposure time points presented obvious temporal response characteristics:
On Day 2 of exposure, the organism’s core was rapid stress response. In terms of molecular function, a large number of genes were enriched in transmembrane receptor activity, ion channel activity, and transporter activity, which are related to signal perception and material transport. Combined with the high gene enrichment in biological processes such as cellular process, metabolic process, and response to stimulus, this indicates that clams rapidly initiated transmembrane signal recognition, material transport, and basic metabolic regulation in the early stage of stress, providing a molecular basis for subsequent energy metabolism conversion.
On Day 7 of exposure, the organism entered the defense and homeostasis maintenance stage. In terms of molecular function, new enrichments in antioxidant activity and protein folding chaperone functions were added, reflecting that the organism initiated defense mechanisms against copper-induced oxidative damage while activating pathways related to protein homeostasis maintenance. Although metabolic and regulatory processes in biological processes continued to be enriched, the number of genes fell back, reflecting the adjustment from rapid initiation of stress response to precise defense adaptation.
On Day 14 of exposure, the organism was in a stage of maintaining homeostasis and multi-dimensional adaptation. The number of enriched genes in each GO category increased significantly. In cellular components, the number of genes related to structural functions such as cell anatomical entities and protein complexes rose sharply. In biological processes, immune system processes and reproductive-related pathways were newly enriched. This indicates that under long-term stress, clams not only continued to strengthen basic metabolism and regulation but also initiated adaptive mechanisms such as immune defense and tissue regeneration to maintain organism homeostasis.
As shown in Figures 9, in the GO functional enrichment analysis of the high concentration copper exposure group (0.568 mg/L), the gene functions at different exposure time points also presented obvious temporal response characteristics:
On Day 2 of exposure, the core was rapid stress defense. In terms of molecular function, functions related to material transport, transcriptional and translational regulation, and signal transduction were rapidly enriched, reflecting that the organism immediately initiated transmembrane material transport and gene expression regulation under high concentration stress. In biological processes, the immune response, interspecies interactions, and positive/negative regulatory pathways were simultaneously activated. Combined with the enrichment of protein complexes and cell anatomical entities in cellular components, this indicates that in the early stage of stress, the organism simultaneously mobilized defense, metabolism, and structural maintenance responses to cope with the acute stress of high concentration copper.
On Day 7 of exposure, the organism entered the stage of metabolic homeostasis and somatic repair. In molecular functions, catalytic activity and binding ability continued to be enriched, and the number of genes increased significantly. The structural support function in cellular components remained stable. In biological processes, cellular process and metabolic process became the core enriched pathways, while immune-related pathways temporarily fell back. This indicates that the organism shifted from acute defense to metabolic adaptation, repairing damage and maintaining energy supply by strengthening basic metabolism and cellular functions to meet the continuous energy demand under high-concentration stress.
On Day 14 of exposure, the organism presented homeostasis maintenance and detailed regulation. Similar to the low concentration group, the number of enriched genes in each GO category further increased. Core functions in molecular functions and cellular components remained stable. In biological processes, immune response pathways were reactivated, and positive/negative regulation pathways were added to achieve more precise metabolic and immune balance. The enrichment of multicellular organismal processes suggested that the organism initiated tissue regeneration and long-term adaptation mechanisms, indicating that Ruditapes philippinarum under long-term high-concentration stress had constructed a system of metabolic-immune defense-protein regulation.
3.6. KEGG Enrichment Analysis
As shown in Figures 10, in the KEGG pathway enrichment analysis of the low concentration copper exposure group (0.0568 mg/L), the activation of pathways at different exposure time points presented a clear temporal evolutionary characteristic:
On Day 2 of exposure (D2), it was in the stage of initial stress and metabolic adjustment. Pathway enrichment was centered on fatty acid metabolism pathways and cellular stress defense, accompanied by the initiation of initial immune responses. This indicates that, in the early stage of stress, clams adjusted energy metabolism and activated basic defense mechanisms to cope with acute copper ion stimulation.
On Day 7 of exposure (D7), it was in the stage of immune defense and metabolic adaptation. Pathway enrichment turned towards immune metabolism related directions. The lysosome pathway continued to be activated, while carbohydrate and lipid metabolism adjustment pathways were newly added. Inflammation-related pathways were also simultaneously activated, reflecting that the organism shifted from initial stress to a response mode with immune defense as the focus accompanied by fine metabolic adaptation.
On Day 14 of exposure (D14), it was in the stage of comprehensive homeostasis remodeling. Core metabolic networks, xenobiotic detoxification systems, and oxidative stress defense pathways were widely activated. Ruditapes philippinarum built a more stable intracellular homeostasis through multi-dimensional metabolic remodeling and protein synthesis remodeling to adapt to the long-term copper stress environment.
As shown in Figures 11, in the KEGG pathway enrichment analysis of the high concentration copper exposure group (0.568 mg/L), the activation of pathways at different exposure time points presented temporal response characteristics different from those of the low concentration group:
On Day 2 of exposure (D2), it was in the stage of neuro-metabolic stress and multi-system disturbance. Exposure to high concentration early on triggered multi-system disturbances related to neuro-metabolic stress. Neuroprotective pathways, amino acid and lipid metabolism adjustment pathways, as well as cell autophagy pathways were simultaneously annotated, indicating that high concentration copper stress in the early stage had a rapid and broad impact on the nervous system and energy metabolism network of clams.
On Day 7 of exposure (D7), it was in the stage of large-scale metabolic activation and strong detoxification. Core metabolic pathways were comprehensively activated, while the enrichment of xenobiotic detoxification systems, cell damage clearance pathways, and ribosomal protein synthesis-related pathways was significantly enhanced. This reflects that under continuous high concentration stress, the organism coped with toxic pressure by strengthening metabolic supply and damage repair.
On Day 14 of exposure (D14), it was in the stage of steady-state defense and oxidative stress response. Core metabolic pathways and detoxification systems remained in an activated state, and new enrichments in oxidative stress defense pathways appeared. This indicates that under long-term high-concentration stress, clams stabilized the intracellular environment by maintaining basic metabolic supply, continuous detoxification, and enhanced defense against oxidative damage.
4. Discussion
4.1. Energy metabolic shift and compensatory mechanism
During the early phase of copper stress (D2), Ruditapes philippinarum relied on a synergistic glycolysis–fatty acid cycle to cope with acute stress. Transcriptomic data showed that the fatty acid metabolism pathway was significantly enriched in the low-concentration group, indicating activation of fatty acid oxidation. Concurrently, activities of HK, PK and LDH peaked on day 2, reflecting an emergency activation of the glycolytic pathway. Such metabolic reprogramming can maintain ATP supply while alleviating copper toxicity 11, compensating for the low efficiency of glycolysis alone. Compared with previous studies on copper-stressed Pacific oysters, our results demonstrate that fatty acid metabolism is activated as early as the initial stage, highlighting differences in metabolic adaptation to heavy metal stress among bivalves 12. By contrast, the high-concentration group (0.568 mg/L) exhibited an early overload response, characterized by concurrent activation of lipid and amino acid metabolism pathways, with lipid metabolism being markedly more pronounced than in the low-concentration group.
In the mid-stage of copper stress (Day 7), the organism maintained homeostasis through oxidative defense mechanisms, with the activities of superoxide dismutase (SOD) and catalase (CAT) peaking on Day 7. Five unique reductase gene families：PNDR, hydroxylamine reductase/Ni-CO dehydrogenase, VDE, PPO, and Rieske-type—were identified in Ruditapes philippinarum, suggesting their role in antioxidant defense against ROS-induced oxidative stress13.In the low-concentration groups, pathways related to autophagy and antigen metabolism were enriched. The autophagy system facilitated the degradation of damaged components to recycle amino acids, thereby providing energy for the synthesis of antioxidant enzymes14. In contrast, the high-concentration groups entered a phase of extensive metabolic activation, during which core metabolic pathways, such as carbohydrate and lipid metabolism, were fully activated to provide sufficient ATP for the catalytic reactions of antioxidant enzymes.
In the late stage of stress (Day 14), core metabolic pathways were annotated in both the high- and low-concentration groups. This indicates that copper toxicity had been effectively alleviated, eliminating the need for the organism to maintain a state of high-intensity metabolic activation; instead, it shifted toward establishing a stable energy supply system. It has been reported that virus-infected Ruditapes philippinarum produce two lipid metabolites, diacylglycerol (DG) and triacylglycerol (TG), which facilitate lipid regulation and provide long-term energy reserves for the shellfish15. This finding is consistent with the variations in TG content observed in the present study. Furthermore, prolonged copper stress disrupted cholesterol homeostasis and impaired cell membrane structure and function. Consequently, the organism upregulated total cholesterol synthesis to repair the damage, thereby further safeguarding intracellular homeostasis.
4.2. Cellular Autophagy and Immune Defense System
In the early stage of stress (Day 2), GO enrichment analysis revealed that both concentration groups prioritized a rapid stress response. The autophagy pathway was annotated in the low-concentration group. The autophagy pathway recognizes specific carbohydrate molecules to tag damaged components or pathogen-associated molecules, providing a basis for targeted degradation by the autophagic machinery16. In terms of energy coordination, amino acids recycled via autophagy entered the gluconeogenesis pathway, while ATP rapidly generated by glycolysis reciprocally fueled autophagic activity, thereby forming an efficient energy feedback loop. The high-concentration group showed activation of immune warning pathways, including mitophagy and atherosclerosis, in addition to autophagy. The response intensity was significantly higher than that in the low-concentration group.
In the mid-stage of stress (Day 7), the lysosomal pathway was annotated in both concentration groups. During this phase, ROS levels peaked, the autophagy system sequesters damaged mitochondria and misfolded proteins into autophagosomes, which subsequently fuse with lysosomes for degradation. This process not only eliminates ROS sources and oxidatively damaged components but also recycles amino acids to replenish energy, thereby maintaining energy homeostasis and stabilizing SOD and CAT activities. Furthermore, experimental findings revealed that the lysosomal pathway serves functions beyond degradation; it facilitates the sequestration of copper complexes to mitigate free copper toxicity and processes antigenic peptides for presentation to the immune system, thereby initiating a specific immune response17.
In the late stage of stress (Day 14), GO enrichment results formally identified the addition of “immune system processes.” Concurrently, genes related to structural components, such as cellular anatomical entities and protein complexes, were significantly upregulated, indicating that the organism had shifted from a passive damage response to active maintenance of homeostasis. Both concentration groups showed persistent enrichment in the lysosome pathway, hematopoietic cell lineage pathway, and glutathione metabolism pathway. During this phase, lysosomes continued to degrade damaged proteins and release amino acids; ribosomes then utilized these recycled amino acids to preferentially synthesize stress-related functional proteins, such as metallothioneins, heat shock proteins, and antioxidant enzymes, thereby re-establishing proteostasis. ACP and AKP, two essential lysosomal hydrolases, play pivotal roles in disease resistance and non-specific immune responses.18. Studies on shellfish subjected to temperature variation reported that the activation of the immune defense system and the peak of ACP activity lagged behind those of SOD and CAT, coinciding with the activation of the lysosome pathway19. This temporal pattern is consistent with the findings of the present study.
4.3. Protein Homeostasis and Heavy Metal Detoxification
In the early stage of copper stress (D2), the prion disease pathway was significantly enriched in the high-concentration group but not in the low-concentration group, likely because high copper concentrations exacerbate protein damage via thiol oxidation20. In contrast, the enrichment of the ABC transporter pathway in the low-concentration group indicated that active efflux of copper ions was the primary detoxification mechanism21. Furthermore, the glycine, serine, and threonine metabolism pathway was enriched in the high-concentration group, suggesting that amino acid chelation of copper ions mitigated toxicity22.
In the mid-stage of stress (Day7), the cytochrome P450-mediated xenobiotic metabolism pathway was newly enriched in the high-concentration group. The CYP450 enzyme system facilitates the enzymatic transformation and detoxification of toxic substances23. In contrast to the low-concentration group, which achieved detoxification primarily through lipid chelation and lysosomal activity, the high-concentration group necessitated the activation of CYP450-mediated enzymatic transformation mechanisms. This strategy prevented the continuous intracellular accumulation of toxic substances that could otherwise induce apoptosis or necrosis, further corroborating the concentration-dependent nature of the detoxification response.
In the late stage of copper stress (Day14), the lysosome, phagosome, ribosome, and glutathione metabolism pathways remained continuously enriched, indicating that the organism maintained homeostasis by degrading damaged proteins to synthesize stress proteins and scavenge ROS24,25. The enrichment of the cholesterol metabolism pathway in the high-concentration group suggested that the organism restored cell membrane fluidity to minimize copper ion influx26. In contrast, the low-concentration group maintained activation of the ABC transporter pathway, resulting in complete elimination of toxicity.
Acknowledgments
The study was supported by Foundation National Key R&D Program of China (2023YFD2400800), Earmarked Fund for CARS (CARS-49), Tianjin Fisheries Green Development Project (YYLSFZ202501).
Authors’ Contribution
Conceptualization: Chen Li (Lead). Methodology: Chen Li (Lead). Investigation: Chen Li (Equal), Yucheng Liu (Equal), Xiaoyu Wang (Equal). Data curation: Chen Li (Equal), Yucheng Liu (Equal). Writing – original draft: Chen Li (Lead). Writing – review & editing: Chen Li (Equal), Shuang Liang (Equal), Yongren Li (Equal), Dexian Gu (Equal). Funding acquisition: Yongjun Guo (Lead). Formal Analysis: Xiaoyu Wang (Lead). Resources: Xiaoyu Wang (Lead).
Competing of Interest – COPE
No competing interests were disclosed.
Ethical Conduct Approval – IACUC
This experiment was strictly conducted in accordance with the regulations and operational guidelines formulated by the Animal Care and Use Committee of Tianjin Agricultural University, and formal approval was obtained from the committee.
Informed Consent Statement
All authors and institutions have confirmed this manuscript for publication.
Data Availability Statement
All are available upon reasonable request.