Assessing the breeding phenology of a threatened frog species using eDNA and automatic acoustic monitoring

Study site and environmental data

Our study site (Round Field Marsh) is a shallow upland marsh encompassing approximately 8,000 m² at the Queen’s University Biological Station (QUBS) in Eastern Ontario (44.5175° N, 76.3883° W) (Figs. 1 and S1). Around 50% of the marsh area comprises thickets and tall shrubs, while the rest is dominated by cattails (genus Typha; Fig. 1). From March 8^th to July 25^th, 2022 we recorded hourly humidity and air temperature using a HOBO MX2303 data logger (Onset Computer Corporation, Bourne, MA, USA), and water temperature using a HOBO UA-001-08 Pendant data logger (Onset Computer Corporation, Bourne, MA, US). We obtained precipitation and wind speed data from an AcuRite Atlas weather station (Chaney Instrument Co., Lake Geneva, WI, USA) located at 44.585° N, 76.304° W, 10 km away from the study site. A DJI Mavic 2 Air (DJI, Shenzhen, China) unmanned aerial vehicle (UAV) was deployed at the field site during each eDNA sampling session to characterize the field site including ice/snow cover.

Automated acoustic monitoring

We placed three Wildlife Acoustics Song Meters SM4 acoustic recorders on March 8^th, 2022, when the marsh was still frozen and covered with snow (Fig. 1). Song Meters were programmed to record the first 3 min each hour. We retrieved the data after the chorus frog breeding season (July 2022) and found one recorder had experienced technical issues with discontinuous recordings and many gaps (Song Meter C; Fig. 1C). We thus primarily used data from two recorders (Figs. 1A and 1b) with the third (Fig. 1C) for supplementary insights if needed. We listened to and visually inspected the spectrograms in Kaleidoscope Pro 5.4.6 (Wildlife Acoustics Inc, Maynard, MA, USA) and estimated the number of calling males. Chorus frogs have distinctive pulsed advertisement calls with dominant frequency at around 3.5 kHz (Fig. S2; Bee et al., 2010; Nityananda, Bee & Coleman, 2011). Inter-individual call variation allowed counting individuals both audibly and visually and up to five individuals were distinguishable with confidence (Fig. S2). If we estimated different numbers of calling males from the two recorders (due to background noise and/or the distance of calling males to the recorders), we used the higher of the two values.

We tested which environmental variables related to chorus frog calling intensity at our study site. We classified calling activity as binary for each hour (1: calling, 0: non-calling) from 0AM on March 31^st (first day of male call), to 11PM on April 12^th, the day that eDNA sampling ended. We explored the association between environmental variables and hourly calling activity using a Generalized Linear Model (GLM) with a quasibinomial family (over-dispersion) in R version 4.0.5 (R Core Team, 2021). Water temperature (°C), relative humidity (%), precipitation rate (mm/h), wind speed (km/h), and calling in the previous hour (1: calling, 0: non-calling) were included as predictor variables. Adding calling in the previous hour as another predictor variable allowed us to (partially) account for non-independence of hourly calling activity. We did not include air temperature in the model because air temperature was strongly correlated with water temperature (Spearman’s correlation = 0.7; ggcor function of the GGally package) and male chorus frogs typically call at the water surface.

eDNA assay design and validation

We designed a Taqman probe-based eDNA assay for chorus frog (P. maculata mitotype) for the cytochrome b (cytb) mitochondrial gene in Geneious v7.1.9 (Kearse et al., 2012) using parameter values recommended by Klymus et al. (2020) (see Table S1 for details). We used 100 sequences produced by Lougheed et al. (2020) and 124 sequences available in GenBank for target and co-occurring species (Table S2). The forward (5′-CGTAGCACACATCTGCCGTG-3′) and reverse primers (5′- CCGACGAAGGCTGTTGCTAT -3′) labeled with the double-quenched probe (5′ 6-FAM- ACGCAACCTGCACGCAAACGGA-Iowa Black FQ -3′) amplified a 196 bp amplicon, with at least nine mismatches with any of the nine co-occurring anuran species and five mismatches with the P. triseriata cytb mitotype sequence in Ontario (Tables S2 and S3). While other assays exist for chorus frogs in Alberta (Booker, 2016), Ontario (Beauclerc et al., 2019), and Québec (Hernandez et al., 2020; Dubois-Gagnon et al., 2022), we designed an assay that is effective for ddPCR and that could reliably distinguish the two chorus frog mitotypes in Ontario from other co-occurring anuran species in the region.

Optimal annealing temperature was 65 °C based on an initial thermal gradient experiment testing temperatures from 62.8 °C to 68.1 °C on extracted chorus frog DNA (P. maculata mitotype) on a Bio-Rad QX200™ AutoDG™ Droplet Digital™ PCR System (ddPCR) (Fig. S3). We assessed the specificity of the ddPCR assay on co-occurring species and the P. triseriata mitotype from Ontario. We prepared four replicates of each sample using the same ddPCR mix and PCR conditions, three for ddPCR detections and one for Sanger sequencing. The PCR products from the fourth replicate were purified and Sanger sequenced at Génome Québec (Montréal, Canada). None of the co-occurring species amplified, and Sanger sequencing revealed that amplifications of the P. triseriata mitotype with our P. maculata primers were due to cross contamination (Fig. S4). Finally, we validated the assay using eDNA water samples collected in 2021 (triplicates for each site) from: (1) four local marshes where chorus frogs are known to be present, and (2) five marshes and two lakes where chorus frogs do not occur (Table S4). All the amplicons from 2021 were confirmed to be the correct species through Sanger DNA sequencing following the same protocol as above.

We used an 11-point 1:2 dilution series of the chorus frog (P. maculata mitotype) synthetic DNA (Integrated DNA Technologies gBlocks™, 216 bp) with six replicates to assess the Limit of Detection (LOD) and Limit of Quantification (LOQ) of the assay. The highest concentration of this series, measured with ddPCR, was 4.64 copies/μL in the reaction mix. We define the LOD as the highest concentration below which false negatives may be expected (Klymus et al., 2020), and this is the lowest standard concentration at which the chorus frog synthetic DNA could be detected (i.e., at least one positive droplet in one replicate) (Brys et al., 2021). If a positive signal was detected in a control sample (i.e., field control, filtration control, extraction control or ddPCR no-template control; all used to detect false positives), LOD was set to three times the number of positive droplets observed in the control sample to account for possible false positive detections from cross-contamination. LOQ was defined as the lowest standard concentration of the linear range of the standard curve, before the measured concentrations plateaued (Brys et al., 2021).

Water sampling, filtration and eDNA extraction

We sampled water every 2 days starting on March 23^rd, 2022 when Round Field Marsh began to thaw (i.e., some patches of open water). No samples were collected on March 29^th and 31^st when the surface of Round Field Marsh re-froze due to cold weather (Fig. S5). We stopped sampling on April 12^th, 2022 after we heard more than three frogs calling during sampling. We collected between three and nine 1L samples of surface water depending on the amount of ice-free open water using sterile 1L bottles (Nalgene HDPE Narrow Mouth) (Fig. 1). For each sampling day, we included a field negative control (bottle filled with distilled water left open during sampling). The sampling bottles were stored individually in Ziploc™ bags in a cooler with ice packs for return to the lab for filtration. All equipment (including waders) was decontaminated with 10% bleach solution and rinsed with distilled water before and after sampling. Sterile nitrile gloves were changed between water samples.

We filtered the water samples within 6 h of sampling in a dedicated room using Polycarbonate (PCTE) membrane filters (pore size: 1.0 μm, diameter: 47 mm; Sterilitech) housed in a 47 mm in-line filter holder and a Waterra portable peristaltic pump. Prior to filtration, the filter holder and pump tubing were submerged in 10% bleach overnight and thoroughly rinsed using distilled water. The filter holder and tubing were also bleached with 1L of 10% bleach and rinsed with 1L of distilled water between each liter of water sampled. After filtration, the filters were folded with tweezers (dipped in 95% ethanol and flamed to disinfect) and stored in 2 mL tubes containing 700 μL 2% (w/v) cetrimonium bromide extraction buffer (CTAB; Turner et al., 2014). All filters were immediately stored at −20 °C until DNA extraction. When several filters were used for the same 1L sample, each filter paper was individually stored in 700 μL CTAB for independent extraction. We assessed the possibility of cross-contamination by filtering 1L distilled water as a filtration negative control.

The eDNA extractions were conducted using a modified chloroform-based method from Turner et al. (2014) by repeating the steps of adding equal volume of 24:1 chloroform: isoamyl ethanol and skipping the step of adding 5 M NaCl. We eluted the DNA in a final volume of 50 μL nuclease free water. For 24 water samples with multiple filters (two to six filters), we applied a serial elution method to maximize yield. We added 25 μL of nuclease free water to the first filter tube and incubated at 56 °C for 5 min. We then transferred the 25 µL eluate from the first filter tube to the second filter tube and incubated at 56 °C for 5 min. This was repeated for all sequential filter tubes. We conducted the serial elution for each sample twice, so that the 25 µL of eluate from the first and second serial elution added up to 50 μL. We included one no-template control (NTC) for each DNA extraction session (three in total).

eDNA detection and quantification in ddPCR

All samples were run in triplicate with NTCs and chorus frog genomic DNA positive controls (PCs) included in each ddPCR run to test for false negatives and positives, respectively. The PCR mastermix was made to a volume of 25 µL: 12.5 µL of Bio-Rad ddPCR Probe no dUTP mix (2X), 2.25 µL of forward primer (10 µM), 2.25 µL of reverse primer (10 µM), 0.625 µL of Taqman probe (10 µM), 0.175 µL of nuclease free water and 7.2 µL of DNA. PCR cycling conditions were as follows: 95 °C for 10 min, 45 cycles of 94 °C for 30 s and 65 °C for 1 min, then 98 °C for 10 min and holding at 4 °C. We used the direct quantification mode in the QX Manager 1.2 Standard Edition software (BioRad, Hercules, CA, USA). The fluorescence threshold was set for each plate (Baker et al., 2018; BioRad Droplet Digital™ PCR Application Guide) using the corresponding NTCs and PCs as reference for positive and negative droplets. Concentrations were reported as number of copies/µL (i.e., ddPCR concentration measurement) and number of copies in the original sample (1L of water, copies/L) calculated as follows:

$${\displaystyle \frac{\left({\displaystyle \frac{ddPCRconcentration\times VolumeddPCRmix}{VolumeofDNAintheddPCRmix}}\right)\times VolumeofDNAextract}{Volumeoffilteredwater}}$$

Inhibition estimation and mitigation strategy

We used a northern map turtle (Graptemys geographica) probe assay (Feng, Bulté & Lougheed, 2019) as an internal positive control (the northern map turtle does not inhabit shallow marshes) to test for inhibition in all eDNA samples in a separate ddPCR run (i.e., no multiplexing). The PCR mastermix was made to a volume of 25 µL: 12.5 µL of Bio-Rad ddPCR Probe no dUTP mix (2X), 1.125 µL of G. geographica forward primer (20 µM), 1.125 µL of G. geographica reverse primer (20 µM), 0.625 µL of G. geographica Taqman probe (10 µM), 0.425 µL of nuclease free water, 2 µL of G. geographica synthetic DNA and 7.2 µL of eDNA samples (or nuclease free water for the no-template controls). PCR cycling conditions were adapted from Feng, Bulté & Lougheed (2019): 95 °C for 10 min, 45 cycles of 94 °C for 30 s and 60 °C for 1 min, then 98 °C for 10 min with a final holding step at 4 °C. All eDNA samples were run in one replicate.

As we observed strong inhibition in some samples, we tested three mitigation methods on two of them (10L9, 12L5): (1) dilution of the eDNA samples to 1:5 and 1:10 ratio with sterile nuclease free water, (2) addition of bovine serum albumin (BSA; 20 mg/mL) to the ddPCR reaction mix (see below) of the 1:10 diluted samples, and (3) purification of the undiluted eDNA samples using the DNA Clean-Up and Concentration Micro Elute kit (Norgen Biotek Corp., Ontario, Canada) which were eluted in a final volume of 20 µL (2X 10 µL). We also purified two eDNA samples with high concentrations of chorus frog eDNA from 2021 from two local marshes in Kingston, Ontario (Cecil Graham Park and Old Mill Road sites, Table S4) to estimate the target eDNA recovery rate from the clean-up kit. We used the same PCR mastermix recipe for the 1^st and 3^rd methods (see the eDNA detection and quantification in ddPCR section). We used the following recipe for the 2^nd method: 12.5 µL of Bio-Rad ddPCR Probe no dUTP mix (2X), 1.125 µL of chorus frog forward primer (20 µM), 1.125 µL of chorus frog reverse primer (20 µM), 0.625 µL of chorus frog Taqman probe (10 µM), 1.925 µL of nuclease free water, 0.5 µL BSA and 7.2 µL of eDNA or nuclease free water for no-template controls. All samples were run in duplicate.

Male calling activities and environmental conditions

We first found open water patches on March 23^rd under consistent above-zero temperatures, but these refroze completely during a cold snap (air temperature down to −10 °C) from March 27^th to March 30^th (Fig. S5). The marsh started to re-thaw on March 31^st when the daily maximum air temperature exceeded 5 °C. After a warm rain in the afternoon of March 31^st, the first calling male was heard at 7 PM (Fig. 2A), despite the fact that over 90% of the marsh was covered with ice that morning (Fig. S5). Calls were evident from recordings from Song Meter B but not A (the data were missing in recordings from C), suggesting that the male called near the east side of the marsh where the ice first started to thaw (Figs. 1 and S5). No additional calling was heard the rest of day nor on April 1^st. After a warm rain in the afternoon of April 2^nd, two males were heard at 6 PM from both Song Meter B and C (but not A) with clearer sonograms evident in C recordings, suggesting that males were closer to C in the east side of the marsh (Fig. 1). One male was heard calling at 8 PM in recordings from all three recorders (Fig. 2A). The number of hours that males called was low during these first few days and males called only briefly on April 4^th and April 5^th (Fig. 2A). On the warm, rainy day of April 6^th, when the maximum air temperature was 18.9 °C, many males emerged and started to call (Fig. 2). The calls were louder and clearer in sonograms from Song Meter A with fewer or even no individuals heard from Song Meter B, suggesting many individuals emerged near the west side of the marsh (Fig. 1). The marsh was mostly thawed by April 8^th (Fig. S5) and the water levels were high due to the rains on April 7^th and 8^th (Fig. 2C). Male calling activities were intermittent despite the presence of many males in the marsh between April 8^th and April 11^th. Beginning on April 12^th, the chorus remained continuous. Calling activity (i.e., number of hour-time slots for which calling was recorded) in the early breeding season increased with warmer water temperatures (partial R² = 0.248, F_1,308 = 72.474, p < 0.001) and was not affected by humidity (p = 0.801), precipitation (p = 0.236), nor wind speed (p = 0.067). Calling in the previous hour was also a significant predictor (partial R² = 0.059, F_1,309 = 128.998, p < 0.001).

Sensitivity of the ddPCR assay and eDNA detections

The LOQ was set as 0.11 copies/µL (the mean concentration of the three positive 1/32 dilution replicates) (Fig. 3). Any eDNA samples with a concentration below 0.11 copies/µL were considered true detections but could not be accurately quantified. The raw reading of the LOD was 0.06 copies/uL (one positive replicate) but because this concentration is below the LOQ, we calculated the LOD (1/512 dilution) by dividing the LOQ (i.e., 1/32 dilution) by 16, which resulted in 0.007 copies/µL.

In total, 56 eDNA samples were collected and run in triplicate using ddPCR: two samples on March 23^rd, four on March 25^th, four on March 27^th, six on April 2^nd, six on April 4^th, seven on April 6^th, and nine on April 8^th, 10^th and 12^th respectively (Fig. 4). In the ddPCR runs, one replicate from locale 7 sample taken on April 6^th and two replicates from locale 5 taken on April 10^th did not generate sufficient droplets for valid interpretation (i.e., fewer than 10,000, potentially because of inhibitors, see below). Those replicates were excluded from subsequent analysis. The total number of generated droplets across all other samples and their replicates was on average 19,216, with a minimum of 10,015 droplets and a maximum of 21,758. All field, filtration, extraction and ddPCR controls were negative, except one field control replicate from April 12^th, which had two positive droplets; thus, we considered eDNA samples from April 12^th positive only if the number of positive droplets was higher than six.

The earliest positive eDNA sample was from locale 5 on April 6^th with only one positive replicate (Figs. 4 and 5E) although the concentration could not be estimated with confidence (<LOQ). No chorus frog eDNA was detected on April 8^th, despite the presence of calling males in the marsh (Fig. 2A). Two eDNA samples on April 10^th showed positive signals (10L9 and 10L3; Figs. 4, 5A and S6A). 10L3 had one positive replicate although again concentration could not be estimated with confidence (<LOQ, Fig. S6A). All three technical replicates from 10L9 showed positive signals with a mean concentration at 1.28 copies/µL (222.5 copies/L in 1L sample), which was underestimated due to PCR inhibition (Fig. 5A). Two samples on April 12^th were also positive for all three replicates (Fig. 4), locale 3 at 0.36 copies/µL (62.7 copies/L, Fig. S6B) and locale 5 at 2.48 copies/µL (431.2 copies/L, although the concentration was underestimated due to PCR inhibition, Fig. S6C). Finally, two samples (locales 4 and 8 on April 12^th) were inconclusive because the number of positive droplets was lower than six (i.e., lower than the number of positive droplets in the no-template control; Figs. 5B and S6D). In summary, one sample had positive eDNA detections on April 6^th, two on April 10^th, and two on April 12^th with two additional inconclusive samples (Fig. 4).

PCR inhibition

Using the G. geographica IPC, all samples generated enough droplets to be interpreted with confidence (>10,000 droplets). Five samples (April 6^th L3 and L7, April 10^th L9, April 12^th L4 and L5; Fig. 4) showed strong inhibition with the collapse of positive droplets that either lowered the concentration estimates or resulted in false negative detections (Figs. 5A–5C) or a higher amplitude of some positive droplets (10L9 and 12L5; Figs. 5A and S6C). Two samples with positive chorus frog detection (10L9, 12L5) and one with inconclusive detection (12L4) were strongly inhibited (Fig. 4).

We also identified 15 samples showing minor inhibition in the IPC assay (April 2^nd L8, April 4^th L1, April 6^th L1, L5, L6 and L8, April 8^th L8, April 10^th L2, L3, L5 and L8, April 12^th L1, L3, L6, L8; Fig. 4). Minor inhibition did not affect positive detections but can increase concentration estimations (over 20% variation of the IPC concentration in the NTC; Shehata et al., 2017) (Figs. 5D and 5F). They showed one or more of the following patterns: (1) elevated fluorescence of some positive droplets (Fig. 5D), (2) elevated fluorescence amplitude of the negative droplet baseline (Figs. 5E and 5F), and/or (3) broader width of the negative droplet baseline (Fig. 5F). Those minor inhibition patterns were not always observed in the corresponding target chorus frog assay. Three samples with positive chorus frog detections (6L5, 10L3, 12L3) and one that was inconclusive (12L8) showed minor inhibition in the IPC assay (Figs. 5 and S6).

Of the three mitigation methods used to assess inhibition of the locale 9 sample from April 10^th (10L9) and the locale 5 sample from April 12^th (12L5), dilution was the most effective (Fig. S7). The 1:5 and 1:10 dilutions of 10L9 had mean concentrations of 1.5 copies/µL and 0.55 copies/µL respectively, corresponding to 1,298 copies/L and 948 copies/L in the original 1L sample. The 1:5 and 1:10 dilutions of 12L5 had mean concentrations of 2.765 copies/µL and 0.615 copies/µL respectively, corresponding to 2,402 copies/L and 1,074 copies/L in the original 1L sample. With dilution, the estimated chorus frog eDNA concentrations were approximately four to six times higher in sample 10L9 and two to five times higher in sample 12L5 (Fig. S8). However, the strong inhibition patterns were still evident because the collapse of positive droplets remained in the diluted 10L9 and 12L5 samples, resulting in lower concentration estimates (Fig. S8). The combination of the 1:10 dilution with BSA generated either lower or no eDNA positive signals in both samples (Figs. S7 and S8). The commercial clean-up kit did not generate any positive eDNA signal in either of the two samples, possibly due to target eDNA being washed away and/or the inhibitors not removed completely (Figs. S7 and S8). The estimated target eDNA recovery rate of the clean-up kit was 61.14% in the Cecil Graham Park sample and 85.43% in Old Mill Road sample (Fig. S8; these two samples were not inhibited and were collected in 2021 for eDNA assay validation).

Assay design and validation

We developed and comprehensively validated a novel ddPCR probe assay for chorus frog (P. maculata mitotype) detection in wetlands in Eastern Ontario. Species-specific detection can be challenging when attempting to distinguish between closely related species, subspecies, or mitotypes (Gorički et al., 2017; Boardman et al., 2021). Insufficiently specific primers and probes between closely related species or mitotypes might result in both false positive (amplification of non-target species) and negative (competition between target and non-target DNA) detections (Wilcox et al., 2013). Our assay was rigorously designed and tested using in silico approaches, in vitro tests with P. maculata genomic DNA and genomic DNA from co-occurring anuran species plus P. triseriata, as well as in situ tests on positive eDNA samples (known presence of P. maculata in the water bodies), and negative eDNA samples (known absence of P. maculata in the water bodies). LOD and LOQ were estimated, and all positive amplicons generated during the validation process were Sanger sequenced (substantial validation at level 4 out 5 on the eDNA assay validation scale from Thalinger et al., 2021). Following further testing of specificity on samples from a chorus frog P. triseriata/P. maculata mitotype contact zone in Southern Ontario, our ddPCR also could be used to improve the geographic delimitation of the two Designable Units currently recognized in Ontario (COSEWIC, 2008).

The promise of eDNA for assessing breeding phenology of chorus frogs

Both acoustic monitoring and eDNA surveys showed that chorus frog started breeding in early spring when air temperature was low and the marsh just started thawing (Fig. 2; Whitaker, 1971). The ddPCR eDNA approach was effective in detecting chorus frogs at low abundance in a shallow marsh even when only a few males were calling. The first eDNA detection was on April 6^th, 6 days after the first male called but only 2 days after males attained higher chorus activity. This temporal lag was not unexpected; the first male called as soon as the marsh started to thaw and one would expect that it would take time for sufficient amounts of DNA to be shed into the water. Moreover, eDNA was detected at more locales and the eDNA concentrations at locale 3 and 5 increased with more calling activities in later days, suggesting that overall eDNA may be a reasonable proxy for calling assemblage size and location and that the detection of P. maculata mitotype was initially limited by the species small body mass (<1 g), its behavior (i.e., sit on floating vegetation while calling) and a lower number of active individuals in the water. Frogs and salamanders with higher body mass and more mobility within breeding water bodies would likely have higher eDNA detection rates using ddPCR (Pilliod et al., 2014; Goldberg et al., 2016; Goldberg, Strickler & Fremier, 2018). Although eDNA revealed the onset of male activities a few days later than acoustic surveys, our results imply that, with sufficient spatial sampling coverage, eDNA may be an effective approach for detecting chorus frogs before oviposition and emergence of tadpoles (Bylemans et al., 2016; Buxton et al., 2017; Dunn et al., 2017; Everts et al., 2021). By extension, eDNA might be an effective means for broad-scale surveys of other anurans that do not vocalize (e.g., newts), or whose calls are barely perceivable (e.g., African clawed frog).

Our eDNA data also revealed interesting fine spatial patterns of chorus frogs that may be difficult to determine using acoustic monitoring, especially with one acoustic recorder, wind, background noise, or overlap of male calls. For example, on the east side of the marsh (locales 3 and 5) we detected chorus frogs eDNA twice while at locales 1, 2, 6 and 7 we never detected chorus frogs. This heterogeneity within the site could be explained by either a difference in the number of active chorus frogs (i.e., more eDNA was shed in some locales and not others), or differences in water physicochemical factors that allowed longer/shorter persistence of eDNA.

Intensive spatial and temporal eDNA sampling is necessary for effective species surveys when abundances are low (Goldberg et al., 2016; Goldberg, Strickler & Fremier, 2018). We collected multiple 1L water samples every 2 days with the number of samples roughly proportional to the marsh area that was free of ice. The percentage of positive samples was less than 25% on April 6^th and 10^th but increased to almost 50% on April 12^th when many males were present and calling day and night. Low biomass of the focal taxon and the lentic wetland properties with natural barriers of dense vegetation separating water patches, reducing dispersion of shed eDNA, might be the important factors underlying the lower detection rates earlier in the breeding season (Goldberg, Strickler & Fremier, 2018; Shackleton et al., 2019). These factors should be considered in eDNA sampling design; for example, conducting a pilot study can be useful to refine sampling schemes (Goldberg et al., 2016). In addition, climate variables are important in eDNA surveys and data interpretation. For example, the breeding activities in some anuran species are triggered by rainfall (Saenz et al., 2006; Ulloa et al., 2019) but our data imply that rainfall can dilute eDNA and reduce detection probabilities. Thus, repeated sampling across time and collection of abiotic variables such as rainfall are important for interpreting eDNA results (Buxton et al., 2017; Akre et al., 2019).

Inhibition in ddPCR

PCR inhibition can be a challenge for the detection and quantification of eDNA (Sidstedt, Rådström & Hedman, 2020), especially in shallow wetlands with potentially high saprotrophic activity, low pH, and proximity to sediments that can increase the concentration of PCR inhibitor substances (Harper et al., 2019). In water samples, inhibitors can be organic (debris, fulmic and humic acids, phenol) or metallic (metal ions) (Kim et al., 1990; Abbaszadegan et al., 1993; Schrader et al., 2012). Several inhibition mechanisms can co-occur in PCR reactions (Opel, Chung & McCord, 2010; Alaeddini, 2012; Schrader et al., 2012): (i) interaction with the polymerase limiting PCR efficiency (e.g., lowering of the rate of extension); (ii) competition between inhibitors and reagents; (iii) interaction with the DNA (e.g., degradation); and (iv) interaction with the fluorescence dye (e.g., a quenching effect). While we used a small elution volume (50 µL) which could have concentrated inhibitors and DNA, we limited the presence of inhibitors by processing our samples using current best practices (multiple filters, CTAB-based storage method, and chloroform-isoamyl extraction; Hunter et al., 2019) and using a ddPCR approach, which is less sensitive to PCR inhibitors than qPCR (Verhaegen et al., 2016).

Nevertheless, our internal positive control revealed strong inhibitions in five samples that caused lower concentration estimations and potential false negatives. Minor inhibition was observed in 15 samples that did not seem to affect detection but sometimes increased concentration estimation, although this might be assay specific as the patterns were not always consistent between IPC and our target chorus frog assay. There was no obvious spatial pattern of inhibition in the marsh (Fig. 4). We found less inhibition when ice was thawing and after rain suggesting a potential dilution effect, followed by more inhibition when the marsh was fully thawed which could be explained by the release of humic substances. Warmer water temperatures could potentially activate humification (plant decay and transformation by microbial community) which would increase concentrations of humic substances in the water. Although ddPCR can be robust to high concentrations of humic acids (up to 1,000 ng/µL; Kolar, 2015), the effects of other humic compounds such as fulvic acids or humins have not been evaluated. We observed patterns of inhibition similar to Maheshwari et al. (2017) (citrus leaf petiole and fruit columella extract), Zhao et al. (2016) (copper-containing bactericides), and Kolar (2015) (metal ions). We did not investigate the source of inhibition nor the underlying mechanisms affecting eDNA detection probabilities in lentic systems, but metals and bactericides are unlikely to be the main source of inhibition in our study as one would expect that metal concentration from natural sources would be spatially homogeneous in the marsh and would not fluctuate daily (Fig. 4); moreover, there is no agricultural treatment around the site. We found samples with minor inhibition patterns sometimes had higher IPC concentrations (over 20% variation of the IPC concentration in the NTC; Shehata et al., 2017), while samples with strong inhibition always had lower IPC concentrations. This could be due to different concentrations of inhibitors in the sample causing different effects. For example, Zhao et al. (2016) found that copper had an enhancing effect at low concentration levels on both ddPCR and qPCR but an inhibitory effect at high levels. Some locales showed consistent minor inhibition (L8), with an increase or a decrease of inhibition intensity over time (e.g., L4, L5 and L3, L7, respectively; Fig. 4). However, after the heavy rains on April 7^th and 8^th, we detected no inhibition in any of the locales save for L8 (minor inhibition). The rain, by increasing water level/volume, could have diluted inhibitors and eDNA. In sum, quantification and detection of chorus frogs (P. maculata mitotype) in our study may have been impeded by inhibition issues, especially early in the season when eDNA was at low abundance in the marsh. Our results highlight the value of IPCs in eDNA studies, even when using a ddPCR approach. Although well established in qPCR studies, the use of an internal positive control to detect inhibition in environmental samples is still rare in ddPCR studies (but see Shehata et al., 2017; Everts et al., 2021).

Common solutions to PCR inhibition include adding BSA that can bind inhibitory molecules, diluting the samples, and purifying the DNA using commercial clean-up kits (Opel, Chung & McCord, 2010; Goldberg et al., 2016; Sidstedt, Rådström & Hedman, 2020; Takasaki et al., 2021). We found that diluting samples was the most effective in mitigating inhibition (i.e., increased eDNA concentrations), but the collapse of droplets was still evident. The drawback of dilution is that it can also dilute the target eDNA concentration to below the Limit of Quantification (LOQ) or Limit of Detection (LOD) and so may only work for samples with higher eDNA concentrations. Although BSA has been found to be effective in qPCR for some inhibitors (Sidstedt et al., 2015), it decreased the eDNA signals in our ddPCR tests. Finally, purification using a commercial clean-up kit resulted in varying eDNA recovery rates among samples (McKee, Spear & Pierson, 2015) and no eDNA detections in the two samples, likely due to the complete loss of target eDNA. In summary, the clean-up kit and diluting samples can reduce inhibition, but can also reduce eDNA concentration below LOD, and one should interpret such results with caution especially as these might result in false negative signals (McKee, Spear & Pierson, 2015; Mauvisseau et al., 2019).