Allometric scaling of eDNA production in stream‐dwelling brook trout ( Salvelinus fontinalis ) inferred from population size structure

Abstract

Environmental DNA (eDNA) concentration exhibits a positive correlation with organism abundance in nature, but modeling this relationship could be substantially improved by incorporating the biology of eDNA production. A recent model (Molecular Ecology, 10.1111/mec.15543) extended models of physiological allometric scaling to eDNA production, hypothesizing that brook trout eDNA production scales nonlinearly with mass as a power function with scaling coefficients <1 in lakes. To validate this hypothesis, we reanalyzed previously published data (Biological Conservation, 10.1016/j.biocon.2015.12.023) that examined the correlation between eDNA concentration and brook trout abundance in streams. We found that allometrically scaled mass (ASM) (e.g., ∑(individual mass^0.36) best described patterns of eDNA concentration across streams (r² = 0.43). ASM exhibited substantially improved model fit relative to biomass (r² = 0.31, ∆AIC = 5.19), indicating that eDNA production did not scale linearly with biomass. However, the explanatory power of ASM was comparable to density (r² = 0.40, ∆AIC = 1.25). Additionally, the optimal scaling coefficient estimated from the data (0.36) was substantially lower than that found in the previous study. Discrepancies between datasets could be attributable to ecological differences between study habitats (streams vs. lakes) or due to the exclusion of juveniles (i.e., individuals <75 mm) that can be abundant in stream environments. Nevertheless, this study adds to the growing body of literature demonstrating that individual eDNA production does not scale linearly with biomass.

INTRODUCTION

Inferring animal abundance from environmental DNA (eDNA) sampling in aquatic ecosystems is a promising approach (Pilliod et al., 2013). However, substantial unexplained variation in eDNA concentration across environments persists, with observed eDNA concentration in natural systems explaining, on average, approximately 50% of the variation in organism abundance (Yates et al., 2019). To improve our power to infer animal abundance from eDNA sampling, we need more sophisticated models that incorporate an understanding of the origins and fate of eDNA (Barnes & Turner, 2016).

The biology of eDNA production remains understudied in comparison with the significant consideration given to the role of degradation and dispersion in eDNA persistence (e.g., Barnes et al., 2014; Goldberg et al., 2018; Harrison et al., 2019; Strickler et al., 2015). The production of eDNA is likely dominated by two major physiological processes in aquatic organisms: shedding and excretion (Stewart, 2019). Metabolic waste excretion rates scale allometrically according to a power function of individual mass, with mean scaling coefficient values between 0.57 and 0.68 (Vanni & McIntyre, 2016). Similarly, eDNA shedding is likely a function of organism surface area, which also increases with body mass in fish according to a power function with scaling coefficient values between 0.59 and 0.65, depending on the species (Shea et al., 2006). As a result, eDNA production is also likely to scale nonlinearly with individual mass approximately according to the following function: [Image Omitted. See PDF]

where I = metabolic rate, I₀ = a normalization constant, M = organism body mass, and b = an allometric scaling coefficient (Brown et al., 2004; Yates et al., 2020). In functional terms, a scaling coefficient (b) value <1 predicts that larger individuals will have a lower eDNA production rate per unit of mass. Traditional metrics of abundance (density or biomass) represent extremes of a continuum of values of b, with a value of 0 corresponding to density/unit area (i.e., ∑(individual biomass⁰)/area), and a value of 1 corresponding to biomass/unit area (i.e., ∑(individual biomass¹)/area) (Post et al., 1999). At a population level, eDNA concentration can therefore be modeled using ∑(individual biomass^b) as an alternative metric of abundance(Yates et al.,2020.

Using brook trout (Salvelinus fontinalis) populations that inhabit lakes in the Rocky Mountains, Canada, Yates et al. (2020) observed that a scaling coefficient of 0.72 best explained patterns of eDNA particle concentration. This was the first study to formally propose and evaluate this hypothesis, and further studies are necessary to examine the extent to which this pattern can be extrapolated to other species, locations, or habitats. Consequently, we conducted a similar analysis on brook trout inhabiting stream systems in western Montana, USA, using data from Wilcox et al. (2016). We applied these data to test the prediction that incorporating allometric scaling coefficients (i.e., ∑(individual biomass^b), or “allometrically scaled mass” (ASM)) would substantially improve models of abundance and eDNA particle concentration relative to traditional metrics of abundance (density and biomass).

MATERIALS AND METHODS

Collection of brook trout density, biomass, and aqueous eDNA

Wilcox et al. (2016) sampled 49 sites distributed across 16 streams in the Blackfoot River and Shields River watersheds in Montana, USA, located in the Middle Rockies ecoregion (elevation between 1,200 and 2,100 m) between July and September 2013. At each site, 5 L of water was filtered at the bottom and top of each sampled stream reach (mean reach length = 108 m); for these calculations, we used estimated eDNA copy numbers from the downstream site. For further details regarding eDNA sample collection, see Wilcox et al. (2016). After eDNA sample collection, stream discharge was also estimated at the downstream site for each reach (midstream method; Hauer and Lambert (2007)), except for three sites where estimates were obtained from nearby streams. For Deep and Buck creek, stream discharge was estimated from the average discharge at three transects. After eDNA sampling, each stream reach was sampled using a backpack electrofisher to quantify the abundance of brook trout ≥75 mm total length. Abundance at 15 sites was estimated using removal methods in which two or three passes were conducted (Otis et al., 1978). At 34 sites, abundance was estimated from a single pass based on capture efficiencies estimated from the multi‐pass sites (see Wilcox et al. (2016) Appendix A for details). Total length and weight were estimated for all brook trout captured, except at eight sites at which weight was estimated from length–weight regressions obtained from the same stream or a nearby stream in the same basin.

Brook trout eDNA concentration from environmental samples was quantified using a species‐specific qPCR assay (Wilcox et al., 2013). For details on sample extraction, qPCR components and sample preparation, cycling conditions, and standard curve preparation, see (Wilcox et al., 2013,2016). Final eDNA concentrations derived from qPCR were converted to a flow‐corrected estimate (sensu Levi et al., 2019) by multiplying estimated sample eDNA concentration (copies/L) by discharge rate (L/s):[Image Omitted. See PDF]which produces an estimate of the eDNA flow rate (copies/s) at a sampling location.

Allometrically scaled mass

Allometrically scaled mass was calculated according to the following formula:[Image Omitted. See PDF]

where ${\sum }_{i=1}^{N_{\text{cap}}}mas{s}_{\text{cap}}$ is the sum of the masses captured during electrofishing surveys, b is the value of the allometric scaling coefficient, N_cap is the number of fish captured during electrofishing surveys, and $\widehat{N}$ is the estimated abundance of brook trout ≥75 mm per stream reach. This function assumes that the size structure of captured fish is representative of the population inhabiting a stream reach. The number of fish caught per reach varied substantially (0–133 individuals). The primary focus of this study is to estimate the effect of population size structure on total eDNA production—sites where no fish were captured from electrofishing efforts were therefore excluded from the analysis because no population size structure could be estimated (12 sites total). Similarly, sites where three or fewer fish were caught were also excluded from the analysis due to sample size (8 sites total), as the calculation of ASM requires a representative sample of individuals to estimate population size structure. ASM values were calculated for the remaining stream reaches (5–133 individuals captured), resulting in data from 27 sites (across 11 streams) for our analysis (Table 1). These stream reaches were 95–330 m long, so brook trout population estimates were standardized to the number of fish per 100 m. Reaches in one stream (Deep) were sampled continuously every 100 m for 1 km.

1 TableCharacteristics for the 27 stream reaches (“Site”) analyzed from Wilcox et al., 2016

Stream	Latitude (N)	Longitude (W)	Site	Captures	N c	Mean Individual Mass (g)	Reach Length (m)	Discharge (L/s)	% Slope	Flow‐Corrected eDNA Copies
Buck	46.17775	110.38330	Buck−000	6	9	22.3	100	70.8	0.031	3,327
Deep			Deep−000	25	63	29.9	100	104.7	0.019	126,268
			Deep−100	40	53	19.0	100	104.7	0.019	115,902
			Deep−200	133	155	15.3	100	104.7	0.030	136,214
			Deep−300	31	78	18.9	100	104.7	0.030	116,321
			Deep−400	39	51	32.5	100	104.7	0.030	83,027
			Deep−500	25	43	21.7	100	104.7	0.030	83,236
			Deep−600	29	73	35.3	100	104.7	0.030	43,869
			Deep−700	29	38	23.7	100	104.7	0.026	77,059
			Deep−800	50	66	20.9	100	104.7	0.026	105,537
	46.16245	110.43764	Deep−900	39	51	55.3	100	104.7	0.026	49,732
Dugout	46.18554	110.37901	Dugout−1	17	21	14.2	103	56.6	0.028	53,826
	46.18628	110.37804	Dugout−2	47	48	21.4	100	56.6	0.028	35,771
	46.18723	110.37628	Dugout−3	13	16	17.6	100	56.6	0.029	39,337
	46.19266	110.37039	Dugout−3.5	17	21	25.9	120	45.3	0.028	23,284
	46.18979	110.37347	Dugout−4	15	18	34.6	95	50.9	0.015	32,626
	46.19632	110.36652	Dugout−4.5	11	14	15.2	125	62.3	0.017	37,753
	46.19682	110.36569	Dugout−5	17	21	21.5	103	70.8	0.017	20,319
	46.19927	110.36117	Dugout−6	5	5	8.2	100	76.4	0.022	2,444
Marshall	47.28065	113.66119	Marshall	50	55	26.3	100	133.0	0.011	150,955
Horn	47.17455	113.65846	Horn	27	41	27.9	100	104.7	0.024	74,651
Placid	47.17980	113.67796	Placid	50	80	15.0	100	25.5	0.006	11,016
Shields	46.18509	110.38397	Shields	6	8	14.2	330	161.3	0.01	23,388
Smith	46.22832	110.52792	Smith	47	62	15.8	100	17.0	0.019	32,980
Unnamed−1	46.19230	110.38380	Unnamed−1	17	17	21.7	112	11.3	0.035	2,802
Unnamed−2	46.19158	110.38459	Unnamed−2	8	8	40.4	120	11.3	0.044	259
Clearwater	47.30331	113.60398	Clearwater	17	17	26.9	100	121.7	0.016	44,663

Note: Captures = the number of fish captured during electrofishing surveys, N_c = abundance of brook trout estimated from multiple capture passes.

To estimate the optimal scaling coefficient for these data, we iteratively estimated ASM using scaling coefficients ranging from 0.00 to 1.00 (by intervals of 0.01) and modeled observed eDNA concentration as a function of these ASM values using linear regressions. For each scaling coefficient value, model fit was then evaluated using AIC (Akaike, 1974). The distribution of AIC values observed for scaling coefficients ranging from 0.00 to 1.00 should follow an approximately upward parabolic shape (as observed in Yates et al. (2020)), with the value of the lowest AIC scaling coefficient representing the best fit.

To evaluate whether ASM improved the strength of the relationship between eDNA particle concentration and abundance relative to traditional estimates of abundance (i.e., density and biomass), AIC values were compared for the optimal scaling coefficient model, the density model (b = 0) and biomass model (b = 1). All analyses were performed in base R (R Development Core Team, 2017).

RESULTS

Density estimates ranged from 2 to 155 individuals/100 m and biomass estimates ranged from 0.03 to 2.81 kg/100 m. Mean individual brook trout mass at each site ranged from 8.2 to 55.3 g, and mean total length ranged from 104 to 158 mm (Figure 1).

1 Figure. Brook trout mass distributions (g) obtained from electrofishing surveys conducted on 27 stream reach study sites in Wilcox et al. (2016)

Both biomass and density exhibited a significant and positive relationship with flow‐corrected eDNA particle concentration (Table 2, Figure 2). The distribution of AIC values calculated for eDNA × ASM regressions for scaling coefficients ranging from 0.00 to 1.00 exhibited an approximately upward parabola shape, with the best‐fit value of the scaling coefficient equal to 0.36 (ASM^0.36, Figure 3). Models with scaling coefficients between 0.00 and 0.73 generated ΔAIC values <2. While ASM^0.36 represented substantial improvement as a predictor for flow‐corrected eDNA concentration relative to biomass (∆AIC = 5.19), it exhibited only marginal improvement relative to density (∆AIC = 1.25); r² values for density and ASM^0.36 exhibited similar explanatory power (Table 2).

2 TableModel results for regressions between eDNA concentration and density (fish/reach), biomass (kg/reach), and allometrically scaled mass (ASM^0.36/reach)

Model	F	p	r2	AIC	ΔAIC
Density	18.38(1,22)	<.001	0.40	644.84	1.25
Biomass	12.49(1,22)	.002	0.31	648.78	5.19
ASM0.36	20.44(1,22)	<.001	0.43	643.59	‐

2 Figure. Linear regressions for flow‐corrected brook trout eDNA concentration and three metrics of abundance in 27 study reaches: (a) density (individuals/reach, r2 = 0.40), (b) biomass (kg/reach, r2 = 0.31), and (c) allometrically scaled mass (ASM0.36/ha, r2 = 0.43) (n = 27)

3 Figure. AIC values for linear regressions for flow‐corrected brook trout eDNA concentration and allometrically scaled mass (ASM), for scaling coefficients ranging from 0.00 (density) to 1.0 (biomass). The dashed line and horizontal black bar denote range of models with ΔAIC < 2 relative to the “optimal” scaling coefficient (0.36)

DISCUSSION

We found that eDNA production scaled allometrically for brook trout in headwater mountain stream systems in the U.S. Midwest, with a scaling coefficient of 0.36 associated with the lowest AIC value. These findings were consistent with Yates et al. (2020), although the scaling coefficient was substantially lower than that found for this species in mountain lakes in the Canadian Rockies (0.72). Additionally, the optimal model for brook trout in streams in this study accounted for approximately half of the variation in eDNA concentrations relative to the lake‐based model (r² = 0.78; Yates et al., 2020). The low values for allometric scaling coefficients observed in this study imply a weak effect of individual biomass on eDNA production—in other words, fish in the stream sites appeared to produce similar quantities of eDNA among individuals, with a tendency to produce slightly more eDNA as individual biomass increased. Further, although ASM^0.36 produced better fitting models than did biomass, model fit was similar to that of simple individual density (∆AIC = 1.25, r² = 0.43 vs. 0.40 for ASM^0.36 and individual density models, respectively). Notably, physiological processes that scale with individual density represent a form of strong allometric scaling because “density” (individuals/area) is equivalent to ∑(individual biomass⁰)/area. Models with scaling coefficients of 0.00 or 0.36 both provided significantly improved model fit relative to biomass; our results provide evidence that eDNA production in wild brook trout exhibits an allometric relationship with biomass.

The low value of b estimated for the stream sites could be due to the 75 mm cutoff for inclusion in brook trout counts. Sites with smaller mean size distributions could, on average, have a larger abundance of fish <75 mm that were not formally counted during the stream surveys. Sites with smaller fish, on average, may therefore have more “uncounted” individuals and, as a result, overestimate the apparent eDNA particle contribution of fish ≥75 mm. The influence of juvenile fish on eDNA concentrations is likely to be more important in streams when compared to lakes because they constitute the typical spawning and rearing habitat for brook trout (Josephson & Youngs, 1996). Future work should be conducted to evaluate the relative contribution of juveniles to eDNA particle production because they are expected to be particularly metabolically active per unit of biomass (Maruyama et al., 2014).

The population size structure of brook trout also differed between streams and lakes. Mean book trout mass ranged from 14 to 55 g between the 27 sampling sites within and among the 11 headwater streams. In the nine lakes in Yates et al. (2020), fish were often substantially larger; mean individual masses ranging from 43 to 405 g. In small cold headwater stream systems, the opportunity for growth can often be relatively limited for brook trout (Hazzard, 1932; Hutchings, 1993). Accounting for size structure might therefore be most important when the magnitude and variation in population size structure is large.

Finally, differences in collection methodology could account for discrepancies between the two study. A myriad of differences between the studies, including (but not limited to) different filter pore sizes, filter preservation methods, DNA extraction methods, qPCR protocols, and timing of sample collection, could have caused unknown biases in results. Future research to optimize eDNA collection, extraction, and analysis is crucial to better standardize and compare results across eDNA studies (Hinlo et al., 2017; Tsuji et al., 2019).

Overall, explanatory power of the stream model was less than in the lake model, regardless of whether allometry was incorporated. Although it should be noted that the lake model was based on a smaller sample size (n = 9 compared to n = 27 sampling sites across 11 streams), we suspect that the discrepancy in model explanatory power observed is due to characteristics of the systems affecting eDNA deposition and retention. Environmental DNA deposition is likely slower and less complex or variable in lentic systems; horizontal transport tends to be limited as eDNA particles largely diffuse downwards as they slowly settle (Ghosal et al., 2018; Goldberg et al., 2018). By contrast, in lotic systems eDNA deposition is an extremely complex process that depends on multiple environmental variables such as discharge, water velocity, eDNA input from individuals upstream, width, depth, and channel roughness (Robinson et al., 2019; Shogren et al., 2017; Wilcox et al., 2016). This complexity likely contributed substantially to the increased variability in eDNA particle concentrations observed across the stream study sites. While flow‐corrected estimates of eDNA production (Levi et al., 2019) can account for some of this variability, eDNA particle transport and deposition are ultimately much more complex in lotic systems relative to lentic systems. Predictive models of organism abundance based on observed eDNA particle concentrations may therefore be more precise in lentic systems compared to lotic environments. Future work could potentially address this complexity by incorporating the impact of site environmental characteristics as covariates on eDNA detection and concentration (Mackenzie et al., 2002; Wilcox et al., 2016).

As shown here, there are likely ample opportunities to test for allometric scaling in eDNA production using existing datasets. Many studies have examined the relationship between organism abundance and the concentration of eDNA in a natural environment, typically quantifying abundance using density and biomass (e.g., Doi et al., 2017; Erickson et al., 2016; Nevers et al., 2018; Pilliod et al., 2013). As long as both density and individual biomass data have been recorded, it would be possible to explore the extent to which allometric scaling coefficients improve the correlation between eDNA particle concentration and abundance. To this end, we also recommend that future studies examining the relationship between abundance and eDNA concentration record individual biomass data to ensure that allometry in eDNA production can be accounted for. Nevertheless, our findings contribute to a growing number of studies demonstrating that individual eDNA production likely does not scale linearly with individual biomass (Maruyama et al., 2014; Takeuchi et al., 2019; Yates et al., 2020).

ACKNOWLEDGMENTS

We thank Brad Shepard, Adam Sepulveda, Winsor Lowe, Stephen Jane, Andrew Whiteley, Ladd Knotek, Will Schreck, and field technicians from Montana Fish, Wildlife, and Parks, the Wildlife Conservation Society, and U.S. Geological Survey who contributed to collecting and interpreting data from the original Wilcox et al. (2016) study. We would also like to thank an associate editor and two anonymous reviewers, whose comments substantially improved the manuscript. MCY was funded by a FRQNT postdoctoral scholarship.

AUTHOR CONTRIBUTIONS

TMW, KSM, MKY, and MKS collected eDNA samples and abundance data. MCY conducted statistical analyses and wrote the first draft of the manuscript, and all authors contributed substantially to subsequent drafts.

DATA AVAILABILITY STATEMENT

eDNA, population size structure, and site habitat data are deposited in the Dryad Digital repository at: Yates, Matthew et al. (2020), Data from: Allometric scaling of eDNA production in stream‐dwelling brook trout (Salvelinus fontinalis) inferred from population, size structure, Dryad Digital Repository, https://doi.org/10.5061/dryad.bvq83bk6v