We collected 120 adult mussels M. galloprovincialis with a length of 5.5 ± 0.3 cm from a pristine environment along Yangma Island in Yantai, China, and transported them to our laboratory immediately. These mussels were acclimated in aerated seawater (20°C ± 2.0°C, 31.1 psu) for 7 days. Subsequently, we randomly sorted the mussels into six groups: five DEHP (99.5%, Dr. Ehrenstorfer, Augsburg, Germany) treatment groups and one blank control. Each group had two replicate glass tanks. Each glass tank contained 10 mussels in 10 L of aerated seawater. Five different concentrations (4, 12, 36, 108, and 324 μg/L) of DEHP were used to study the effect of dose-dependent exposure on the mussels. Of these five experimental concentrations, the three experimental concentrations of 4, 12, and 36 μg/L were environmentally relevant (Gao et al., 2015). We used dimethyl sulfoxide (DMSO) as a solvent to dissolve DEHP. Studies on the mussel M. galloprovincialis have reported no significant biological difference between the solvent control group and the seawater control group (Ji et al., 2013). Therefore, we only had one solvent control group, and the concentrations of DMSO (0.00025%, v/v) in the solvent control and DEHP-treated groups were identical. In this work, Chlorella vulgaris was used to feed the mussels (approximately 2% of the tissue per dry weight daily) during the acclimation and exposure periods. After 14 days of exposure, eight mussels were randomly selected from each treatment, and the entire soft tissue was dissected. Subsequently, we quickly froze each sample in liquid nitrogen. Finally, each sample was divided into three parts before conducting enzyme activity determination, RNA extraction, and metabolite extraction.

The activities of antioxidant enzymes in mussel tissues (n = 6) were assayed by using commercially available enzyme kits (Jiancheng, Nanjing, China). Furthermore, we measured the levels of three necessary antioxidant enzymes: catalase (CAT; EC 1.11.1.6), glutathione S-transferase (GST; EC 2.5.1.18), and superoxide dismutase (SOD; EC 1.15.1.1). A bicinchoninic acid assay was employed to determine the protein content in mussel tissues.

Total RNA was extracted from the whole soft tissue samples (n = 6) according to the manufacturer’s instructions. First-strand cDNA was synthesized based on information in the M-MLV reverse transcription kit (Promega, Madison, WI, United States). For the analysis of gene expression levels, six stress responsive genes, including SOD, CAT, GST, nuclear factor-k-gene binding (NF-κB), goose-type lysozymes (MgGLYZ), and heat shock protein 70 (HSP70), were selected, and the beta-actin gene (β-actin) was used as the reference gene. All the primers of these genes are listed in Table 1.

The sample cDNA was mixed, and fourfold serial dilution was performed. The diluted cDNA was used as a template for the fluorescent real-time quantitative polymerase chain reaction (qRT-PCR), with three replicates for each template. The fluorescent qRT-PCR reaction system was operated according to the manufacturer’s instructions (Liu et al., 2019). The value during the amplification process was continuously measured in real time to establish a standard curve. The obtained data were analyzed using the ABI 7500 SDS software (Applied Biosystems). The 2^–ΔΔCT method was used to analyze the gene expression levels (Livak and Schmittgen, 2001). The qRT-PCR process is detailed in the Supplementary Material.

Metabolites were extracted using the modified extraction protocol described by Ji et al. (2020). Briefly, the metabolites in each M. galloprovincialis sample (100 mg wet weight) were extracted by 485 μL of methanol/water (4/0.85, v/v), 200 μL of chloroform, and 225 μL of water. The mixture was first vortexed thoroughly and left to stand on ice for 10 min, and it was then centrifuged at 3,000 r/min for 10 min. Finally, the upper polar phase was absorbed and vacuum-dried for testing. We prepared 600 μL of phosphate buffer solution (100 mmol/L, pH = 7.0) in heavy water before adding the solution to a numbered EP sample tube. The tube was vortexed thoroughly for its contents to dissolve before being centrifuged at 3,000 r/min for 10 min, and 550 μL of supernatant was pipetted into a nuclear magnetic resonance (NMR) tube (Φ = 5 mm) for further NMR analysis.

¹H NMR spectral analysis was conducted using a Bruker AV 500 NMR spectrometer with a 5-mm reverse probe (Wu et al., 2013a,b). The corresponding proton resonance frequency was 500.18 MHz. The other parameters were as follows: spectrum width = 6009.6 Hz, delay time = 0.3 s, mixing time = 0.1 s, line broadening factor = 1, and number of scans = 128. Moreover, the water peak was suppressed by presaturation irradiation during the relaxation delay and mixing period when the ¹H NMR spectrum was obtained (Wu et al., 2013a).

TopSpin 2.1 was used to perform phase correction and baseline adjustment on the spectrum, and trimethylsilyl propanoic acid was used as the internal standard of the chemical shift reference peak (0.0 ppm) (Viant et al., 2009). The processing and analysis of NMR data are detailed in the Supplementary Material. In addition, Chenomx (Evaluation Version, Chenomx, Edmonton, Alberta, Canada) was used to identify major metabolites, specifically by analyzing the chemical shifts of metabolites (Viant et al., 2009; Wu et al., 2013c; Li et al., 2015).

Data analysis was performed using Minitab (Version 15, Minitab, United States). The Shapiro–Wilk test was used to verify the normality of the distribution for all variables. These variables were antioxidant enzyme activities and gene expression levels in the mussel samples, which were subject to Ryan–Joiner’s test and Bartlett’s test, respectively (Liu et al., 2019). The values are presented in terms of the mean ± standard deviation for continuous variables. One-way analysis of variance (ANOVA) with Tukey’s test was performed to determine the differences between the control and DEHP treatments. A P-value less than 0.05 was regarded as statistically significant.

The dose-dependent effects of environmental pollutants depend on the testing endpoints (Agathokleous et al., 2018). First, the activities of three typical antioxidant enzymes, namely CAT, SOD, and GST, were measured. As summarized in Figure 1, both SOD and CAT activities in mussel samples demonstrated typical inverted U-shaped curves with increasing DEHP concentration, which indicated the presence of DEHP-induced hormetic effects in mussels. Specifically, the environmentally realistic concentration of 12 μg/L induced significant (P < 0.05) increases in CAT and SOD activities, whereas the higher concentrations of 108 and 324 μg/L significantly (P < 0.05) inhibited the activities of these two antioxidant enzymes.

As an antioxidant enzyme, SOD can catalyze the dismutation of superoxide radicals into oxygen (O₂) or hydrogen peroxide (H₂O₂) in coordination with CAT, consequently preventing reactive oxygen radicals from damaging organisms. Therefore, these two enzymes play critical roles in protecting cells from the oxidative stress induced by reactive oxygen species. In this study, the increased activities of SOD and CAT revealed that the mussels enhanced the capacity of these two enzymes against the potential oxidative stress induced by DEHP. However, the inhibited activities of SOD and CAT indicated interference to antioxidant enzymes induced by high concentrations (108 and 324 μg/L) of DEHP, which was similarly observed in the DEHP-treated yellow catfish P. fulvidraco (Yuan et al., 2017). Notably, the hormetic effects of cadmium were also detected in both earthworm Eisenia fetida and polychaete P. aibuhitensis when CAT and SOD activities were observed (Zhang et al., 2009; Liu et al., 2019). As a known phase II metabolic isozyme, GST is the key antioxidant enzyme that catalyzes the conjugation of glutathione to several electrophilic substrates in organisms. Basically, it functions in the coupling of the electrophilic group of harmful substances with the sulfhydryl group of the reduced form of glutathione for the detoxification of toxicants (Chen et al., 2018). We noted no significant DEHP-induced hormetic or monophasic response in GST activity in the mussels. However, GST activity was significantly (P < 0.01) inhibited in all DEHP-treated groups, suggesting that DEHP might produce damage to proteins and subsequent enzymatic inactivation of GST and reduced the detoxification capacity of mussels.

Figure 1 shows the dose-dependent effects of DEHP at the mRNA level when we evaluated the levels of six stress-responsive genes (CAT, SOD, GST, NF-kB, MgGLYZ, and HSP70). In general, these genes were responsive to oxidative stress (SOD, CAT, and GST) and immune stress (NF-kB, MgGLYZ, and HSP70) (Cellura et al., 2006; Wang et al., 2012; Sang et al., 2020). After exposure to 36 μg/L of DEHP treatment for 14 days, the expression level of SOD was significantly (P < 0.05) downregulated. For both CAT and GST, their expression profiles exhibited U-shaped curves (Figure 1), which were intriguingly the inverse of the response profiles (inverted U-shaped curves) of CAT and GST activities (Figure 2). This constituted a major finding in understanding how mRNA expression is encoded protein’s tendency, which does not always happen due to post-translational or post-transcriptional modifications. Additionally, the enzyme activities contained the total activities of the entire enzyme family. Therefore, the disparity of the gene expressions of CAT and GST with the corresponding CAT and GST activities was not unexpected (Wu et al., 2013b).

NF-κB is an evolutionarily conserved family of regulatory molecule that is widely present in various tissues and plays a vital role in regulating the expression of diverse immune-related genes that respond to pathogens and the immune system (Hayden and Ghosh, 2008). Lysozymes are a non-specific immune factor that can hydrolyze the peptidoglycan of the bacterial cell wall, resulting in cell wall rupture and content escape. Furthermore, lysozyme scan inactivates pathogens by binding to the negatively charged viral protein, and thus, they play critical roles in the defense and elimination of pathogens (Wang et al., 2012). HSP70 is a molecular chaperone belonging to the HSP family and with standard features in the aggregation of unfolded or misfolded proteins and in the degradation of denatured protein (Wang et al., 2012). In this study, NF-κB expression levels were substantially and significantly (P < 0.01) upregulated in the groups treated with 108 and 324 μg/L of DEHP, suggesting that these two high doses of DEHP induced immune stress in mussels. Furthermore, the expression profile of MgGLYZ exhibited an approximately U-shaped curve with increasing dose, and we noted a significantly decreased expression level of MgGLYZ in the group treated with 12 μg/L of DEHP. This finding is consistent with our previous finding on CAT and GST expression profiles and on the presence of DEHP-induced hormesis in mussels at the mRNA level. In all DEHP treatments, HSP70 expression levels were significantly (P < 0.01) downregulated without linear or hormetic dose responses. Future studies should conduct transcriptomic analysis to globally profile the dose-dependent responses induced by DEHP at the mRNA level.

Metabolomic analysis can characterize the alteration of biological samples, allowing the researcher to identify perturbed metabolic pathways and obtain vital information pertaining to biochemical networks (Viant et al., 2009; Wu et al., 2012). Therefore, metabolomics has been extensively applied in marine ecotoxicology due to its capability to globally characterize the metabolic responses in organisms stressed by environmental pollutants, which can elucidate the toxicological effects and select biomarkers (Cappello et al., 2013, 2015; Ji et al., 2015; Yu et al., 2016; Liu et al., 2019). In our study, metabolomic analysis based on the NMR technique was employed to characterize the dose-dependent effects of DEHP in mussels at the metabolite level. As shown in Figure 3, we identified 34 metabolites in mussel tissue extracts, including energy storage compounds (glycogen, ATP, ADP, and glucose), amino acids (arginine, leucine, alanine, isoleucine, valine, glutamate, glycine, and tyrosine), organic osmolytes (betaine, hypotaurine, taurine, dimethylglycine, and homarine), and intermediates in the Krebs cycle (fumarate and succinate). Principal component analysis was conducted on the NMR data from the control group and all DEHP treatment groups. In general, the average scores along PC1 (46.38% variation) decreased linearly (Figure 4A), presenting a monophasic dose–response relationship along the PC1 axis. Although the lowest concentration (4 μg/L) of DEHP treatment did not significantly differ from the control, the three groups exposed to higher DEHP concentrations (36, 108, and 324 μg/L) significantly (P < 0.05) differed from the control group along the PC1 axis. However, the average scores along the PC4 (4.48% variation) axis had an inverted U-shaped curve with two moderate concentrations (12 and 36 μg/L) of DEHP treatments located at the top (Figure 4B), which was consistent with the dose–response curves for SOD and CAT activities. These findings confirmed that DEHP mainly induced a linear dose-dependent response at the metabolite level, accounting for 46.38% of the total variation induced by DEHP treatments. However, we also found significant hormetic responses, as indicated by the inverted U-shaped curve, verifying that at environmentally relevant concentrations, DEHP induced hormesis in mussels.

As detailed in Figure 5, orthogonal projection to latent structures discriminant analysis (O-PLS-DA) was conducted. At the concentrations of 12 and 324 μg/L, DEHP exhibited a typical biphasic dose–response relationship in mussels. In mussel samples treated with 12 μg/L of DEHP, 12 metabolites, including 10 amino acids (such as threonine, asparagine, and lysine) and two energy storage metabolites (glycogen and glucose), were considerably (P < 0.05) altered. Of note, in the group treated with 12 μg/L of DEHP, these amino acids exhibited obviously reduced levels, whereas glucose and glycogen had markedly increased levels. In DEHP treatment at the highest concentration (324 μg/L), both arginine and lysine levels decreased, which was similar to what happened under the 12 μg/L DEHP treatment. However, glucose, glycogen, and tyrosine were, conversely, altered in the highest concentration of the DEHP-treated group. In addition, levels of glutamate, aspartate, homarine, and inosine were significantly elevated, whereas taurine was depleted in the group receiving the highest concentration (324 μg/L) of DEHP. Among these metabolites, the average concentration of homarine gradually increased with the increasing concentration of DEHP (from 4 to 324 μg/L), exhibiting a monophasic response. Other metabolites (especially amino acids, including isoleucine, threonine, alanine, and asparagine) recovered to normal levels in the control group, exhibiting DEHP-induced hormetic responses as well.

In marine invertebrates, amino acids play essential roles in energy metabolism and osmotic regulation (Viant et al., 2003; Bonnefille et al., 2018; Dumas et al., 2020). In the 12 μg/L DEHP treatment, the levels of amino acids were approximately negatively correlated with the contents of glucose and glycogen. These results illustrated that DEHP exposure (12 μg/L) enhanced the consumption of amino acids and reduced the consumption of glucose for energy supply. This suggested that mussels elevated their energy storage to cope with the DEHP-induced stress, indicating an adaptive mechanism to DEHP exposure. However, in the group exposed to the highest concentration (324 μg/L) of DEHP, depleted glucose and glycogen, as well as the recovered amino acids, meant that the mussels consumed more energy to survive under this high concentration of DEHP. This hormetic response could be explained by the overcompensation of energy metabolism to enhance energy storage, ultimately leading to excessive energy consumption to maintain homeostasis (Calabrese, 2001). This overcompensation response has been considered to be the main mechanism accounting for the low-dose stimulation of hormesis (Kim et al., 2018; Agathokleous et al., 2019; Costantini and Borremans, 2019). Homarine and taurine are known as organic osmolytes in mussels (Ji et al., 2014). The increased homarine in DEHP treatments from 4 to 324 μg/L verified that DEHP induced monophasic osmotic stress in mussels, which means that DEHP could affect cell integrity in mussels. The decreased taurine level in the group exposed to the highest concentration (324 μg/L) of DEHP might be due to taurine being used to compensate for the increased homarine. Furthermore, the increased inosine level implied that at 324 μg/L, DEHP disturbed muscle movement in mussels.