Accurate knowledge on spatiotemporal distributions of marine species and their association with surrounding habitats is crucial to inform adaptive management actions responding to coastal degradation across the globe. Here, we investigate the potential use of environmental DNA (eDNA) to detect species-habitat associations in a patchy coastal area of the Baltic Sea. We directly compare species-specific qPCR analysis of eDNA with baited remote underwater video systems (BRUVS), two non-invasive methods widely used to monitor marine habitats. Four focal species (cod Gadus morhua, flounder Platichthys flesus, plaice Pleuronectes platessa and goldsinny wrasse Ctenolabrus rupestris) were selected based on contrasting habitat associations (reef- vs. sand-associated species), as well as differential levels of mobility and residency, to investigate whether these factors affected the detection of species-habitat associations from eDNA. To this end, a species-specific qPCR assay for goldsinny wrasse is developed and made available herein. In addition, potential correlations between eDNA signals and abundance counts (MaxN) from videos were assessed. Results from Bayesian multi-level models revealed strong evidence for a sand association for sedentary flounder (98% posterior probability) and a reef association for highly resident wrasse (99% posterior probability) using eDNA, in agreement with BRUVS. However, contrary to BRUVS, eDNA sampling did not detect habitat associations for cod or plaice. We found a positive correlation between eDNA detection and MaxN for wrasse (posterior probability 95%), but not for the remaining species and explanatory power of all relationships was generally limited. Our results indicate that eDNA sampling can detect species-habitat associations on a fine spatial scale, yet this ability likely depends on the mobility and residency of the target organism, with associations for sedentary or resident species most likely to be detected. Combined sampling with conventional non-invasive methods is advised to improve detection of habitat associations for mobile and transient species, or for species with low eDNA concentrations. 

Water samples were collected directly from the surface water and extracted using a disposable 60 mL sterile syringe and injected into an enclosed Sterivex-GP capsule filter (0.22 µm pore size, SVGPL10RC, Millipore, CA, USA). A total of 1000 mL seawater was filtered for each sample, however on one occasion, due to clogging, only 850 mL could be filtered. After filtration, the sterivex filters were closed using sterile luer lock caps (Cole Palmer, Vernon Hills, IL, USA), put in individual zip-locked bags and stored in a cooling box on ice before subsequently transferring them to a -20˚C freezer upon return from the field. A maximum of two water samples were taken on a given sampling day, for which site selection was based on accessibility of the site (depending on weather conditions) or otherwise balanced among sites as much as possible. One field blank was taken on 18 June 2018 at 12:50h by filtering 900 mL nuclease free water at site to control for potential DNA contamination related to the field sampling procedure. Further details on the sampling scheme are provided in the Supporting Information (Table S1).

Immediately following the eDNA water sampling, we deployed two BRUVS at the same field site. BRUVS consisted of GoPro’s Hero 3 or 4 (GoPro, USA), attached to a steel pole (1 m height; 3 cm diameter) at 20 cm above the seabed and positioned firmly on a concrete tile with the camera field of view parallel to the seabed. BRUVS deployments were separated by at least 100 m to minimize the risk of double counting individuals between the two BRUVS while staying within the designated habitat type. We used chopped Atlantic herring (Clupea harengus) as bait, packed tightly in a mesh bait bag. The bait bag was positioned 15 cm above the seabed and attached via an 80 cm bait arm in the lower center of the field of view. Soak time of the BRUVS varied between 1h – 2h 40min depending on the camera’s battery life, with a small number of BRUVS (< 5) marked unsuccessful due to battery failure or loss of the bait bag. Accordingly, a total of 41 seawater samples were taken for eDNA extraction and coupled with 78 BRUVS deployments across the four different field sites (Fig. 1; Table S1).  

We analyzed the video recordings using VLC media player (https://www.videolan.org/vlc/). The four focal species were identified and counted by trained video observers making use of reference images and conspicuous morphological features for identification. Individuals that were challenging to identify to species level (e.g. due to distance from the camera or poor visibility) were instead labelled on a genus or family level and omitted from the analysis in the present study. We expressed the abundances of focal species as the maximum number of individuals per species in a single video frame (MaxN), a metric widely used to avoid duplicate counts of individuals that are visible in multiple frames throughout the video (Cappo, Speare, & De’Ath, 2004). In addition, video observers estimated the functional visibility in each video clip from taped markers on the bait pole and by comparing the field of view to reference images. 

All DNA extractions and setup of quantitative PCR (qPCR) reactions were conducted in a dedicated clean laboratory facility. DNA extraction was based on a modified version of a previously published eDNA extraction protocol (Spens et al., 2017). This protocol uses the DNeasy blood and tissue kit (Qiagen, Hilden, Germany) to extract eDNA directly from the Sterivex filters (Supporting information, Note S1). Each extraction batch included a negative control (extraction blank) to test for exogenous DNA contamination through the used reagents or from the laboratory. Final DNA concentration was measured for every sample using a QubitTM flourometer and the QubitTM dsDNA high sensitivity kit (ThermoFisher Scientific, Waltham, MA, USA).

Four different qPCR assays targeting the mitochondrial cytochrome b gene (cytb) were used to analyze eDNA from the collected filters (Table 1). Assays for cod, flounder and plaice were based on published species-specific primers and probes (Knudsen et al., 2019) while the assay for goldsinny wrasse was designed for the study. This assay was developed based on thorough in silico and in vitro tests to ensure specificity. In order to ensure a broader use of the assay, the in silico analysis included all closely related species of wrasse with an available genetic cytb sequence in Genbank (https://www.ncbi.nlm.nih.gov/genbank/). After database establishment, all primers and probes were aligned to the reference sequences using Geneious Prime® 2019.0.4 (https://www.geneious.com/) to verify specificity (Supporting Information, Table S2). Specificity was further tested in vitro via qPCR on 1 ng of extracted DNA from the target species, as well as Ballan wrasse (Labrus bergylta) and Corkwing wrasse (Symphodus melops), which are the only other wrasse species potentially occurring in the sampling area (Muus, Nielsen, Dahlstrom, & Nystrom, 1999; Peter Rask Møller, Natural history Museum of Denmark, personal communication). Finally, the assay targeting goldsinny wrasse was optimized for sensitivity by testing all possible combinations of different concentrations of primers (200, 400, 600 and 800 nM) and probes (200, 300 and 400 nM). The combination showing the lowest Ct-value (cycle-threshold) was chosen for final qPCR analysis.

Quantitative PCR reactions were performed in 20 µL reaction volumes containing the optimal concentration of primers and probe (Table 1) with 8 µL of TaqMan Environmental Master Mix 2.0 (ThermoFisher Scientific, Waltham, MA, USA) and 4 µL DNA template. All qPCR reactions were analyzed on a StepOne Plus Real-time PCR instrument (Life Technologies, Carlsbad, CA, USA) using the following thermal cycling profile: 5 min at 50 °C and 10 min initial denaturation at 95°C, followed by 50 cycles of 95°C for 15 seconds and 60°C for 1 min. Initially, all samples were analyzed using a positive control (IPC, Applied Biosystems) to ensure that potential PCR inhibitors did not interfere with the PCR reaction. Subsequently, all samples were run in triplicates. These included one field blank, two extraction blanks, and six NTCs per analyzed species to test for potential exogenous DNA contamination from the field or laboratory. The number of DNA copies per reaction was estimated using a dilution series. DNA Amplicons were generated by conventional PCR using the qPCR primers and were purified using the Nucleospin Gel and PCR clean-up kit (Macherey-nagel, Düren, Germany). Amplicon concentrations were measured using a QubitTM flourometer and the QubitTM dsDNA High Sensitivity Kit to enable estimation of copy numbers.

All qPCR results were compared to a standard dilution series and categorized in terms of limit of detection (LOD) and limit of quantification (LOQ). We defined the LOD as the lowest copy number detected in at least one of the triplicates in the standard dilution series, and LOQ as the lowest copy number detected in all triplicates in the standard dilution series (Ellison, English, Burns, & Keer, 2006). The sensitivity of all assays was evaluated before analyzing the extracted samples based on analysis of a standard dilution series of 1-1⨯105 copies of DNA targets per reaction. Estimation of eDNA copy numbers in the collected samples was based on a slightly higher dilution series of 4-4⨯105 copies of DNA targets. Evaluation of LOD/LOQ was based on the lowest values observed from these tests. However, most qPCR estimates of copy numbers were below the LOQ due to low DNA concentrations of the four target species in the collected samples (see Results), which precluded a quantitative data treatment (Klymus et al., 2020). Thus, as a proxy for DNA concentration, we used the proportions of positive amplifications of the individual samples (similar to Biggs et al., 2015), as it is reasonable to expect that amplification rates increase with the relative concentration of target DNA molecules (Furlan, Gleeson, Hardy, & Duncan, 2016).

Explanation of the different variables used in the two datasets are provided in the ReadMe file.