Abstract
The bacteria of shellfish bioaccumulate mainly in the gill and digestive gland tissues, which can affect shellfish health status and disease susceptibility under stressful conditions or environmental effects. Ruditapes philippinarum is a filter-feeding shellfish with both ecological and economic significance, and it is classified as a Mollusca (phylum), Bivalvia (class), Veneridae (family), and Ruditapes (genus). In this study, 16S high-throughput sequencing was used to explore the microbiomics of digestive glands and gills of R. philippinarum under temporary rearing without feeding (including the purification process) in different water environments (natural vs. artificial seawater) from 0–7 days. The results revealed that the digestive glands and gills of R. philippinarum had their own unique bacterial community structures. Tissue microorganisms under the overlap of different water environments and starvation factors showed different dynamic effects within 0–7 days. The sand spitting and purification steps (posttransport rehydration stage of 24 hours) effectively reduced microorganism abundance. There were different biomarkers in the prerearing and postrearing periods, and R. philippinarum may be more susceptible to the enrichment of opportunistic pathogenic bacteria during the postrearing period in artificial seawater.
1. Introduction
The Manila clam Ruditapes philippinarum is an economically valuable bivalve that is widely distributed in the intertidal zone of the Bohai and Yellow Seas in China [1]. Compared to other marine invertebrates, bivalves are more likely to accumulate and concentrate pathogenic microorganisms (including human pathogens and viruses) from seawater because of their filter feeding [2]. Host-microbiota relationships are of broad interest to researchers in marine invertebrates, and hosts can establish microbiota associated with healthy homeostasis through vertical and lateral transmission [3]. These host microorganisms significantly impact the host’s digestive function, metabolism, and infection immunity [4, 5], as well as assist the host in adapting to a changing environment and acting as a barrier against pathogens. Microbial composition structures are known to serve as health status markers in bivalve shellfish, and potential pathogens may lead to potential disease outbreaks during ecological dysregulation [3, 6].
Research on replacing natural seawater with artificial seawater has received attention in the fields of aquaculture and algae culture because of factors such as year-round variations in seawater quality and regional restrictions [7]. Artificial seawater primarily relies on artificial sea salt dissolution in inland aquaculture enterprises and commercial retail stores far from the coast. Artificial sea salt is a chemical product of seawater and saline, which is produced through a series of processes such as evaporation, centrifugation, and the concentration of seawater or brine, while others are prepared directly from mineral salts. To simulate the composition of seawater, trace elements, such as iron, manganese, zinc, and copper, are added during processing. Artificial seawater can maintain the osmotic balance and metabolism of seafood, keeping it alive in the short term. Rearing water, composed of dissolved artificial sea salt, has become the base water environment used in various studies of seafood such as Litopenaeus vannamei and Corbicula japonica [8, 9].
Intensive culture of bivalves is generally carried out on coastal mudflats [10]. After being harvested, the shellfish still has unspitted sand and various microorganisms accumulated in its body due to filter feeding. Commercial shellfish are often subjected to the sand spitting and purification steps process, where they are placed in a controlled aquatic environment that relies on the filtering of clean seawater to remove gastrointestinal contents thereby eliminating shellfish contaminants such as microbes and sands. Gills are susceptible sites exposed to the aquatic environment and serve as an important barrier against pathogen invasion. Digestive glands play a role in digestion, nutrient storage, and detoxification. Therefore, the abovementioned tissues were selected as the target organs for 16S high-throughput sequencing in this study, combined with alpha and beta diversity analyses, species taxonomic composition and relative abundance assessment, and linear discriminant analysis effect size analysis. The dynamic changes in the bacterial community structure of Manila clams’ digestive glands and gills in different water environments were investigated, and the enrichment patterns of opportunistic pathogenic bacteria were further compared and analyzed. This work provides fundamental data on controlling opportunistic pathogens in shellfish temporary rearing and guidance for the current potential of artificial seawater applications in the seafood industry.
2. Materials and Methods
2.1. Sample Collection
R. philippinarum were collected from Donggang City, China (39° 50′, 123° 45′) in December 2021 (Figure 1). The clams were tied in nylon mesh bags and transported to the laboratory by a cold chain for two days. Healthy, vigorous, and similar sized clams were selected and their biometric data were collected (shell length, 37 ± 2 mm; shell width, 23 ± 2 mm; and shell height, 16 ± 1 mm).

2.2. Seawater Preparation
Natural seawater was extracted from the Zhoushan Sea area at a depth of 20 m and transported to the laboratory after the purification process of sedimentation and sand filtration. Artificial seawater was created using a combination of tap water and sea salt. Before utilization, natural and artificial seawater was sterilized by UV. The salinity was 25 ± 1‰, the temperature was maintained at 23 ± 1°C, and dissolved oxygen was ≥7 mg/L. The contents of the main elements in natural and artificial seawater at the same salinity of 25% are shown in Table S1.
2.3. Temporary Rearing without Feeding
Temporary rearing was conducted using natural seawater and artificial seawater, and the processes were carried out in six tanks with a volume of 50 L. R. philippinarum were kept 20 cm from the bottom of the tank to prevent secondary absorption of spitted sand and reduce fecal and pseudofecal contamination. The seawater was replaced every 12 h with continuous aeration (without feeding). The sample names were set as follows: C: initial R. philippinarum group (sand-unspitted group); A: natural seawater group; B: artificial seawater group; H: digestive gland tissue; G: gill tissue; 0, 1, 3, 5, and 7: number of days of temporary rearing without feeding.
The digestive glands and gills from the original R. philippinarum samples were labeled as C0G and C0H, respectively. 21 clams (after selection in part 2.1) were randomly chosen and dissected with sterile blades and forceps. The two sample pools (digestive glands and gills) were separately and evenly divided into three parts before being placed into lyophilization tubes (i.e., C0H-1, C0H-2, C0H-3, C0G-1, C0G-2, and C0G-3). Then, they were snap-frozen in liquid nitrogen and stored at −80°C for subsequent microbiome analysis.
The remaining clams were divided into the natural seawater group (group A) and the artificial seawater group (group B). On the first day, 42 clams were randomly selected from 6 tanks (i.e., 7 from each tank), and the tissue collection and subsequent operations in groups A and B were repeated in group C. Tissues from group A were labeled as A1H and A1G, and those from group B were labeled as B1H and B1G. The same experimental operations and grouping arrangements were used for the next 3, 5, and 7 days.
2.4. Microbiological Analysis
DNA samples were extracted from the digestive gland and gill microorganisms according to the method described by Pan et al. [11], with some minor modifications. The V3-V4 region of the bacterial 16S rRNA gene was amplified by PCR (98°C for 30 s, 26–27 cycles at 98°C for 15 s, 50°C for 30 s, and 72°C for 30 s, and a final extension at 72°C for 5 min) using the primers 338 F (5′-ACTCCTACGGGAGGCAGCA-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′), and sample-specific 7-bp barcodes were incorporated into the primers for multiplex sequencing. The amplification results were subjected to 2% agarose gel electrophoresis, and the target fragments were excised and recovered using the Axygen Gel Recovery Kit. The PCR products were quantified on a microplate reader (BioTek, FLx800) using the Quant-iT PicoGreen dsDNA Assay Kit. Libraries were built using the Illumina TruSeq Nano DNA LT Library Prep Kit, and 2 × 250 bp double-end sequencing was performed on the qualified libraries using the MiSeq Reagent Kit V3 (600 cycles).
2.5. Sequencing Data Processing and Analysis
The microbiome biological information was analyzed according to QIIME2 (version 2019.4) and according to a process modified and improved by the official tutorial (https://docs.qiime2.org/2019.4/tutorials/). The raw sequence data were decoded and processed using the demux plug-in, first calling Qiime cutadapt trim-paired excised primer fragments of sequences. Unmatched primer sequences were discarded, and DADA2 was used for quality control, denoising, splicing, and chimera removal [12]. The obtained sequences were grouped at 100% similarity to generate tables of characteristic amplicon sequence variants (ASVs) and abundance data for subsequent analyses. The taxonomic information corresponding to each ASV was obtained using the Greengenes database, and the ASV characteristic sequences were compared with the reference sequences in the database. ASVs with abundance values lower than 0.001% (1 in 100,000) of the total number of sequenced samples were removed, and the abundance matrix of the removed rare ASVs was used for subsequent analysis.
The following data were processed using QIIME2 and R software. First, the total number of sequences in each sample in the ASV abundance matrix was randomly sampled at different depths to plot the sparse curves. Second, the ASV abundance matrix was randomly sampled with 95% of the sequence volume of the sample with the lowest sequence volume among all samples to correct for diversity differences between samples caused by sequencing depth. Diversity indices (Chao1, Shannon, and Good’s coverage indices) were calculated, and box plots were drawn. The nonmetric multidimensional scaling (NMDS) approach using Jaccard distance was visualized. Histograms were plotted for different classification levels based on ASV delineation and classification status identification results. LEfSe, which combines nonparametric Kruskal F02D Wallis and Wilcoxon rank-sum tests with the linear discriminant analysis (LDA) effect size, was used to detect categorical units that were rich in differences between groups.
3. Results
3.1. Sequencing Analysis
The samples yielded a total of 6078245 original sequence readings (H generated 2991654 and G generated 3086591 original sequence readings). Following quality filtering and denoising, an average of 101128 readings remained, and the proportion of effective sequences in all samples exceeded 84%, indicating that the valid sequences obtained from the experiment were able to fulfill the requirements of subsequent analysis. The overall sequence length of H and G was 2025543898 bp (Tables S2 and S3). The Venn diagram (Figures 2(a) and 2(b)) showed that H and G had 82 and 123 similar ASVs, respectively, indicating that the microbial composition of different tissue samples varied during purification. The rarefaction curves approached the saturation plateau, indicating that the sequencing depth was appropriate for all samples (Figure S1).

(a)

(b)

(c)

(d)
3.2. Alpha Diversity Analysis
Alpha diversity is also known as the within-habitat diversity. In this study, the Chao1 index was selected to characterize richness, the Shannon index to characterize diversity, and the Good’s coverage index to characterize coverage. The Good’s coverage values in this study ranged from 0.988% to 100%, indicating that the sample pool coverage was high, and the sequencing depth was sufficient to truly reflect the microbial community.
The various indicators of diversity within the AH, BH, AG, and BG are as follows (Figure S2): Chao1, , and 0.39; Shannon: , and 0.044; and Good’s Coverage: , and 0.44, respectively.
Figures 2(c) and 2(d) showed that the Chao1 indices of the control group (C0H and C0G) were higher than those of the temporary rearing group (AH, BH, AG, and BG). The Chao1 indices of AH and BH decreased continuously from 0 to 3 days and maintained a stable trend for 3–7 days. There was no significant difference within the diversity of H, whereas there were significant differences in the diversity of G. The Shannon in BG consistently showed an increasing trend from 1–7 days. In conclusion, temporary rearing of R. philippinarum was effective in reducing the abundance of the tissue microorganism community, especially at day 1. At 3–7 days, group A was more effective in controlling tissue microbial diversity than group B.
3.3. Beta Diversity Analysis
Beta diversity analysis demonstrated variability among samples by dimensionality reduction of multilatitude data. In this study, NMDS analysis was chosen to rank the sample distances by simplifying the data structure such that the sample ranking matched the distances of similar distances to each other as much as possible. NMDS results are generally considered more reliable when the stress value is less than 0.2. Figure 3 showed that intergroup differences were greater in the digestive glands than in the gills. The C0H, A1H, A3H, and B1H distances were more dispersed, and the remaining digestive gland sample groups were closer to each other. In the gills, the most significant differences were found between the C0G and BG samples (3–7 days). In conclusion, temporary rearing affected the digestive gland community more than the gills, and digestive gland microorganisms gradually became similar at 3–7 days. In different seawater samples, the respective microbial structures were observed in the gills.

(a)

(b)
3.4. Taxonomic Composition of Microorganisms
The species composition in the digestive glands and in the gills of R. philippinarum were characterized by phylum, class, family, and genus, respectively. As shown in Figure S3a, the dominant phyla in H were Proteobacteria, Tenericutes, Chlamydiae, Firmicutes, and Actinobacteria. Proteobacteria and Chlamydiae dominated C0H, and Tenericutes and Proteobacteria were the most abundant phyla during temporary rearing. In 3–7 days, BH had a higher percentage of Tenericutes than AH. The major phyla in group G were Proteobacteria, Firmicutes, Actinobacteria, Bacteroidetes, and Chlamydiae (Figure S3b). Proteobacteria maintained its highest abundance at 74.62–95.25%.
At the class level (Figure 4(a)), the percentage of Gammaproteobacteria increased significantly in A1H and B1H, with A1H and B1H group mean being 2.2 and 2.8 times that of C0H, respectively. Interestingly (Figure 4(b)), the relative content of Alphaproteobacteria in BG was higher than in AG, whereas the ratio of Gammaproteobacteria to Alphaproteobacteria was negatively correlated.

(a)

(b)

(c)

(d)
At the family and genus level (Figures 4(c), 4(d), S3c, and S3d), Pseudomonas (family Pseudomonadaceae) had the highest proportion in C0H and A1H, whereas the most dominant group in B1H was Pseudoalteromonas (family Pseudoalteromonadaceae). In the late rearing period (3–7 days), Mycoplasma (family Mycoplasmataceae) overtook Pseudomonas to become the most abundant component, and the relative proportions of BH were both greater than AH. The other subcomponents of the H group were Vibrio (family Vibrionaceae), Neptunomonas (family Oceanospirillaceae), Shewanella (family Shewanellaceae), and Rhodococcus (family Nocardiaceae). Rhodobacteraceae made up a comparatively larger portion of the BG than the AG in the G group, whereas Pseudomonadaceae showed the opposite trend. With increased rearing time, the proportion of Pseudomonas (family Pseudomonadaceae) in the BG decreased gradually, but there was no significant change in the AG.
3.5. Sample Variance Analysis
As shown in Figure 5, the LEfSe results (LDA = 3) revealed that there were 25 taxa available to distinguish gill bacterial communities, and there were relatively more differentials in the digestive gland. C0H was found to be abundant in 4 phyla, 13 orders, 11 families, and 5 genera. Verrucomicrobiaceae (phylum to family level), Clostridium (family to genus level), and Desulfobactorales (order level) were enriched in C0H, while Synechococcophycideae were enriched in both C0H and C0G. Pseudoalteromonas (order to genus level) were enriched in B1H, and orders Vibrionales and Shewanella (family to genus level) were enriched in A1G and B1G, respectively. In the later stages, Rhodobacterales (Sulfitobacter, Nautella, and Anaerospora), Rhodobacteraceae (Phaeobacter), Campylobacteraceae (Arcobacter), Aurantimonadaceae (Aurantimonas), and Oceanospirillaceae (Neptunomonas) had higher relative proportions in BH and BG than in AH and AG, which may be the main differences between artificial seawater and natural seawater.

(a)

(b)
4. Discussion
To date, few pieces of research studies have been recorded about the composition of the microbial community of the Manila clam (R. philippinarum). Meisterhans et al. used capillary electrophoresis DNA fingerprints to characterize the microbiota of organs (i.e. gut, gills, and remaining tissues) of Manila clams from different habitats [3]. Milan et al. employed 16s rRNA gene amplicon sequencing to characterize the hepatopancreas microbiota of the Manila clam, while analyzing the effects of seasonal fluctuations and chemical contamination on the clam microbiota [13]. In this study, 16S high-throughput sequencing was utilized for the first time to explore the tissue microbiomics (digestive glands and gills) of R. philippinarum in different water environments (natural vs. artificial seawater).
Statistical plots of alpha diversity, beta diversity, and taxonomic composition revealed that the digestive glands and gills had their own unique bacterial community structure. The diversity of the digestive gland was more stable than that of the gills. This was probably because the host had a strong filtering capacity for external microorganisms, and the digestive gland was potentially a more stringent environmental filter than the gill. The initial microbiota of shellfish came from the sea where the clams were grown and the fishery facilities where they were harvested and processed. The abundance of some microorganisms significantly enriched in the initial R. philippinarum tissues decreased after temporary rearing, probably because Clam partially changed its community structure by pumping clean seawater (up to 10 L/h) during environmental fluctuations, especially during the sand purification phase. Verrucomicrobia and Synechococcophycideae have been reported successively in various aquatic environments, marine plants, soil, and animal intestines [14–17]. The genus Synechococcus (phylum cyanobacteria) was considered a threat to the sustainability of marine ecosystems and can systematically accumulate in multiple tissues in bivalve mollusks (such as mussels, clams, and scallops), especially in the digestive glands [18, 19]. Synechococcus has been reported to produce toxigenic strains, which cause human diseases such as paralytic shellfish poisoning (PSP) through the traditional food chain pathway of filter-feeding mollusks (especially bivalves) [20].
The phylum Proteobacteria was consistently dominant in digestive glands and gills before and after purification. Proteobacteria were mainly identified as Alphaproteobacteria and Gammaproteobacteria, which have previously been shown to dominate in bivalves [21, 22]. The dominant genera in the digestive glands were Mycoplasma, Pseudomonas (family Pseudomonadaceae), Pseudoalteromonas, and Vibrio (family Vibrionaceae), and the most dominant genus in the gills was Pseudomonas. The class Mollicutes in digestive glands can be mainly identified as Mycoplasma, and Mycoplasma has the same dominance in Manila clams [13]. Mycoplasma spp. have parthenogenetic anaerobic properties and can ferment glucose or hydrolyze arginine, so they can adapt to the digestive tract and are therefore present at low levels (0–0.4%) in the gills [23]. The genus Pseudomonas (Class Pseudomonadaceae) was highly abundant in all tissues, especially in the gill tissues. Pseudomonas had strong environmental adaptability and was highly represented in skin microorganisms such as Fundulus grandis and rainbow trout [24, 25]. Pseudoalteromonas has been identified as the dominant genus in the fish gut [26]. Some species in this genus have competitive advantages in antibacterial, lysis, and antifouling abilities against other microorganisms and also contain species that enhance mortality in seafood [27]. The genus Vibrio is widespread in estuaries and seawater, and many species of this genus have been detected in bivalves, crustaceans, and fish [28].
From the temporal dimension, the microbial composition structures of AH and BH were significantly different from 0–1 days and 3–7 days, while the difference between BG and AG was mainly reflected in 3–7 days of BG. The dominant proportion of Gammaproteobacteria increased rapidly in the digestive gland at day 1 of decontamination. This may be because 0–1 days was the stressful phase of clams rehydrating after transport, during which the transport stress was lifted and viability began to recover. The relative content of Gammaproteobacteria in the BG continued to decline in 0–7 days, whereas the ratio in the AG remained more stable. Zhang et al. suggested that the Gammaproteobacteria spectrum contains the branched-chain amino acid (BCAA) catabolic pathway [29]. A reduced abundance of Gammaproteobacteria may lead to impaired neurotransmitter production.
The genus-level opportunistic pathogenic bacteria identified in this study were Psychromonas, Psychrobium, Shewanella, Vibrio, Nautella, and Arcobacter. Simultaneously, the relative abundance of tissues in the artificial seawater environment was generally greater than that in natural seawater (Figure 6). It is believed that the microbiota in stressed seafood is dominated by potentially pathogenic bacteria. It can be observed that artificial seawater has a greater impact on the stress of clams. The relative abundances of Psychromonas, Psychrobium, Shewanella, and Vibrio were high in the early stage and gradually decreased with the decontamination, especially during 1–3 days. The relative abundance of Vibrio showed a decreasing trend in the late period, whereas it remained the main potential pathogenic bacteria with Pseudomonas, Nautella, and Arcobacter in the late period. Nautella isolated from diseased samples of Litopenaeus vannamei were thought to be positively associated with Oryzias melastigma larvae, which showed a positive correlation with dysregulated brain genes [29, 30]. The genus Arcobacter had a relative abundance in B7H of 0.5%, which is 17 times higher than that of A7H. Some species of Arcobacter were associated with human gastroenteritis diseases [31]. The current species have been isolated mainly from seawater or seafood (especially shellfish), and their extremely high abundance has been repeatedly observed in unhealthy or dying marine organisms (e.g., oysters and abalone) [32]. However, their pathogenic mechanisms in shellfish remain unknown.

(a)

(b)
The physiological activity, pumping rate, and behavioral responses of the shellfish varied with changes in the seawater environment [33]. The results infer that 0–24 h was the peak period for sand removal and purification in R. philippinarum and 3–7 days was the stable period for temporary starvation rearing. In this study, R. philippinarum was subjected to stresses such as transportation, rehydration, and starvation from harvesting, and the specific mechanism of the effect of oxidative stress on clam microorganisms still requires further study. The elements in natural seawater differed from those in artificial seawater. At the same salinity, the macronutrients (Na, Mg, Ca, and K) were similar in artificial seawater than in natural seawater. Fe and Co were more abundant in the artificial seawater, whereas As was more abundant in the natural seawater. Compared to natural seawater, the elements in artificial seawater are less stable. In general, the various elements in artificial seawater are far less complete than those in natural seawater, and long-term culture may injure some of the metabolic functions of seafood and reduce enzyme activity and ion permeability.
The seawater was sand-filtered, and the possibility of phytoplankton or other large detrital materials supplementing energy requirements was excluded. The sample clams were continuously starved, and the substrate of energy metabolism in shellfish was altered, generally through endogenous energy (major nutrients such as fat, protein, and carbohydrates) or possibly through alternative energy sources to increase nutrition. The energy consumed primarily varies among the different species. As an alternative energy source, dissolved organic matter (DOM) in seawater can be absorbed by clams as a nutrient [34]. The ability of marine invertebrates to absorb DOM has been studied for more than a century, and DOM is considered to be an important source of energy for adaptation to periods of starvation and replenishment [35]. Different DOM levels simultaneously affect the composition of bacterial communities dispersed in an aqueous environment [36].
5. Conclusion
In conclusion, the bacterial communities of R. philippinarum were affected differently by different water and starvation conditions during the pretemporary and posttemporary rearing periods. Rehydration and sand splitting for 24–48 h can effectively reduce the tissue microbial diversity. The opportunistic bacterial community of artificial seawater compared to that of natural seawater may continue to increase with an extended duration of temporary rearing. Artificial seawater can still be an option for short-term cultivation of seafood if the resources of natural seawater are difficult to obtain or if there are contaminants such as heavy metals. The effect of artificial and natural seawater on the metabolism and flavor of R. philippinarum will be investigated in future studies.
Data Availability
The datasets used and analyzed during the current study are available from the corresponding author on reasonable request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Authors’ Contributions
Feixia Zhu conceptualized the study, designed the methodology, and wrote of the original draft. Xiangyu Zhang curated and validated the data. Jiao Liu and Xiaohan Ma conducted the investigation. Jing Yang collected the resources. Yue Wang designed the software. Yongjiang Lou acquired the funds. Yongyong Li supervised the study and reviewed and edited the manuscript. All authors read and approved the final manuscript.
Acknowledgments
This work was supported by National Key R&D Program of China (2020YFD0900903).
Supplementary Materials
Figure S1. The rarefaction curves of each sample based on Chao1 (a) and Shannon (b) diversity. Figure S2. Alpha diversity indices (AH(a), BH(b), AG(c), AG(d)) Figure S3. Relative abundance of the microbial community. The digestive gland organisms at the phylum (a) and genus (c) level; the gill organisms at the phylum (b) and genus (d) level. Table S1. Main elemental contents (mg/L) in seawater samples. Table S2. Statistic results of sequence data of digestive gland. Table S3. Statistic results of sequence data of gill. (Supplementary Materials)