Numerical Study on Evaluation of Environmental DNA Approach for Estimating Fish Abundance and Distribution in Semi-Enclosed Bay

Abstract

Despite efforts to use environmental DNA (eDNA), accurately quantifying fish populations remains a challenge. A recent eDNA approach provided reliable estimates of coastal fish population abundance, but it was not as effective for assessing spatial distribution due to a lack of eDNA samples relative to the study area. Therefore, we conducted a numerical case study to evaluate the ability of the eDNA approach to estimate fish (Jack mackerel) abundance and distribution based on the number of eDNA samples in a semi-enclosed bay (Jinhae Bay). Our study revealed that the eDNA approach can provide reliable estimates of fish abundance, even with knowledge of the eDNA concentration in just 1% of the study area. However, for estimating spatial distribution and fish school, significant estimates were obtained only when the eDNA concentration was identified in more than 70% of the study area. Our results confirm that the eDNA approach can reflect fish abundance but has limitations in estimating fish distribution.

1. Introduction

Estimating the abundance and distribution of fish is crucial not only for the sustainable management and use of fishery resources but also for monitoring ecosystem health, supporting species conservation efforts, and guiding effective fisheries management policies [1,2]. Overfishing has emerged as a global challenge in recent decades, driven by population growth and the development of civilization, which has led to an increase in fishing activities [3,4]. Since the early 1970s, the proportion of sustainably harvested stocks has gradually decreased, with a recent estimate suggesting that only 67% of stocks are harvested sustainably [5]. To effectively manage overfishing, it is essential to conduct quantitative assessments, such as estimating fish abundance and distribution, as overfishing can cause ecological imbalances and habitat changes in the coastal system [6,7].

There exist conventional methods to estimate fish abundance, such as gill-netting, bottom-trawl, mark-recapture, echo-sounder, vessel logbooks, and multibeam sonar, which can be classified as either fishery-dependent or fishery-independent methods [8,9]. Fishery-dependent methods statistically estimate fish abundance based on fishery logs, such as vessel logbooks. These methods are less costly in terms of financial and human resources, but they involve several biases, including gear selectivity and variable fishing efforts, which can affect the quantitative evaluation of fishery resources, such as the catch per unit effort [10]. On the other hand, fishery-independent methods, such as mark-recapture and echo-sounder, are not influenced by gear selectivity since they use similar gears [11]. Moreover, these methods mainly serve scientific sampling and provide reliable data for quantitative assessment [8]. However, they often require expensive equipment and are not useful for broad-scale applications. Mark-recapture is associated with high costs because of the need for repeated capture-count-mark-release processes in terms of human and time resources [12]. Accounting for migration, mortality, and recruitment is also an ongoing challenge for these types of studies [13]. An echo-sounder requires the target strength parameter of the target species, which changes depending on the characteristics of the target fish (e.g., body size and shape), necessitating individual estimation [14]. Conventional methods may have limitations such as high costs, small-scale applicability, biases, habitat disturbance, and mortality. In particular, they may have an undetected probability for rare species, such as endangered and/or protected species.

The environmental DNA (eDNA) approach, which is an emerging method for the investigation of aquatic organisms, is cost-effective, noninvasive, and has been proposed as an alternative to conventional methods [15,16]. The eDNA methodology is less affected by investigational circumstances (e.g., accessibility and uneven distribution) and could reduce the costs related to data collection [17]. Furthermore, it has been evaluated as a way to minimize habitat disturbance because it requires only water samples for analysis [18]. In recent years, the eDNA approach has shown the ability to reliably quantify aquatic organisms. Specifically, studies have been conducted pertaining to the release and/or degradation rate of eDNA [19,20,21], biodiversity [22,23], detection [24,25], abundance [26,27,28], distribution [29,30,31,32], and comparison with conventional methods [33,34]. Early studies mainly examined the relationship of eDNA with biodiversity and the presence of species. Although more recent studies have focused on the abundance and/or distribution of aquatic species, their analyses remain a challenge because of unclear processes such as shedding, degradation, transport, and exogenous input of eDNA in the natural environment [30]. Fukaya et al. [30] proposed a novel approach while considering these processes for estimating the abundance and distribution of Jack mackerel (Trachurus japonicus) in a semi-enclosed bay. They showed that the approach could reliably estimate the abundance of Jack mackerel, but the spatial distribution was not as clear. They envisaged that the lower number of eDNA samples relative to the number of grid cells could have led to the unclear spatial distribution. Because of the absence of related studies, it is uncertain whether insufficient eDNA samples caused this disagreement. Most studies pertaining to eDNA have only been conducted since the early 2000s, and those on the estimation of abundance and distribution are very rare. While recent efforts have been made to estimate fish populations using eDNA [35,36,37], Fukaya et al.’s research [30] is the only study to estimate its distribution. A quantitative evaluation of whether eDNA can effectively estimate both abundance and distribution is now needed, and our study seeks to address this.

Herein, we evaluated the eDNA approach to estimate the abundance and distribution of Jack mackerel in a semi-enclosed bay based on a numerical case study using a number of eDNA samples relative to the study area as a simulation condition. Based on the results of the numerical case study, we evaluated the amount of eDNA samples required to significantly estimate fish abundance and distribution.

2. Materials and Methods

2.1. eDNA Approach for Estimating Fish Abundance and Distribution

The eDNA approach proposed by Fukaya et al. [30] consists of forward and backward inferences. First, the forward inference was used to calculate the eDNA concentration by inputting the current field, rate parameters, and fish density. This generates a design matrix A that can be used to calculate the fish density. Backward inference was then defined as a process to estimate fish density by utilizing matrix A and eDNA concentrations. Specifically, the fish density was estimated by multiplying A⁻¹ with an eDNA concentration vector that is interpolated to the whole cell, with limited known values of eDNA. More details regarding the eDNA approach have been reported by Fukaya et al. [30].

The tracer model required the current field, rate parameters, and fish density as inputs. The rate parameters were the eDNA shedding rate of fish and the degradation rate of eDNA. In this study, the shedding rate (9.88 × 10⁴ copies individual⁻¹ h⁻¹) and degradation rate (0.044 h⁻¹) of Jack mackerel, as introduced in the studies by Fukaya et al. [30] and Jo et al. [38], were used. In other words, the shedding and degradation rates of eDNA are constant (homogeneity). Additionally, the abundance of the target species remains stationary throughout the simulation period (demographic closure assumption [39]). For the construction of the current field, the Princeton Ocean Model (POM) was used.

2.2. Simulation of Current Field

The current field simulated in this study was aimed at Jinhae Bay, South Korea (Figure 1). The current field includes hydrodynamic processes (e.g., three-dimensional flow velocity, temperature, salinity, and diffusion coefficient) that determine the transport of eDNA in the field, such as advection and diffusion. The current field was produced using POM within the bay for approximately one month. Jinhae Bay is a tide-dominated area, and the current field does not show significant variation across seasons [40]. However, when conducting this study in other areas, seasonal changes in the current field should be taken into consideration. Specifically, the model grid was discretized using 74 × 87 horizontal grid cells with a resolution of 500 m, and the sigma (σ) coordinate was adopted for the vertical grid with 10 σ layers. The total number of grid cells was 64,380, with 24,480 aquatic cells. We then verified the tide level and tidal flow using a time series and tidal ellipse, respectively. The phase lag and amplitude of the tide level showed small errors of 0.5–7.4° and ±1.0 cm, respectively. The calculated tidal current showed good agreement between the calculated and observed values for tidal ellipse and phase. We verified temperature and salinity using the objective functions of determination coefficient (R²) and skill score (SS). The SS ranges from 0 to 1, and the closer it is to 1, the better the agreement [41]. The R² and SS values for water temperature were 0.96 and 0.99, respectively, and those for salinity were 0.74 and 0.93, respectively. More details regarding the current field have been provided by Park et al. [42].

2.3. Latent Fish Density for Simulation

In this study, randomized latent fish densities were used to represent actual fish density in the simulation. Three cases of latent fish density were considered, namely Case1, Case2, and Case3. In Case1, all values were randomized to 0–10 individuals m⁻³ (ind. m⁻³). In Case2 and Case3, high fish density (15 ind. m⁻³) was assigned to a specific point in the surface layer and bottom layer, respectively, to represent fish schools. A total of 15 latent fish densities were obtained, with five densities for each case, as shown in Figure 2.

2.4. Evaluation of the eDNA Approach

After stabilizing the tracer model for approximately one month, the eDNA concentrations and design matrix were calculated using the latent fish density, current field, and rate parameters. The eDNA concentrations collected from the field may be underestimated or overestimated, so different known value ratios (KVR) of eDNA concentrations from the top and bottom 1%, 3%, 5%, 7%, 10%, 20%, 30%, and so on, up to 90% were considered (i.e., KVR represents the number of eDNA samples relative to the number of aquatic grid cells). To prevent any bias in the known value selection, we divided the study area and selected known values from each section. Next, the selected eDNA concentrations were interpolated to the entire cells and estimated fish density by multiplying the eDNA concentration vector with A⁻¹. Any fish density estimates below zero were set to zero, and outliers were eliminated using the generalized extreme studentized deviation (GESD) method [43], with a maximum number of outliers set at 10% of the total grid cell. Finally, the estimated fish abundance and distribution were evaluated by comparing them with those of the latent condition. To assess the underestimation and overestimation of fish abundance, the estimated fish abundances were evaluated as a relative ratio (reproducibility) to latent fish abundance using the following equation:

$$Reproducibility={\displaystyle \frac{Estimated fish abundance}{Latent fish abundance}}$$

(1)

In the analysis, a reproducibility value of one indicated a perfect match, whereas values greater or less than one indicated overestimation or underestimation, respectively. The estimated fish distributions (i.e., the spatial distribution of fish densities) were evaluated by visually inspecting them and calculating the correlation coefficient (R) with the latent fish densities. We also compared the histograms and scatter plots of the estimated fish densities with those of the latent fish densities.

3. Results

3.1. Comparison of Estimation Accuracies Between Cases

The results of the evaluation for fish abundance (reproducibility) and distribution (R) with the top and bottom 1%, 50%, and 90% KVR conditions for all cases are presented in

Table 1. Evaluation of estimated fish abundance (reproducibility) and distribution (R) in top and bottom 1%, 50%, and 90% KVR conditions for all cases.

	Reproducibility						R Between Latent and Estimated Fish Densities
	Top KVR			Bottom KVR			Top KVR			Bottom KVR
	1%	50%	90%	1%	50%	90%	1%	50%	90%	1%	50%	90%
Case1-1	0.96	0.94	0.98	0.86	0.96	0.99	0.08	0.44	0.85	0.03	0.49	0.87
Case1-2	1.04	0.95	0.97	0.87	0.94	0.99	0.10	0.45	0.86	0.04	0.52	0.88
Case1-3	0.92	0.93	0.97	0.87	0.93	0.99	0.09	0.43	0.83	0.00	0.49	0.84
Case1-4	0.93	0.96	0.97	0.84	0.93	0.99	0.10	0.42	0.85	0.00	0.50	0.87
Case1-5	1.01	0.93	0.97	0.86	0.93	0.99	0.07	0.43	0.85	0.06	0.47	0.87
Case2-1	0.96	0.93	0.97	0.85	0.93	0.99	0.10	0.48	0.84	0.03	0.47	0.87
Case2-2	0.99	0.95	0.97	0.84	0.94	0.99	0.10	0.46	0.86	0.02	0.52	0.87
Case2-3	0.99	0.96	0.97	0.87	0.93	0.99	0.08	0.39	0.82	0.01	0.46	0.86
Case2-4	0.94	0.93	0.98	0.83	0.94	0.98	0.09	0.48	0.84	0.02	0.48	0.88
Case2-5	1.02	0.95	0.98	0.91	0.94	0.99	0.11	0.48	0.86	0.02	0.48	0.89
Case3-1	0.94	0.94	0.97	0.94	0.91	0.99	0.06	0.45	0.82	0.02	0.52	0.87
Case3-2	0.99	0.94	0.97	0.93	0.92	0.98	0.08	0.45	0.84	0.05	0.50	0.87
Case3-3	0.92	0.96	0.97	0.86	0.91	0.99	0.09	0.42	0.84	0.04	0.43	0.85
Case3-4	1.02	0.94	0.97	1.05	0.94	0.99	0.11	0.48	0.84	0.03	0.47	0.87
Case3-5	0.98	0.95	0.98	0.82	0.93	0.98	0.07	0.43	0.83	0.04	0.47	0.84
Mean	0.97	0.94	0.97	0.88	0.93	0.99	0.09	0.45	0.84	0.03	0.48	0.87
Std. Dev.	0.04	0.01	0.00	0.06	0.01	0.00	0.02	0.03	0.01	0.02	0.03	0.01

. The mean reproducibility (averaged across all cases) for the top and bottom 1% KVR was 0.973 and 0.877, with a standard deviation of 0.039 and 0.059, respectively. The standard deviation decreased as the KVR increased, with values being 0.004 and 0.003 at the top and bottom 90% KVR, respectively. The reproducibility for the top and bottom 1% KVR for Case3-4 was 1.020 and 1.047, respectively, which was 0.047 and 0.170 higher than the mean values. This may have been overestimated due to insufficient eDNA samples (i.e., low KVR). The difference in reproducibility between the cases was up to 0.120 and 0.229 at the top and bottom 1% KVR, respectively. It is expected that the reproducibility of the estimated fish abundance may vary depending on the fish distribution if the eDNA sample is insufficient. The R values between the latent and estimated fish densities showed a small standard deviation of 0.01–0.03, regardless of the KVR. The R values showed a range of 0.07–0.11, 0.39–0.48, and 0.82–0.86 at the top 1%, 50%, and 90% KVR, respectively. There was no significant case-dependent (i.e., fish distribution-dependent) difference between the R values.

3.2. Evaluation of an eDNA Approach for Estimating Fish Abundance

The evaluation results for fish abundance (reproducibility) with all KVR conditions for Case1-1, Case2-1, and Case3-1 are presented in

Table 2. Evaluation of fish abundance (reproducibility) and distribution (R) under all KVR conditions for Case1-1, Case2-1, and Case3-1.

		Case1-1		Case2-1		Case3-1
		Abundance (Reproducibility)	R	Abundance (Reproducibility)	R	Abundance (Reproducibility)	R
Top KVR (%)	1	2.83 × 1010 (0.96)	0.08	2.89 × 1010 (0.96)	0.10	2.81${ \times 10}^{10}$ (0.94)	0.06
	3	2.92${ \times 10}^{10}$ (0.99)	0.10	2.95${ \times 10}^{10}$ (0.97)	0.13	2.85${ \times 10}^{10}$ (0.95)	0.10
	5	2.82${ \times 10}^{10}$ (0.96)	0.13	2.91${ \times 10}^{10}$ (0.96)	0.14	2.72${ \times 10}^{10}$ (0.91)	0.14
	7	2.82${ \times 10}^{10}$ (0.96)	0.14	2.83${ \times 10}^{10}$ (0.94)	0.16	2.70${ \times 10}^{10}$ (0.90)	0.16
	10	2.86${ \times 10}^{10}$ (0.97)	0.16	2.85${ \times 10}^{10}$ (0.94)	0.20	2.73${ \times 10}^{10}$ (0.91)	0.19
	20	2.83${ \times 10}^{10}$ (0.96)	0.23	2.74${ \times 10}^{10}$ (0.91)	0.28	2.75${ \times 10}^{10}$ (0.92)	0.24
	30	2.76${ \times 10}^{10}$ (0.94)	0.32	2.75${ \times 10}^{10}$ (0.91)	0.34	2.79${ \times 10}^{10}$ (0.93)	0.30
	40	2.78${ \times 10}^{10}$ (0.94)	0.38	2.80${ \times 10}^{10}$ (0.93)	0.40	2.81${ \times 10}^{10}$ (0.94)	0.37
	50	2.78${ \times 10}^{10}$ (0.94)	0.44	2.81${ \times 10}^{10}$ (0.93)	0.48	2.82${ \times 10}^{10}$ (0.94)	0.45
	60	2.80${ \times 10}^{10}$ (0.95)	0.54	2.78${ \times 10}^{10}$ (0.92)	0.58	2.84${ \times 10}^{10}$ (0.95)	0.52
	70	2.86${ \times 10}^{10}$ (0.97)	0.64	2.87${ \times 10}^{10}$ (0.95)	0.66	2.86${ \times 10}^{10}$ (0.96)	0.62
	80	2.80${ \times 10}^{10}$ (0.95)	0.72	2.87${ \times 10}^{10}$ (0.95)	0.73	2.84${ \times 10}^{10}$ (0.95)	0.70
	90	2.88${ \times 10}^{10}$ (0.98)	0.85	2.94${ \times 10}^{10}$ (0.97)	0.84	2.91${ \times 10}^{10}$ (0.97)	0.82
Bottom KVR (%)	1	2.52${ \times 10}^{10}$ (0.86)	0.03	2.56${ \times 10}^{10}$ (0.85)	0.03	2.80${ \times 10}^{10}$ (0.94)	0.02
	3	2.79${ \times 10}^{10}$ (0.95)	0.01	2.47${ \times 10}^{10}$ (0.82)	0.05	2.70${ \times 10}^{10}$ (0.91)	0.05
	5	2.71${ \times 10}^{10}$ (0.92)	0.07	2.64${ \times 10}^{10}$ (0.87)	0.08	2.78${ \times 10}^{10}$ (0.93)	0.05
	7	2.68${ \times 10}^{10}$ (0.91)	0.10	2.71${ \times 10}^{10}$ (0.90)	0.05	2.74${ \times 10}^{10}$ (0.92)	0.09
	10	2.74${ \times 10}^{10}$ (0.93)	0.12	2.80${ \times 10}^{10}$ (0.93)	0.07	2.68${ \times 10}^{10}$ (0.90)	0.11
	20	2.68${ \times 10}^{10}$ (0.91)	0.21	2.72${ \times 10}^{10}$ (0.90)	0.17	2.72${ \times 10}^{10}$ (0.91)	0.23
	30	2.76${ \times 10}^{10}$ (0.93)	0.30	2.73${ \times 10}^{10}$ (0.90)	0.27	2.65${ \times 10}^{10}$ (0.89)	0.33
	40	2.76${ \times 10}^{10}$ (0.94)	0.41	2.74${ \times 10}^{10}$ (0.90)	0.37	2.70${ \times 10}^{10}$ (0.90)	0.43
	50	2.82${ \times 10}^{10}$ (0.96)	0.49	2.81${ \times 10}^{10}$ (0.93)	0.47	2.73${ \times 10}^{10}$ (0.91)	0.52
	60	2.83${ \times 10}^{10}$ (0.96)	0.59	2.87${ \times 10}^{10}$ (0.95)	0.60	2.78${ \times 10}^{10}$ (0.93)	0.61
	70	2.87${ \times 10}^{10}$ (0.97)	0.68	2.91${ \times 10}^{10}$ (0.96)	0.71	2.84${ \times 10}^{10}$ (0.95)	0.70
	80	2.87${ \times 10}^{10}$ (0.97)	0.79	2.96${ \times 10}^{10}$ (0.98)	0.79	2.89${ \times 10}^{10}$ (0.97)	0.79
	90	2.92${ \times 10}^{10}$ (0.99)	0.87	3.00${ \times 10}^{10}$ (0.99)	0.87	2.95${ \times 10}^{10}$ (0.99)	0.87

. The latent fish abundances for Case1-1, Case2-1, and Case3-1 were 2.95 × 10¹⁰, 3.03 × 10¹⁰, and 2.99 × 10¹⁰ individuals, respectively, which were disproportionate to the bay size of 612 km². This was caused by the high initialization of latent fish density at 0–10 ind. m⁻³. It is expected that the estimation accuracy would not differ by changing the fish density scale because our process is linear. The top KVR conditions showing reproducibility closest to the latent fish density were 3% for Case1-1 and Case2-1 and 90% for Case3-1, with the reproducibility being 0.989, 0.974, and 0.973, respectively. In the bottom KVR condition, the 90% KVR showed reproducibility closest to the latent for all cases, with the reproducibility being 0.989, 0.991, and 0.986 for Case1-1, Case2-1, and Case3-1, respectively. The approach showed high reproducibility (>0.800) not only in the 90% KVR but also in low KVR conditions. In the bottom 1% KVR condition, the reproducibility was 0.856, 0.846, and 0.938 for Case1-1, Case2-1, and Case3-1, respectively, and it gradually increased with increasing KVR. The underestimation in both KVR conditions was caused by the elimination of outliers. The reproducibility of all cases before outlier elimination was overestimated to 1.652–1.869 and 1.194–1.453 in the top and bottom 50% KVR, respectively. The eDNA approach could estimate the fish abundance with reproducibility of over 0.800 (i.e., error under 0.200), regardless of KVR and fish distribution.

3.3. Evaluation of an eDNA Approach for Estimating Fish Distribution

The estimated fish distributions by KVR for Case2-1 and Case3-1 are shown in Figure 3. The R value for the top 1% KVR of Case2-1 was 0.10, which was 0.07 higher than the bottom 1% KVR (Table 2). The difference was whether they were able to estimate the fish school. The top 1% KVR of Case2-1 was partially capable of materializing a fish school, whereas the bottom 1% KVR was not able to do so. The eDNA selection process may be the reason for this difference because the eDNA copies shed from a fish school are more likely to be selected in the top 1% KVR condition than in the bottom 1% KVR. The failure to select eDNA shed from a fish school in the bottom 1% KVR is likely what caused the materialization to fail. Unlike Case2-1, the top 1% KVR in Case3-1 could not materialize a fish school. Even under the top 5% KVR, the fish schools materialized in Case2-1, Case2-2, Case2-3, and Case2-5, but not in Case3. It depends on whether the fish school is located in the surface or bottom layer of the bay. The tidal residual current in the bottom layer of Jinhae Bay is slower than that at the surface [42]. This indicates that the transport of eDNA copies shed from the fish school in Case3 was slower than that in Case2, and there were more eDNA copies in Case3. The fish densities around the fish school in Case3 were overestimated and treated as outliers. Case3 was, therefore, not capable of materializing a fish school. We further checked that the scatter points of latent fish densities of 10–15 ind. m⁻³ (i.e., fish school) in the top and bottom 1% KVR of Case3-1 were eliminated as outliers (Figure 4). The R values increased with increasing KVR and were 0.85, 0.84, and 0.82 in the top 90% KVR of Case1-1, Case2-1, and Case3-1, respectively, and it was 0.87 in the bottom 90% KVR for all three cases. The improvement in the estimation of fish distribution according to the increase in KVR is clearly shown in Figure 3.

At 70% KVR conditions, fish schools materialized in both Case 2-1 and Case 3-1 (Figure 3). To verify if the materialization of fish schools depended on their location, the fish distribution estimates under the Top 70% KVR condition for other cases were presented in Figure 5. Regardless of the fish school’s location, they all successfully materialized. Additionally, the R values between the latent fish densities and the estimates under the Top 70% KVR condition were 0.63, 0.62, 0.65, and 0.65 for Case 2-2, Case 2-3, Case 2-4, and Case 2-5, respectively, and 0.62, 0.63, 0.62, and 0.60 for Case 3-2, Case 3-3, Case 3-4, and Case 3-5, respectively. This indicates that the number of eDNA samples is the most critical factor for estimating fish density distribution and materializing fish schools, and fish schools successfully materialized at and above the 70% KVR condition in all cases.

The comparison results of the estimated fish density histogram and scatter plot with those of the latent condition are shown in Figure 4. The histograms have been plotted as a bar graph with an interval of 0.5 ind. m⁻³. In the top 1% KVR for Case2-1, the histogram counts for 0–0.5 ind. m⁻³ were 3320 higher than those observed for the latent condition. This was because the fish densities were underestimated below 0 ind. m⁻³ and were set to 0 ind. m⁻³. The count between 2 and 5 ind. m⁻³ in the top 1% KVR for Case2-1 was 5170, while the count for the latent condition was 9213. The count estimated over 15 ind. m⁻³, which did not exist in the latent condition, was 1004. As a result, this approach may underestimate and/or overestimate fish density under low KVR conditions. These issues were resolved as the KVR increased. In the top 90% KVR of Case2-1, the count for 0–0.5 ind. m⁻³ and over 15 ind. m⁻³ decreased to 1082 and 49, respectively, and the count between 2 and 5 ind. m⁻³ increased to 8318. We also confirmed the improvement in agreement between estimated and latent fish densities according to the increase in KVR using the scatter plots. These above-mentioned results were confirmed for other cases as well.

4. Discussion

In recent years, eDNA approaches have been applied to studies on aquatic ecology and have shown their potential for quantifying aquatic organisms. To the best of our knowledge, the eDNA approach proposed by Fukaya et al. [30] is the most advanced method to estimate fish abundance in a semi-enclosed bay. However, this approach has several limitations (e.g., stationarity of fish population and homogeneity of eDNA shedding rate) that may introduce bias in the estimation of fish distribution. In their study, the fish distribution was not reasonably estimated, and one of the reasons discussed was that the number of eDNA samples was relatively smaller than the number of model grid cells. That means the number of eDNA samples is one of the major factors in determining the reasonability of estimating fish distribution. A quantitative evaluation of the relationship between the number of eDNA samples and the reasonability of the eDNA approach is needed, but the eDNA research is in its early stages, and research materials are insufficient. Therefore, we conducted a numerical study to evaluate the eDNA approach according to the number of eDNA samples. Specifically, this study followed these steps: (1) randomizing latent fish density and calculating eDNA concentration from the tracer model and latent fish density, which represents the in situ value; (2) selecting cells and assuming that we only know the eDNA concentration in the selected cells, representing the eDNA sampling process; (3) estimating fish abundance and spatial distribution using the eDNA approach; (4) finally, evaluating the eDNA approach by comparing the estimated fish densities with the latent fish densities.

From the numerical evaluation, we revealed that the eDNA approach can reasonably estimate the fish abundance regardless of the number of eDNA samples if outliers of the fish density estimates are eliminated. The reproducibility of all cases before outlier elimination was overestimated to be between 1.194 and 1.869 under 50% KVR condition; those became close to 1 after outlier elimination, with values being 0.913–0.960. This overestimation of fish abundance was consistent with the findings of Fukaya et al. [30]. They also observed an improvement in reproducibility, from 2.108 to 1.420, after omitting the cells with extremely high fish density near the fish market. We expect that further improvements in reproducibility can be achieved by eliminating outliers. These results suggest that outlier elimination is a key factor in improving the reproducibility of abundance estimation. Further studies are needed to quantitatively evaluate the relationship between outlier elimination methods and reproducibility.

We figured out that a reasonable estimation of fish distribution requires the identification of 70% or more eDNA concentrations relative to the study area, but this is practically impossible. Extensive eDNA sampling requires more costs than conventional methods, and it is worthless in engineering. Estimation of fish distribution is hard-pressed because it should clearly account for physicochemical processes of eDNA (e.g., shedding, degradation, advection, diffusion, settling, and resuspension) in the whole study area [30,44]. Those processes have been studied experimentally [19,45,46,47,48,49]; however, applying these small-scale studies to a large geographic range is very difficult. Thus, in the application of the eDNA approach to estimate a fish distribution, additional models, which could expand a few eDNA samples to inference over the whole study area, may be required. Meanwhile, Shelton et al. [50] demonstrated to expand a few eDNA samples to the whole study area using a Bayesian state-space model for modeling eDNA concentration in the coastal ocean. In the study, the eDNA concentration was defined as a function of spatial coordinates and sample depth. A combination of the Bayesian state-space model with the eDNA approach proposed by Fukaya et al. [30] may improve the reasonability of fish distribution estimation.

This study not only evaluates the performance of the eDNA approach proposed by Fukaya et al. [30] but also provides a more rigorous and quantitative assessment through a controlled numerical case study. By varying the number of eDNA samples under simulated conditions, we were able to identify key factors that affect the accuracy of fish abundance and spatial distribution estimates. The findings are crucial for extending the application of eDNA to other aquatic environments, offering a replicable framework that can be adapted to different bay types, coastal areas, or open oceans. This enhanced understanding of sampling density and spatial variability significantly contributes to the refinement of eDNA methodologies for future research.

5. Conclusions

The eDNA approach proposed by Fukaya et al. [30] successfully estimated fish abundance in a semi-enclosed bay, but it was unable to estimate spatial distribution. While they speculated that the cause was the relatively small number of eDNA samples compared to the area of the bay, the exact reason remains unclear. Therefore, to evaluate this eDNA approach, we conducted a numerical case study using varying numbers of eDNA samples as a simulation condition. The estimated abundances showed high reproducibility between 0.818 and 1.047 (if perfectly matched, reproducibility is 1), regardless of the number of eDNA samples. The approach reliably estimated the abundance, even with a small number of eDNA samples, if outliers of the fish density estimates were eliminated; however, this was not the case for the estimation of fish distribution. To obtain a correlation of over 0.6 and to materialize the fish schools regardless of fish distribution and eDNA sampling bias, it is necessary to know 70% or more eDNA concentrations relative to the study area. Therefore, it is necessary to explore other methodologies for broadly estimating eDNA concentrations with fewer samples or for estimating fish distributions in a point-to-point manner (i.e., estimating fish density at eDNA sample points). Nevertheless, this eDNA approach is useful for enhancing our ability to estimate fish abundance in semi-enclosed bays. We expect our results to provide researchers with insights into the estimation of the abundance and spatial distribution of fish using eDNA.