TITLE PAGE

TITLE: Red fish, blue fish, native fish, new fish: eDNA as a tool to monitor fish in estuarine systems

RUNNING TITLE eDNA as a tool to monitor fish in estuarine systems

AUTHORS

Corresponding author: Alison W. Watts, University of New Hampshire, Alison.watts@unh.edu

Contributing authors:

Laura C. Crane, Wells National Estuarine Research Reserve, lcrane@wellsnerr.org

Jason Garwood, Apalachicola National Estuarine Research Reserve, Jason.garwood@FloridaDEP.gov

Jason S. Goldstein, Wells National Estuarine Research Reserve, jgoldstein@wellsnerr.org

Megan S. Lamb, Apalachicola National Estuarine Research Reserve, Megan.Lamb@FloridaDEP.gov

Christopher R. Peter, Great Bay National Estuarine Research Reserve, christopher.r.peter@wildlife.nh.gov

Yoshimi M. Rii, Hawaiʻi Institute of Marine Biology, University of Hawaiʻi at Mānoa, shimi@hawaii.edu

Shon Schooler, South Slough National Estuarine Research Reserve, shon.schooler@dsl.oregon.gov

Devin W. Thomas, University of New Hampshire, Devin.thomas@unh.edu

W. Kelley Thomas, University of New Hampshire, Kelly.thomas@unh.edu

Bree K. Yednock, South Slough National Estuarine Research Reserve, bree.yednock@dsl.oregon.gov

ABSTRACT

Environmental DNA (eDNA) is emerging as a potentially powerful tool to assess aquatic ecosystems, and to monitor fish assemblages. We conducted a pilot eDNA water sampling program at 5 Reserves within the National Estuarine Research Reserve System from temperate and tropical estuaries to better understand how eDNA data can support detection of fish species within an standardized monitoring network. Samples were collected in coordination with an existing long term monitoring program, and the sampling plan was designed to minimize additional work for field teams. We found that results varied; in some estuaries the number and type of fish species was consistent with expected occurrence, and eDNA analysis detected fish that often eluded traditional sampling methods. In warm turbid waters, however, we detected fewer species than expected, suggesting that additional processing is required at these sites. Managers interested in applying standardized eDNA monitoring across a range of sites may want to consider a decision tree protocol, where a baseline methodology is developed for all sites, enhanced by additional laboratory or analysis steps when specific conditions are present.

1. INTRODUCTION

Estuaries are among the most prolific habitats in the world, comprising complex and diverse ecosystems that are home to a wide range of species. They serve as both primary and nursery habitat for numerous fish including commercially important, rare and endangered, and invasive species (Kelley, 1988; Lenanton et al., 1987). Estuaries are also areas of high human activity that can influence fish communities and habitats, through shoreline development, global transportation, commercial and recreational fishing, and other industries (Kennish, 2002; Schulz et al., 2020). Management, conservation, and restoration in these regions rely on accurate assessments of the species present and how they vary over space and time (Fleishman et al., 2011). Although comprehensive species and habitat assessments with conventional monitoring methods provide a strong basis for understanding estuarine structure and function (Neckles et al., 2015), they are often time-consuming, expensive, and subject to technical and resource limitations. For example, traditional capture-and-release methods are labor-intensive and require extensive capacity, limiting the frequency of assessments and the ability to track temporal dynamics of fish communities. Many fish and invertebrates may also evade capture; species that are larger and faster may escape, while very small organisms may be able to move through nets, resulting in biased data. In addition, capture-and-release techniques are not universally applicable within all estuarine environments, due to bank and bottom topography, water depth, limitations on capture of sensitive species, and other environmental factors (Bulleri & Chapman, 2010; Kennish, 2002).

Because of these constraints, traditional monitoring methods may miss early detection of newly-arrived invasive species or losses of rare native species that occur in low densities or at locations that are difficult to survey. These factors highlight the importance of developing alternative methods and tools that may individually, or in combination with traditional methods, support more comprehensive monitoring in diverse environments. Underwater video, passive and active acoustic monitoring and environmental DNA (eDNA) are potentially powerful methods of augmenting traditional monitoring. Here we will explore the use of eDNA-based monitoring for fish species at five estuaries within the National Estuarine Research Reserve System (NERRS).

Advances in DNA methods and rapid reductions in analytical costs present an opportunity to harness this new technology and fundamentally improve our capacity to monitor biological communities and individual species (Baillie et al., 2019; Gilbey et al., 2021; Ip et al., 2021; Pawlowski et al., 2018; Thomsen & Willerslev, 2015). DNA present in a water sample is comprised of whole microorganisms (microalgae, bacteria, etc.), fragments of tissue, reproductive and waste products, and other cellular material. eDNA analysis allow resource managers to identify species in an aquatic ecosystem without having to actually capture and physically identify individual organisms. eDNA monitoring protocols are being developed for monitoring invasive species targeted for control, and for native species that may be in decline for special management attention (e.g. Crane et al., 2021; Darling & Mahon, 2011; Goldberg et al., 2016). eDNA methods have limitations, including the potential for false negative or positive results, vulnerability to contamination, misinterpretation of transported DNA and database limitations that may preclude identification of place-specific species (Jerde, 2019). In addition, eDNA results yield relative abundance of species within the sample, making data analysis and comparison beyond presence/absence difficult. Estuaries present particular challenges to eDNA studies; turbid water, transport from rivers, tidal forces, and shipping can complicate interpretation of results, and the use of eDNA to support coastal management has been limited to date. In this study, we evaluate the use of eDNA for fish surveys in estuarine settings as an additional tool for resource managers to support monitoring and assessment of coastal ecosystems.

1.1. National Estuarine Research Reserve System (NERRS) The NERRS is a network of estuarine research sites founded in 1972 under the United States Coastal Zone Management Act (https://coast.noaa.gov/czm/). At present, a total of 29 NERRs are designated and funded through a partnership between the National Oceanic and Atmospheric Administration (NOAA) and a state agency, university, or non-profit organization. Each NERR is subdivided into four distinct, but integrated sections: Stewardship, Research, Coastal Training, and Education sectors are each focused on the overarching goal of coastal management. The research and monitoring sector of the NERRS includes the System-Wide Monitoring Program (SWMP), which encompasses monitoring for water quality, biological species, and watershed and land-use/land-cover characterizations to examine short- and long-term change (Kennish, 2004; Sanger, 2002). The SWMP program is a required element of all NERRs and serves as a system-wide network of ecosystem-based monitoring with a large quantitative database produced through standardized operating procedures and data quality and data assurance metrics. The SWMP program thus results in initiatives to prioritize research in estuarine systems for understanding the structure and function of coastal wetlands for sustainability and future resilience. The SWMP initiative enables inter-and intra-NERR comparability of environmental characteristics of these systems; however, NERRs span multiple biogeographic regions representing multiple unique ecological habitats. Therefore, monitoring programs focused on biotic factors are more difficult to standardize. These constraints amplify the need to formulate methods to link biological research initiatives in a systematic manner among the NERRs. The addition of standardized methods and recommendations for eDNA sampling to the existing, long-term monitoring program at these sites would be crucial to examining changes in biological diversity and understanding the underlying mechanisms driving such changes across the Reserve system.

The use of eDNA has a tremendous range of research applicability from tracking rare or invasive species to assessing biodiversity at all trophic levels (Andruszkiewicz et al., 2017; Baird & Hajibabaei, 2012; Giakoumi et al., 2016; Gilbey et al., 2021). For management areas focusing on fish communities, eDNA is a promising monitoring tool to help identify species composition by augmenting or potentially replacing traditional surveys to track communities over time affected by local or global driving factors (e.g., eutrophication, dam removal, climate change), and further improve habitat management for fish production (Aglieri et al., 2020; Carim et al., 2017; Gilbey et al., 2021; Nguyen et al., 2020). However, as with any new method, eDNA must be field tested and evaluated to better understand where application is appropriate and useful, and where other methods might better serve management needs.

NERRS staff have identified management needs that generally center around identifying and mitigating impacts from climate change, eutrophication and other forms of contamination, advancing our understanding of ecosystem function, and investigating traditional and novel options for habitat restoration, all of which have the potential for advancement through eDNA monitoring. However, eDNA is a relatively new tool, with limitations associated with species identification, analytic and bioinformatic capacity, and interpretation of results. As part of the process of understanding the strengths and weakness of using eDNA to support estuarine management, we designed and implemented an eDNA sampling program at 5 NERR sites in Apalachicola FL, Great Bay NH, Heʻeia HI, South Slough OR, and Wells ME (Figure 1). The Hudson Bay NERR in NY also participated, but samples were collected in a tributary stream rather than an estuary, so the results are not included in this analysis. The overall goal was to deploy eDNA-based biological assessments as a complement to existing, long-term monitoring programs in a range of coastal settings to support estuarine management. There are many viable approaches to eDNA sampling, and project specific methods were selected based on project goals, end user needs, and the capacity of the science team, including field, lab, and analytic resources. We attempted to identify a practical, time efficient approach with a level of detail and specificity that complement existing NERRS SWMP and resources. The selected sampling methods are readily applied by a field team with access to a clean workspace and filtration apparatus.

2. MATERIALS and METHODS 2.1. Field sampling Samples were collected from January-November 2019 in 5 estuaries. Between 4 and 15 locations were sampled in each estuary with most of the samples collected at existing SWMP sites (Table 1). Replicate samples were collected from just below the water surface in new, sterile Whirl-pak bags or acid-washed or bleached sample bottles, then transported to a clean laboratory for filtering. Larger volume (5-L) samples were collected at South Slough in June and August 2019 and field-filtered using a portable peristaltic pump. Smaller volumes (500ml) were collected at Wells. A negative field control, consisting of lab or bottled water was opened at the site, then closed and processed with field samples. Collected samples were kept in a cooler during transport, and either frozen or filtered within 24 hours. Disposable nitrile gloves were worn throughout this process and all sampling equipment was either new/single use or sterilized with bleach or 10% w/v hydrochloric acid solution. Samples were filtered through 1.5 µm glass fiber filters (VWR 691) using a vacuum or peristaltic pump and stored at -20°C or -80°C. Frozen filters were packed in ice and delivered or shipped overnight to the lab. Prior to extraction, sample filters were cut in half, with one half retained for storage and one half processed for DNA extraction.

2.2. Extraction and PCR All samples were extracted at the Uviversity of New Hampshire Hubbard Center for Genome Studies in a dedicated lab space separate from PCR and raw sample handling. Extraction was performed following the protocols given in the QIAamp DNA Mini Hanbook (Qiagen, 2016); filters were placed in a lyse and spin basket with 400 mL of buffer ATL and 20 µL of proteinase K and incubated at 56 °C for 1 hr, then centrifuged. The remainder of the filter extraction was performed on a QIAcube Connect system (QIAGEN®, Hilden, Germany) following the QIAamp DNA Mini Kit program. Tissue collected for voucher specimens was extracted with DNeasy Blood and Tissue kits. One negative control sample was extracted with every 11 field samples. DNA concentration was measured using an Invitrogen Qubit 2.0 (Thermo Fisher®, Waltham, MA). The manufacturer’s protocol for the Qubit dsDNA HS Assay kit was followed using 1 uL of sample. Samples were kept at -20°C until processed, then transferred to -80°C storage.

A 12S rRNA metabarcoding assay designed for fish (Miya et al., 2015) was applied in a two-step process. After the first round of thermocycling, the presence of target DNA was confirmed with gel electrophoresis, then TrueSeq adapters were added, and a second round of PCR was conducted. A negative and positive control was included with each PCR set. Most of the samples were processed using the following thermocycling profile: Step 1) 95 °C for 3 min, 2) 98 °C for 30 sec, 3) 55 °C for 30 sec, 4) 72 °C for 30 sec, 5) repeat steps 2 to 4 43 times, 6) 72 °C for 5 min, and finally 7) hold at 4°C. Samples collected during the summer in South Slough contained bacteria sequences that masked the 12S response. These samples were processed using a Touchdown procedure (Pitz et al., 2020) for step 1, then processed in the same manner as other samples for step 2. Samples were then sequenced on an Illumina HiSeq 2500 (Illumina Inc., San Diego, CA), demultiplexed, and returned for bioinformatic analysis.

All samples were diluted 1:10 to reduce inhibition, and amplification in each sample was confirmed by reviewing gels after the first PCR application. The extraction method used (Qiagen QIAmp DNA) includes inhibition removal, and was generally very effective. In rare cases where inhibition was apparent, samples were re-extracted (from the reserved ½ filter) using PowerSoil Pro (Qiagen) which has enhanced inhibition removal. Bands of competing sequences could also be observed on gels, particularly 16S bacterial sequences, which in some samples completely overwhelmed the fainter 12S bands. In these samples, gel bands were excised to isolate target sequences. These efforts improved the number of reads associated with each detection, but did not increase the number of fish species detected. We also applied a cytochrome c oxidase subunit I (COI) primer (Leray et al., 2013) to a subset of samples, without substantial increase in number of species detected.

2.3. Bioinformatics and data analysis The demultiplexed sequence data was processed through QIIME2 dada2 (Caporaso et al., 2018) with a truncation length of 150 bp for both the forward and reverse reads, and with a trim length of either 21 bp forward, 27 bp reverse or 27 bp forward and 33 bp reverse depending on whether the MiFish primers used included a 6 bp spacer. Each lane was denoised separately, and the results were merged together. An initial approximate taxonomy was done comparing the sequences to a reference set of 12S sequences downloaded from NCBI https://github.com/dwthomas/MiFish-Reference-Database. Unassigned and off target sequences are filtered from the samples, and the finalized taxonomy was curated by hand using the web interface to BLAST to manually assign taxonomy, based on similar reference sequences, and knowledge of the geographic range of place-specific species. This draft species ass In a few cases, improbable species identified in a sample were traced to lab contamination. In these cases the suspect species were flagged and removed from the sample list. Non-fish chordates, such as mammals and birds, were segregated from the fish results.

Species accumulation curves were generated for each site; incidence data was uploaded to iNEXT online, which computes seamless interpolation and extrapolation using Hill numbers.  An order hyperparameter of q = 0 was selected, corresponding to the species richness.  The hill number with incidence data counts the effective number of equally frequent species in the assemblage (Chao et al., 2014). Species accumulation curves were calculated using iNEXT online and plotted with matplotlib (Hsieh et al., 2016).

2.4. Traditional fish surveys The approximate number of fish species expected to be present was provided by resource managers at each site based on long term surveys and local knowledge. Long term survey methods vary, but included seine surveys (Great Bay, Apalachicola, South Slough), larval fish tows (Wells), cast net and visual surveys (He’eia). Managers augmented these data with information on species likely present, but not always captured by common survey methods. For example American Eel (Anguilla rostrata ) is present in the Great Bay region, but rarely caught in estuarine seine surveys. They used this information to estimate the number of fish species likely present in the estuaries at each site to allow comparison with total number of species detected with eDNA. In Great Bay a more detailed species comparison was conducted by pairing eDNA sampling with seine surveys performed by New Hampshire Fish and Game during the same time frame; seining was conducted using a bag seine (30.5 m long and 1.8 m high, 6.4 mm mesh). Seine hauls were performed during daylight within 2 hours of low tide and were set by boat 15–25 m from shoreline. In September 2019, water samples were collected for DNA analysis in coordination with seine sampling at four locations. Water samples were collected a day after the seining to avoid cross contamination and species disturbance associated from the net and boat.

3. RESULTS Overall, a total of 124 unique fish species were detected in over 370 environmental samples (Figure 2). Most of the species detected represent fish known to be in the region, with some exceptions discussed below. Taxa from the same family which are clearly different, but cannot be resolved to the species level are indicated by ‘sp.’ Many fish were only detected in one sample, while the more abundant species are generally detected more frequently.

In cases where sequence data identifiedd a species that was not likely to be present in the system, we re-reviewed all lab and bioinformatics QA steps to ensure that neither contamination nor incorrect sequence assignment were likely. We then considered the local estuarine environment to identify potential DNA transport pathways. In several cases, known commercial or recreational fisheries were noted as likely sources: In South Slough both hagfish (Eptatretus sp. ) and albacore (Thunnus alalunga ) were identified in eDNA samples during this study. Although not likely to be living in the region, both species are handled at local commercial facilities. In Apalachicola, an estuarine grouper (Epinephelus sp .) was detected in a freshwater tributary, probably transported by recreational fishing staging from a nearby dock. At all sites, DNA from freshwater species known to be present in the watershed was detected occasionally in estuarine samples. Detections of freshwater fish were most common near the entrance to tributary rivers, indicating transport from upstream reaches. In contrast, at Heʻeia, DNA from 3 different native goby species was detected at the specific watershed sites where they have been detected through other methods, indicating that DNA transport may be limited in that system. This is consistent with studies finding that degradation and transport of DNA fragments in aquatic environments will vary depending on characteristics such as tidal flow, water chemistry and temperature (Caza-Allard et al., 2021). Additional species which do not have a clear source were occasionally detected; some of these are fish that are known to be present in the region, but further offshore such as Atlantic flyingfish (Cheilopogon melanurus ).

The number of fish species detected in each sample is shown in Figure 3. The number of species ranged from 0 to 7 in individual samples. In the cases where few or no fish species were detected the samples were reviewed to ensure that inhibition or other issues were not interfering with amplification.

At most locations triplicate samples were collected. Replicate samples provide additional information, without adding additional field travel time. Figure 4 shows the number of fish species detected in triplicate vs single 1-liter samples in Great Bay in June 2019. Triplicate samples approximately doubled the number of detected species, although costs associated with analyses also increase.

4. DISCUSSION

The estuaries sampled in this study represent very different site conditions (e.g., latitude, geomorphology, temperature, salinity, etc.) and fish communities. By applying the same sampling and analysis method at each location we are able to compare the results from a standardized method to better understand the practical use of eDNA monitoring in estuarine systems. We considered two questions commonly asked by resource managers: How many species are detected by eDNA sampling in a given estuary? How does eDNA-based monitoring compare to traditional fish surveys? And what is the relative cost of these methods?