Molecular Diet Analysis of Asian Clam for Estuarine Biodiversity Monitoring: A Case Study of Nakdong River Estuary

1. Introduction

Estuaries are highly complex and dynamic ecosystems that provide diverse habitats, such as tidal mudflats, sandbars, marshes, and transition zones for numerous organisms [1,2,3]. However, estuaries are often subjected to serious threats from anthropogenic impacts, including overexploitation, reclamation, pollution, and barrage construction, leading to rapid declines in habitat and biodiversity [4,5]. Thus, monitoring and detecting significant changes in estuarine ecosystems and biodiversity are important [6]. Traditional aquatic species monitoring methods that rely on direct detection or capture are time-consuming and labour-intensive [7,8]. Morphological identification using the naked eye or a microscope is skill-dependent and can result in misidentification [9]. Furthermore, the dynamic environment and high biodiversity of estuarine ecosystems complicate monitoring [10].

Environmental DNA (eDNA) metabarcoding is a promising alternative to traditional monitoring methods [11]. This molecular technique enables the identification of the entire community from a single environmental sample (for example, water, soil, air, faeces, or gut contents) without directly observing or capturing the organisms. Through polymerase chain reaction (PCR) using universal or group specific primer and next-generation sequencing (NGS) techniques, researchers can identify the presence of various organisms, including rare, elusive, or endangered species [11,12]. This method provides a noninvasive and efficient way to assess biodiversity in different habitats [13]. Consequently, studies employing eDNA metabarcoding for biodiversity assessment have gained significant attention in recent years and have explored various potential eDNA sources, such as biofilms, faeces, or gut contents [14,15,16].

Filter feeders, which filter water and ingest organic particles, may accumulate eDNA within their bodies without an artificial filtering process and can be used to identify biodiversity in aquatic ecosystems [15,17,18]. Their unique and effective feeding habit makes them potential eDNA samplers [15,19]. In particular, bivalves such as clams and mussels are widely distributed across aquatic ecosystems, including lakes, rivers, and estuaries, demonstrating broader applicability than other filter feeders [20]. However, few studies have extracted eDNA from bivalves, and its potential has rarely been investigated, especially in specific ecosystems, such as closed estuaries.

In this study, we conducted a first-of-its-kind investigation applying bivalve eDNA metabarcoding for diversity monitoring in a closed estuary, using the Nakdong River Estuary as a case study. This estuary exemplifies an artificially regulated ecosystem due to the presence of an estuarine barrage, which was constructed in 1987 and reopened in 2022 to restore a brackish ecosystem [21]. The bivalve C. fluminea is widely distributed in the Nakdong River Estuary. We hypothesised that eDNA analysis of C. fluminea gut contents could reveal the biodiversity of the Nakdong River Estuary.

To test this hypothesis, (1) we extracted eDNA from the gut contents of C. fluminea collected from an estuary and investigated its biodiversity using metabarcoding analysis. Then (2) we compared the results of the three study sites in the Nakdong River Estuary. Finally, the potentials and limitations of this case study are discussed.

2. Materials and Methods

2.1. Study sites

The Nakdong River is the second-largest river system in South Korea and maintains a well-developed estuarine system (35° 05' N, 128° 55' E). It is also recognised as an important biodiversity conservation area, including winter bird habitats and stopover sites on the East Asia-Australasian Flyway [22]. Another important feature of the estuary is flood control activities, mainly through the estuarine barrage built in 1987. This structure divides the brackish areas into distinct freshwater and saline zones. This division is believed to influence biodiversity changes [21], necessitating more detailed and efficient monitoring methods to address growing concerns regarding regional biodiversity protection. In the present study, we selected three sites within a brackish area, each with different salinity levels (Figure 1). Site 1 is closest from the barrage), with Site 3 situated between the two (approximately 2.0 and 3.9 km from the barrage, respectively) and Site 2 situated between the two (approximately 2.7 km from the barrage).

2.2. Water quality survey

A field survey was conducted in September 2021 concurrently with the sampling of C. fluminea. Water samples were collected from the surface layer (at a depth of approximately 0.5 m) using a 10 L polypropylene bucket. Dissolved oxygen (D.O., mg L^-1, %) and water temperature (°C) were measured with a YSI 550A Dissolved Oxygen Instrument (YSI, Ohio, USA). The pH levels were determined using a YSI Model 60 handheld pH–temperature system (YSI, Ohio, USA). Electrical conductivity (μS cm^-1) and salinity (ppt) were assessed using a YSI Pro30 Conductivity Meter (YSI, Ohio, USA). Following the field survey, the water samples were transported to the laboratory in refrigerated storage for turbidity and alkalinity analyses. Turbidity (NTU) was measured using an APERA TN500 Portable White Light Turbidity Meter (APERA, Ohio, USA), and alkalinity (mg L^-1) was determined using the neutralisation method in accordance with standard procedures [23].

2.3. C. fluminea sampling and pretreatment

C. fluminea samples were collected using a fishing dredge with a width of 123 cm and height of 22 cm (Figure A1). The dredge net was made of polyethylene and was 320 cm long with an 11 × 11 mm mesh size. Collected C. fluminea samples were placed in separate polyethylene bags at each site and transported to the laboratory for refrigerated storage. The samples were stored at -80 °C until further analysis.

From the collected C. fluminea samples, ten individuals per site were randomly selected (totalling 30 individuals; shell length (cm), 2.0–3.0; shell height (cm), 1.9–2.4; shell width (cm), 1.1–1.5; and total weight (g), 2.990–6.051). The gut was eviscerated and dissected to obtain the contents. During the dissection process, scalpels, tweezers, and scissors were flame-sterilised between samples to minimise contamination. The extracted gut contents were placed in 1.5 mL microtubes separately (n=30) and stored at -20 °C until further analysis.

2.4. DNA extraction and amplification

The gut contents from the C. fluminea samples were homogenised with sterilised homogeniser pestle, and gDNA was extracted using ‘DNeasy Blood & Tissue Kit’ (Qiagen, Hilden, Germany), according to manufacturer’s instructions. The extracted gDNA samples were stored at -20 °C.

We implemented two consecutive PCR steps for next-generation sequencing (NGS) process. The 1st PCR was performed using primer sets targeting universal eukaryotes to amplify the V9 regions of 18S rRNA (18S V9 primer). The forward primer sequence is 5’-TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGCCCTGCCHTTTGTACACAC-3’, and the reverse is 5‘-GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAGCCTTCYGCAGGTTCACCTAC-3’. We used ‘AccuPower HotStart PCR PreMix’ (Bioneer, Deajeon, Korea), and the volume of the PCR reaction solution was total 20 uL (DNA template 1 uL, forward primer 1 uL, reverse primer 1 uL, and distilled water [D.W.] 17 uL). The PCR condition consist of 1 cycle of initial denaturation (94 °C, 10 min) and 35 cycles of Denaturation (94 °C, 1 min), Annealing (50 °C, 1.5 min), extension (72 °C, 1 min) and 1 cycle of final extension (72 °C, 10 min). After 1st PCR, we confirmed the size of the products using 1.5% agarose gel electrophoresis and stored them at -20 °C.

The 2nd PCR was performed using ‘KAPA HiFi HotStart ReadyMix’ (KAPA Biosystems, Wilmington, MA, USA) and ‘Nextera XT Index Kit v2’ (Illumina, San Diego, CA, USA). The volume of the PCR reaction solution was total 25 uL (1st PCR product 2.5 uL, Forward index 2.5 uL, Reverse index 2.5 uL, KAPA mix 12.5 uL, D.W. 5 uL). The PCR condition consist of 1 cycle of initial denaturation (95 °C, 3 min) and ten cycles of denaturation (95 °C, 30 s), annealing (55 °C, 30 s), extension (72 °C, 30 s), and 1 cycle of final extension (72 °C, 5 min). Subsequently, we confirmed the size of the products by 1.5% agarose gel electrophoresis and stored them at -20 °C.

The PCR products were purified by beads clean-up process using ‘AMPure XP Reagent’ (Beckman Coulter, Indianapolis, IN, USA) and then pooled in equal concentration (10nM) using ‘DeNovix QFX Fluorometer’ and ‘Denovix dsDNA Ultra High Sensitivity Assay’ (Denovix, Wilmington, DE, USA) according to the manufacturer's protocol. The generated library was stored at -20 °C until DNA sequencing.

2.5. DNA sequence analysis and statistics

We sequenced library samples using NGS and performed taxonomic identification. The library was sequenced on an Illumina iSeq platform (Illumina, San Diego, CA, USA) and data processing was performed using USEARCH (v11.0.667) [24]. Demultiplexed raw sequences (FASTQ files) were merged into one sequence, allowing a maximum of ten mismatches. Merged reads with expected errors>3.0 were discarded after quality filtering. The remaining sequences were dereplicated and clustered into operational taxonomic units (OTUs) at a 97% OTU cutoff value, removing chimeric and singleton sequences.

The resulting OTUs sequences were searched against the NCBI database (Release 250.0; July 2022) using BLASTn [25]. We obtained a list of the top 100 taxa for each OTU with the highest identity percentage. Considering the known characteristics and distribution of each taxon, we excluded taxa that were not present at the study sites. If all taxa found in the search results had low scores (query cover<80%, identity <90%) or were regarded as not inhabited, the OTUs were deleted. Finally, the one with the highest identity was assigned to each OTU. OTUs with 97% or higher identity were identified at the species level, and the rest (93–97%) were identified at the gene level. OTUs identified as C. fluminea were considered ‘self-DNA’ and excluded from the subsequent process. Obtained sequences were deposited in the NCBI repository under accession number SAMN35796656.

Next, we categorised the OTUs identified through taxonomic classification and general habitat based on the 2022 World Register of Marine Species (WoRMS) database (https://www.marinespecies.org/) and additional references (Table A1). The number of OTUs was visualised as pie charts and bar graphs using Microsoft Excel 2019 (Microsoft, Redmond, WA, USA) based on the phylum and habitat categories. The taxon composition at each study site was visualised using a Venn diagram. The number of OTUs and taxa composition at each study site were analysed based on the phylum and habitat categories. Rarefaction analysis was performed using Past 3.01 (Natural History Museum – University of Oslo, USAOSLO Norway) to determine the sample size required to obtain a sufficient number of OTUs.

3. Results

3.1. Water parameters

The three study sites at the Nakdong River Estuary exhibited salinity values ranging 10.0–14.2 ppt, indicative of a brackish area (Table 1). Salinity increased in the following order: Site 2, Site 3, and Site 1, regardless of the distance from the estuarine barrage. Electrical conductivity showed a trend similar to that of salinity, as it was also affected by dissolved salts. The other parameters (that is, dissolved oxygen, pH, water temperature, alkalinity, and turbidity) were not significantly different among the study sites.

3.2. eDNA from the gut contents of C. fluminea

DNA metabarcoding analysis generated 17,272 paired-end reads from 30 samples. After quality filtering, 16,980 (98.3 %) sequences were obtained, comprising 21 OTUs. OTUs were identified as 21 eukaryotic taxa (ten phyla, 14 classes, 17 orders, 19 families, and 20 genera) belonging to Animalia, Chromista, Fungi, and Plantae. We identified 12 OTUs at the species level (that is, identity ≥ 97%) and nine OTUs at the genus level (94.41–97%) (Table 2). Cyclops sp., Neocercomonas sp., and Nannochloris sp. were identified at the genus level despite having an identity higher than 97%, as the search results from the database were already at the genus level.

The OTUs identified in the gut contents of C. fluminea recovered taxonomically diverse species and genera. Among the results, Animalia accounted for the highest percentage at 38.1% (eight OTUs), followed by Chromista (28.6%; six OTUs), Plantae (23.8%; five OTUs), and Fungi (9.5%; two OTUs) (Figure 2). Animalia comprises the phyla Arthropoda, Chordata, and Platyhelminthes and includes Mytilicola orientalis, a species known to be parasitic in the intestines of bivalves. Fish accounted for 9.52% (two OTUs) of the total, and copepods accounted for 14.29% (three OTUs). The OTUs in the kingdom Chromista spanned four phyla, with Ochrophyta representing the highest proportion. Fungi contained the lowest number of OTUs detected. All Plantae OTUs belonged to the phylum Chlorophyta (green algae).

The OTUs recovered from the various habitat environments were marine, brackish, freshwater, and terrestrial. The freshwater habitat group had the highest percentage of 52.4% (11 OTUs), followed by the terrestrial (19.0%; four OTUs), marine (19.0%; four OTUs), and brackish (9.5%; two OTUs) groups (Figure 2). The freshwater group comprised six Animalia OTUs and four Plantae OTUs (green algae). Marine and brackish groups primarily included Chromista species and genera. Terrestrial species and genera, including Liposcelis sp. (booklice), Aspergillus penicillioides, Trichosporon asahii (fungi), and Neocercomonas sp. (Cercozoa), were also detected.

3.3. Comparison between study sites

We compared the number and composition of OTUs among the three study sites representing different salinity levels. Of the 21 OTUs, 16 OTUs were detected at Site 1, nine OTUs at Site 2, and 14 OTUs at Site 3 (Figure 3), indicating that the number of OTUs at Sites 1 and 3 was 25% higher than that at Site 2.

We compared the number of OTUs based on taxonomic classification (Phylum). Sites 1 and 3 recovered eight and nine of the ten phyla, respectively (Figure 3). At Site 2, where the number of OTUs was the lowest, only five phyla were recovered, and the kingdom Chromista was undetectable. Although Site 3 had fewer OTUs than Site 1, it had the largest number of phyla among the sites.

In addition, we compared the number of OTUs based on the general habitat. Sites 1 and 3 included species and genera inhabiting terrestrial, freshwater, brackish, and marine environments, whereas Site 2 did not include a brackish group (Figure 3). At all sites, the freshwater group had the highest percentage (56.3%, 77.8%, and 64.3% at the three sites, respectively) and the percentage of brackish water was the lowest (6.3%, 0.0%, and 7.1%, respectively).

We also compared the OTU compositions of the three sites. Of the 21 detected OTUs, we identified 11 “common OTUs” (OTUs detected in more than one site; 52.4% of the total) and ten “site-specific OTUs” (OTUs detected in only one site; 47.6% of the total) (Figure 4). Sites 1 and 3 recovered five site-specific OTUs, whereas Site 2 did not contain any site-specific OTUs. Thus, all OTUs detected at Site 2 were common OTUs.

3.4. Relationship between sample size and number of OTUs

A rarefaction analysis was conducted to determine the sample size required to detect sufficient OTUs for biodiversity monitoring. After analysing 30 samples, the expected number of OTUs did not converge, although the sample size reached its maximum value (Figure 5). However, the curve exhibited a gradually decreasing slope, suggesting an approximate trend in the number of OTUs.

4. Discussion

In this study, we detected 21 OTUs in the gut contents of 30 C. fluminea samples collected from the Nakdong River Estuary. These results suggest that species and genera are present at the study sites and support the applicability of C. fluminea as an eDNA sampler for estuarine diversity monitoring. The gut contents of C. fluminea contained taxonomically diverse organisms, including the phyla Arthropoda, Chordata, Ochrophyta, and Chlorophyta, and successfully recovered diverse habitats in the estuaries. Moreover, this is the first study to apply bivalve eDNA to a closed estuary, providing a perspective on eDNA metabarcoding in regulated estuarine ecosystems.

To examine biodiversity through faeces or gut contents, the characteristics of both the habitat and the sample species should be carefully considered. Our sample species, C. fluminea, has an efficient filter-feeding mechanism and has potential as a biological indicator of microplastic or heavy metal pollution and for environmental applications, such as eutrophication alleviation [26,27]. They do not represent strong selectivity for food [28], enabling the detection of diverse species, and are widely distributed in aquatic ecosystems worldwide [29], indicating broad application. Importantly, as C. fluminea can accumulate eDNA in its body without an artificial filtration process, they may serve a complementary role for filtered water samples, particularly in environments where it is difficult to use water samples due to high sediment loads and elevated turbidity, such as estuaries [30].

However, C. fluminea is a typical benthic filter feeder [31], and the possibility of the dominant detection of benthic organisms cannot be ruled out. In addition, as freshwater and seawater are mixed, seawater with a high-density flows forms a vertical salinity gradient in the brackish zone [32,33]. Salinity is an important factor that determines the distribution of aquatic organisms [34,35]; therefore, salinity stratification can influence the biological community, suggesting that species that prefer saltwater may be detected dominantly in eDNA analysis. However, in the present study, the OTUs did not tend to be biased toward benthic or marine organisms. The exact reason for this is unknown, but the biological distribution by depth or salinity does not seem to have a significant impact on eDNA detection using C. fluminea. These results suggest a positive effect of the availability of C. fluminea eDNA in estuarine ecosystems. However, a more detailed investigation of the relationship among water depth, salinity gradient, and eDNA detection is required.

Our study site, the Nakdong River Estuary, represented the ecosystem of a typical closed estuary [21]. This artificial regulatory system may exhibit different singularities from those of natural estuaries (open estuaries) in biodiversity monitoring using eDNA metabarcoding. The major distinction is that the estuary barrage separates freshwater from seawater and blocks the water flow. These features make it more difficult to filter water for eDNA analysis because of the accumulation of more floating organic matter and sediments than in natural environments [30]. However, eDNA transport and dispersal from the origin may also be prevented, thereby enabling more efficient eDNA sampling. Our study has yet to address this issue; therefore, it is impossible to suggest an effect on eDNA metabarcoding in C. fluminea. However, it will be possible to investigate and compare these with natural estuaries in future studies.

In the present study, 18S V9 primers were used to target universal eukaryotes. The 18S rRNA V9 barcode is suited to DNA metabarcoding for diet analysis because of potential to amplify degraded DNA and detect a relatively broad range of eukaryotic organisms [36,37]. However, because the targeting region is relatively short, the 18S rRNA region is considered limited for distinguishing taxonomically close species [38,39]. To compensate for this, we investigated the known characteristics and geographical distribution of each taxon identified using BLASTn, and directly excluded taxa that did not inhabit South Korea. This process minimises the risk of false positives that can be caused by short amplicons, compared to a simple program that selects the taxa with the most similar sequences from the BLASTn results. For efficiency, it is possible to use AI or program codes to make this process more convenient.

We also detected various species and genera inhabiting estuarine habitats, but obtained only 21 OTUs, which was less than expected. This is probably a limitation of eDNA sampling through direct gut content extraction. According to Heo, et al. [40], the gut content of Asian clams (C. fluminea) recovered relatively fewer OTUs than the pseudo-faeces. Jeunen, et al. [18] reported that eDNA could not be obtained from the gills of blue mussels (Mytilus galloprovincialis). Since there are few cases of eDNA extraction from bivalves and no established protocols exist, various sample processing methods are needed to examine their potential.

DNA metabarcoding using faeces or gut contents usually poses the risk of excessive detection of self-DNA (that is, DNA of the sample species itself) or the DNA of the species interacting with the sample species (for example, parasites and symbionts) [41,42,43]. It can be excessively amplified during DNA processing and detected in much larger amounts than pre-DNA, leading to non-informative results [42,44]. Our study also identified ten OTUs as C. fluminea and the OTU of Mytilicola orientalis, which is known to be parasitic in the intestines of bivalves. Therefore, to avoid these biases, it is necessary to consider alternatives, such as blocking oligonucleotides that prevent self-DNA amplification, in further studies.

5. Conclusions

In this study, we investigated the potential use of C. fluminea as an eDNA sampler to assess the biodiversity of the Nakdong River Estuary. We extracted eDNA from the gut contents of C. fluminea and used the 18S V9 primer to examine overall biodiversity at the study sites. The metabarcoding results recovered 21 taxa belonging to different kingdoms, including Animalia, Chromista, Fungi, and Plantae, and represented a wide range of environments, including marine, brackish, freshwater, and terrestrial environments. These findings indicate the potential use of C. fluminea in eDNA metabarcoding for biodiversity investigations. This is the first case study to apply bivalve eDNA metabarcoding in a closed estuarine system, supporting the potential application of bivalves as eDNA samplers in diverse aquatic ecosystems.