Multidimensional Drivers of Fish Community Assembly Across Seasonal and Hydrographic Gradients in the Yangtze River Estuary and Adjacent East China Sea: Insights from eDNA Analyses

TL;DR

This study uses environmental DNA to explore how seasonal and environmental factors shape fish communities in the Yangtze River Estuary and East China Sea.

Contribution

The study reveals how multidimensional environmental and seasonal gradients jointly influence fish community assembly patterns in a dynamic estuarine-coastal system.

Findings

Surface and bottom fish communities differ in species composition and assembly patterns.

The Mixed Transitional Water zone supports the most diverse and interactive fish assemblages.

Spring communities show higher diversity and complex trophic networks compared to autumn.

Abstract

Fish communities in estuaries and adjacent coastal seas are shaped by strong environmental gradients and seasonal variability, yet how these factors jointly influence community organization remains incompletely understood. In this study, we investigated fish communities in the Yangtze River Estuary and Adjacent East China Sea using environmental DNA, a non-invasive approach that allows species to be detected from water samples without direct capture. Fish communities were compared across water depths, hydrographic zones, and seasons. We found that surface and bottom communities differed not only in species composition but also in their assembly patterns. In addition, pronounced spatial variation was observed across estuarine-offshore transition areas, which supported higher diversity and more complex species associations. Seasonal changes were also evident, with spring communities showing greater complexity than those observed in autumn. Together, these results indicate that water depth, hydrological setting, and seasonal dynamics interact to structure fish communities in estuarine-coastal systems. This study demonstrates the value of environmental DNA for studying fish community patterns across multiple spatial and temporal scales and provides ecological insights that are relevant for understanding biodiversity organization in highly dynamic coastal environments.

Marine fish communities in the Yangtze River Estuary and Adjacent East China Sea (YRE-ECS) are subject to complex environmental gradients; however, their multidimensional assembly mechanisms remain insufficiently resolved. Here, we integrated environmental DNA (eDNA) metabarcoding, co-occurrence network analysis, and environmental profiling to examine fish community structure across vertical layers, hydrographic zones, and seasons. Vertically, surface communities dominated by pelagic-associated Perciformes and Clupeiformes showed more variable assembly patterns, whereas bottom communities enriched in Gobiiformes and Pleuronectiformes were more strongly associated with temperature and dissolved oxygen. Horizontally, among three zones delineated by salinity and hydrographic characteristics, the Mixed Transitional Water (MTW) supported the most diverse and interactive assemblages and functioned as an ecological connector between estuarine (EHSW) and offshore (OWSW) waters. Seasonally, community structure shifted markedly: spring communities exhibited higher diversity and denser trophic networks supported by zooplankton-rich, phototrophic plankton (e.g., Arthropoda, Bacillariophyta), whereas autumn communities were simpler, dominated by Chlorophyta and microbial taxa, with fish assemblages showing increased modularity and reliance on fewer planktonic groups. This seasonal pattern suggests a transition from diversified energy pathways to more constrained trophic coupling. βNTI and Mantel analyses jointly revealed a stratified environment-response-feedback framework driving community differentiation through combined stochastic and deterministic mechanisms. These findings highlight the importance of integrated spatial-temporal monitoring and suggest that protecting transitional zones and spring food-web integrity is critical for ecosystem resilience in the YRE-ECS.

Figures

Full text

1. Introduction

Global climate change and human activities have substantially altered marine fish community structure and dynamics, highlighting the importance of spatiotemporal analyses for assessing ecosystem stability and functional resilience [1]. The Yangtze River Estuary (YRE), located at the land–sea interface of the East China Sea (ECS), is a highly dynamic estuarine system characterized by strong seasonal hydrodynamics, pronounced salinity gradients, and substantial freshwater and nutrient inputs from the Yangtze River [2,3]. Regulated by the East Asian monsoon system, these environmental features generate marked spatial and temporal heterogeneity in fish community composition and trophic structure [4,5]. Beyond climate-induced variability, fish communities in the YRE have been strongly shaped by long-term anthropogenic disturbances, particularly intensive fishing pressure and associated trophic alterations [6]. For instance, high fishing pressure has led to the depletion of large, high-trophic-level fishes and their replacement by smaller, lower-trophic-level species. This pattern is reflected in the well-documented decline in the mean trophic level of fishery catches in the East China Sea. Such changes are commonly described as the “fishing down” effect, which simplifies food-web structure and may reduce ecosystem stability by diminishing the functional role of top predators [7]. The implementation of the Yangtze River Basin Ten-Year Fishing Ban in 2021 represents a landmark management measure that has substantially reduced fishing pressure, thereby creating unique opportunity to assess how fish community diversity, richness, and structure respond to the combined effects of environmental variability and fisheries management [8].

However, studies on estuarine fish communities in the Yangtze River Estuary and Adjacent East China Sea (YRE-ECS) have emphasized diversity indices and abiotic drivers, often overlooking species interactions and indirect resource-mediated effects [9,10]. Such limitations constrain process-oriented interpretations of community-level responses to both climate change and conservation policies. Consequently, an integrative assessment that explicitly considers spatiotemporal heterogeneity and ecological interactions is essential for evaluating ecosystem responses to ongoing environmental change and management interventions in the YRE-ECS. Recent advances in environmental DNA (eDNA) technology have provided an efficient, non-invasive approach for assessing fish community composition with high taxonomic resolution [11]. Compared with traditional fisheries surveys and gut content analyses, which provide critical information on biomass and trophic relationships but are often constrained by sampling efficiency, depth accessibility, and seasonal logistics [12], eDNA metabarcoding enables standardized and repeatable sampling across water layers, hydrographic zones, and seasons [13]. This advantage is particularly relevant in estuarine–coastal environments characterized by strong vertical stratification, pronounced horizontal gradients, and marked seasonal variability, where conventional sampling can be difficult to implement consistently. It is important to note that correlation-based networks inferred from eDNA primarily describe co-occurrence patterns and do not provide direct evidence of trophic interactions [14]. eDNA-based surveys offer clear advantages in detecting rare taxa, early life stages, and spatially heterogeneous assemblages compared to traditional methods [15]. Moreover, when integrated with network-based approaches, eDNA data can be used to characterize statistical co-occurrence structure (e.g., modularity and connectivity patterns) and to identify highly connected taxa (network hubs) that may serve as indicators of community organization [16,17]. Despite these methodological advances, how fish co-occurrence networks respond to multidimensional environmental gradients, especially spatial replacement and seasonal turnover, remains poorly quantified in dynamic coastal systems.

Notably, phytoplankton, bacteria, and other planktonic primary producers serve as fundamental components of marine food webs. Their seasonal shifts in abundance and composition may be associated with fish community dynamics in ways consistent with potential bottom–up influences, particularly in nutrient-limited or production-fluctuating environments [18,19]. For example, phytoplankton blooms in spring typically enhance primary productivity and trophic connectivity, while autumn often features increased heterotrophic bacterial abundance and simplified food web structures [20,21]. Although key aspects of trophic dynamics and planktonic community variability relevant to fish assemblages in Chinese coastal waters have been examined in previous studies [22], the ecological effects of recent large-scale fisheries management interventions are still emerging. Following the implementation of the Yangtze River Basin Ten-Year Fishing Ban, the Yangtze River Estuary and adjacent coastal waters appear to be undergoing an early stage of ecological reorganization, characterized by changes in fish community diversity and species composition, together with concurrent variability in environmental conditions and planktonic community structure [23,24]. Under this new management context, how seasonal shifts in basal planktonic resources interact with multidimensional environmental gradients to shape fish community structure across multiple spatial dimensions has not yet been systematically evaluated.

Against this background, this study focuses on the nearshore waters of the Yangtze River Estuary and Adjacent East China Sea, integrating eDNA-based biodiversity surveys, co-occurrence network analysis, and multidimensional community structure assessment. We aimed to characterize the spatial (vertical and horizontal) and temporal (seasonal) variation in fish community diversity and composition, and to examine their associations with planktonic basal resources. Specifically, we address the following questions: (1) How does fish community structure vary along multidimensional environmental gradients? (2) Are seasonal shifts in plankton communities associated with seasonal responses in fish community structure in ways consistent with potential bottom–up influences? Our study contributes to a deeper understanding of marine community assembly patterns and their environmental associations in dynamic shelf environments and offers a scientific basis for coastal biodiversity conservation and fishery management strategies.

2. Materials and Methods

2.1. Sampling Areas and Sample Preparation

Water samples were collected from 23 sampling stations in the Yangtze River Estuary and Adjacent East China Sea in April and November 2024 (30°45′–31°44′ N, 122°16′–123°44′ E). The detailed information on the sampling stations is provided in Table S1. At each station, surface (approximately 0.5 m) and bottom (20–50 m) water samples were collected in three biological replicates. Each replicate corresponded to one filter membrane, and the standard volume filtered per membrane was 1–2 L. A field blank (filtered deionized water under identical conditions) was included at each site to monitor potential contamination. Filtration was conducted immediately after sampling using 0.45 μm pore-size glass fiber filters (47 mm diameter) and a vacuum filtration system. Each sample was filtered until visible particulate matter accumulated or until 1–2 L of water had passed through. The filters were folded, placed into sterile 2–10 mL centrifuge tubes, and temporarily stored at −20 °C before DNA extraction.

2.2. Horizontal Water Mass Classification

Based on salinity and hydrographic characteristics [25], sampling stations were horizontally classified into three distinct water masses (Figure 1): (1) Estuarine Halocline Stratified Water (EHSW), located inside the estuarine front and characterized by a pronounced halocline with a surface and bottom salinity difference of (ΔS ≈ 0.27) (ΔS represents the difference in salinity between the bottom and surface layers); (2) Mixed Transitional Water (MTW), situated between two frontal zones and showing vertically mixed salinity with (ΔS ≈ 0.03); and (3) Outer Warm Current Saline Water (OWSW), distributed outside the plume front, dominated by the warm current with high and homogeneous salinity (ΔS ≈ 0.04). To assess temporal variation, samples from corresponding stations were collected in both April and November for comparative analyses of seasonal differences in fish community composition.

2.3. Environmental Parameter Measurements

Water temperature and salinity were measured in situ using a conductivity-temperature-depth (CTD) multi-parameter profiler in strict compliance with national and industry standards GB/T 12763.2-2007 [26] and GB 17378.4-2007 [27]. The pH value, dissolved oxygen (DO) concentration, chemical oxygen demand (COD), and suspended solids (SS) were determined in accordance with GB 17378.4-2007. Total nitrogen (TN) and total phosphorus (TP) concentrations were determined by flow injection analysis in accordance with HY/T 147.1-2013 [28]. Water samples were digested before automated analysis for absorbance detection and concentration calculation. All measurements were conducted under rigorous quality control, including the analysis of procedural blanks and duplicate samples. These measures were implemented to ensure the reliability, precision, and comparability of data across different sampling stations and periods.

2.4. Environmental DNA Extraction

Environmental DNA (eDNA) was extracted using the DNeasy Blood & Tissue Kit (Qiagen, Hilden, Germany) following the manufacturer’s protocol with slight modifications. The extracted eDNA was homogenized and aliquoted into 20 μL subsamples. DNA quality was examined using agarose gel electrophoresis and a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA), and both DNA concentration and A260/280 ratio were recorded. All qualified eDNA extracts were stored at −80 °C for subsequent molecular analyses.

2.5. PCR Amplification and High-Throughput Sequencing

Fish eDNA was amplified using the universal primer pair MiFish-U-F and MiFish-U-R, which targets a hypervariable region of the mitochondrial 12S rRNA gene with an expected fragment length of approximately 170 bp. The primer sequences were MiFish-U-F (5′-GTCGGTAAAACTCGTGCCAGC-3′) and MiFish-U-R (5′-CATAGTGGGGTATCTAATCCCAGTTTG-3′). To simultaneously detect non-fish taxa commonly present in coastal assemblages, including algae, arthropods, mollusks, and other metazoans, an additional COI metabarcoding primer set was employed. The primer pair mlCOIintF (5′-GGWACWGGWTGAACWGTWTAYCCYCC-3′) and jgHCO2198 (5′-TAIACYTCIGGRTGICCRAARAAYCA-3′) amplifies a 313 bp fragment of the mitochondrial COI gene. PCR reactions were performed in a total volume of 25 μL, containing 12.5 μL of 2× San Taq PCR Mix (Sangon Biotech, Shanghai, China), 1 μL of each primer (5 μmol L^−1^) with unique barcodes, 1 μL of template DNA (10 ng μL^−1^), and nuclease-free water to volume. The thermal cycling program consisted of an initial denaturation at 94 °C for 4 min, followed by 34 cycles of 94 °C for 30 s, 65 °C for 30 s, and 72 °C for 1 min, and a final extension at 72 °C for 10 min. Amplification products were visualized on 2% agarose gels, and positive amplicons were purified using the MiniBEST Agarose Gel DNA Extraction Kit Ver. 4.0 (Takara, Dalian, China). Negative controls were included in all PCR assays to detect potential contamination from reagents or the environment. Purified PCR products were pooled in equimolar amounts to achieve the desired sequencing depth and were subjected to paired-end sequencing on the Illumina MiSeq platform (Shanghai Meiji Biomedical Technology Co., Ltd., Shanghai, China).

2.6. Data Processing and Analysis

All statistical analyses of environmental variables and biodiversity data were conducted using R version 4.2.2. The Shapiro–Wilk test and Levene’s test were used to assess normality and homogeneity of variance, respectively. To evaluate differences in α-diversity indices (Shannon diversity, Pielou’s evenness, and ACE richness) and environmental factors across groups (e.g., vertical layers, water masses, and seasons), the non-parametric Kruskal–Wallis test was used. When significant effects were detected (p < 0.05), post hoc pairwise comparisons were conducted using Dunn’s test with Benjamini–Hochberg false discovery rate (FDR) adjustment for multiple testing, and adjusted p-values < 0.05 were considered significant. Community β-diversity was assessed using Bray–Curtis distance matrices and visualized through Principal Coordinates Analysis (PCoA) and Non-metric Multidimensional Scaling (NMDS). The significance of community structural differences was tested using Analysis of Similarities (ANOSIM) and PERMANOVA (adonis2, vegan).

To infer mechanisms of community assembly, the β Nearest Taxon Index (βNTI) was calculated to distinguish between deterministic and stochastic processes. Co-occurrence networks were constructed using pairwise Spearman’s rank correlation coefficients (ρ) among the top 50 most abundant fish species across all samples. An edge was retained only when the correlation exceeded a stringent threshold (|ρ| > 0.6) and was statistically significant (p < 0.01). Given the large number of pairwise tests, we applied these strict thresholds to reduce spurious edges; however, p-values were not further adjusted for multiple comparisons. Therefore, the inferred co-occurrence relationships should be interpreted cautiously, correlation-based networks derived from eDNA data were interpreted as patterns of co-occurrence rather than direct ecological interactions. Key topological parameters including the number of nodes, average degree, graph density, modularity, and clustering coefficient were extracted for comparison. All statistical plots were generated using OriginPro 2022b, and figures were finalized and assembled using Adobe Illustrator 2023 for publication.

To balance ecological interpretability and taxonomic robustness, prey community composition was summarized at the phylum level to characterize broad trophic structure and basal resource regimes. This approach was adopted because fine-resolution taxonomic assignments based on COI metabarcoding can be uneven across diverse planktonic groups in coastal waters, and phylum-level patterns provide a more stable representation of dominant resource environments across vertical layers, water masses, and seasons. Although trophic values can vary within a phylum, summarizing at the phylum level helps capture the general functional roles of plankton groups (e.g., primary producers, grazers) without over-interpreting fine-scale trophic relationships. In contrast, genus-level analyses were conducted to identify key prey taxa exhibiting significant seasonal shifts, providing complementary insights into finer-scale trophic signals. Pearson correlation coefficients (two-tailed) were calculated to assess associations between prey phyla and fish orders for each sampling month (April and November). For each month, all pairwise prey–fish correlations were tested (10 prey phyla × 15 fish orders = 150 tests). To reduce the risk of Type I errors due to multiple comparisons, p-values were adjusted using the Benjamini–Hochberg false discovery rate (FDR) procedure. FDR-adjusted q-values were reported, and statistical significance was evaluated at q < 0.05.

3. Results

3.1. Spatial and Temporal Patterns of Fish Community Diversity and Composition

An overview of fish community diversity and abundance metrics, including ASV number, total ASV reads, and α-diversity indices (Shannon diversity, Pielou’s evenness, and ACE richness), across spatial and temporal classifications is provided in Table S2, indicating broadly comparable values between vertical layers and seasons, but more pronounced variability among water masses. At the vertical scale (Figure 2a), the Shannon index was slightly higher in bottom samples compared to surface samples, indicating a marginally greater evenness and diversity in bottom waters. This interpretation was further supported by Pielou’s evenness index, which showed a consistent pattern with Shannon diversity, confirming that the observed differences primarily reflected changes in community evenness rather than richness alone. In contrast, the ACE index, a proxy for species richness, was higher in bottom waters than in surface waters, but with lower variability, suggesting that while species richness was slightly higher in the bottom waters, surface samples exhibited greater heterogeneity, as also evident from greater dispersion and multiple outliers. These patterns suggest that although surface waters are often associated with higher productivity, the bottom layer may exhibit a more even distribution of detected species, potentially reflecting differences in habitat structure or sampling context. The Bray–Curtis-based β-diversity analysis (PCoA) clearly separated surface and bottom communities, and ANOSIM analysis confirmed this difference was statistically significant (R = 0.1710, p = 0.003).

At the horizontal scale (Figure 2b), distinct gradients in α-diversity were observed among water masses. MTW exhibited the highest Shannon diversity, indicating a more even community structure, consistent with higher Pielou’s evenness, while EHSW showed the largest variability in ACE, suggesting greater spatial heterogeneity within the estuarine front. OWSW showed relatively low and narrowly distributed diversity values, indicating a comparatively homogeneous community structure. PCoA results indicated partial separation among the three water masses, although substantial overlap was evident. ANOSIM analysis further confirmed that differences among water masses were statistically significant but of small effect size (R = 0.0738, p = 0.04). Pairwise comparisons further revealed that community composition in EHSW differed significantly from both MTW (p < 0.05) and OWSW (p < 0.01), whereas no significant difference was detected between MTW and OWSW, indicating a higher degree of similarity between the two offshore-influenced water masses. These results suggest that contrasts involving the inner estuarine water mass contributed more strongly to the observed horizontal differences, while transitional and outer water masses exhibited a higher degree of overlap in community composition.

At the seasonal scale (Figure 2c), both Shannon and ACE indices were higher in April than in November, indicating greater richness and evenness in spring, a pattern consistent with Pielou’s evenness results. PCoA ordination showed a pronounced seasonal separation, with April and November samples forming distinct clusters (R = 0.4041, p = 0.001). The broader dispersion of April samples suggests greater environmental heterogeneity in spring, whereas the more compact distribution of November samples reflects a more stable community composition during autumn. Overall, seasonal variation had the strongest influence on fish community structure when compared with vertical stratification and differences in water masses.

3.2. Integrated Analysis of Community Composition, Species Overlap, and Network Structure

Table 1 provides an overview of the taxonomic composition and ecological traits of fish species identified by eDNA sequencing, with 21 orders, 34 families, 62 genera, and 72 species. At the order level (Figure 3a–c), fish community composition exhibited clear differentiation across vertical, horizontal, and seasonal dimensions, with Gobiiformes, Perciformes, Pleuronectiformes, and Clupeiformes consistently representing the four most dominant orders (each accounting for >10% in at least one grouping). Vertically, Gobiiformes and Pleuronectiformes were more abundant in bottom waters (47.6% and 15.2%, respectively) compared with surface waters (37.6% and 12.2%) (p < 0.05). In contrast, Perciformes and Clupeiformes showed higher proportions in the surface waters (27.1% and 13.8%, respectively) than those (16.7% and 13%, respectively) in the bottom waters (p > 0.05). These differences demonstrate a clear separation in the relative dominance of benthic-associated and pelagic-associated orders between vertical layers. Horizontally, Gobiiformes were most abundant in MTW (43.4%), followed by EHSW (37.9%) and OWSW (37.7%) (p > 0.05). Perciformes decreased from EHSW (34.8%) to MTW (24.9%) and reached the lowest level in OWSW (16.4%) (p < 0.05). In contrast, Pleuronectiformes increased from EHSW (7.8%) to MTW (13.6%) and were most abundant in OWSW (22.8%), while Clupeiformes were lowest in MTW (11.2%) and higher in EHSW (13.3%) and OWSW (16.0%) (p > 0.05). This pattern suggests a gradual taxonomic transition from estuarine to offshore waters, with MTW representing an intermediate assemblage influenced by both estuarine and offshore processes. Seasonally, Gobiiformes (54.3%) and Pleuronectiformes (18.5%) were more abundant in April than in November (23.8% and 6.3%, respectively), whereas Perciformes dominated in November (47.9%) relative to April (8.2%) (p < 0.05). Clupeiformes showed comparable contributions in April (14.1%) and November (12.7%) (p > 0.05). These changes reveal pronounced seasonal turnover in the dominance structure of major fish orders.

At the species level, Venn diagrams based on the most abundant taxa revealed distinct patterns of overlap and uniqueness across spatial and temporal dimensions. Vertically (Figure 3a), surface and bottom layers shared 36 species, accounting for 50.0% of the total species detected (n = 72). The surface layer contained 26 unique species (36.1%), whereas only 10 species (13.9%) were exclusive to the bottom layer, indicating a higher degree of species exclusivity in surface waters. This higher proportion of surface-specific species likely reflects greater environmental variability and transient occurrences of species in upper waters. Across water masses (Figure 3b), 27 species (37.5%) were shared among EHSW, MTW, and OWSW. In contrast, 20 species (27.8%) were unique to EHSW, substantially exceeding the numbers unique to MTW (3 species, 4.2%) and OWSW (6 species, 8.3%), highlighting stronger species differentiation in estuarine-influenced waters. Seasonally (Figure 3c), 31 species (43.1%) were shared between April and November. However, November exhibited more unique species (27 species, 37.5%) than April (14 species, 19.4%), indicating pronounced seasonal turnover in species composition. These results demonstrate that species-level assemblages exhibit substantial spatial heterogeneity and seasonal variability, with differences arising from both habitat-related separation and temporal turnover.

The co-occurrence network of the top 50 fish species consisted of 50 nodes and 47 edges, with an average degree of 1.88 and a graph density of 0.038, indicating weak overall connectivity and a highly sparse network (i.e., most species were linked to fewer than two others and only ~3.8% of all possible links were present), and was fragmented into multiple connected components. The network showed an average clustering coefficient of 0.671 and a modularity value of 0.860, reflecting a modular network structure with limited overall connectivity. Figure 4 shows that the temporal co-occurrence network of the top 50 fish species was primarily dominated by Perciformes (15%), Gobiiformes (14%), and Scombriformes (12%), which together accounted for over 40% of all nodes in the network. Pleuronectiformes and Clupeiformes each contributed an additional 10%, indicating that a limited number of major orders consistently dominated the composition of nodes within the temporal co-occurrence network. In addition, the temporal co-occurrence network of the top 50 fish species was organized into several well-defined modules, with Modules I (16%), II (14%), and III (12%) representing the three largest components. This modular organization indicates that species interactions at the temporal scale were structured around a limited number of dominant network modules rather than being evenly distributed across the assemblage. Within Module I, Parachaeturichthys polynema and Amblychaeturichthys hexanema exhibited strong associations with multiple co-occurring species, as reflected by their high degree values. Similarly, Scomber japonicus in Module II and Nibea albiflora in Module III showed elevated connectivity within their respective modules, identifying them as highly connected taxa that contribute to within-module connectivity in the co-occurrence network, indicating central positions within the observed co-occurrence structure at the temporal scale. Overall, the temporal co-occurrence network displayed a clearly modular structure, with species associations unevenly distributed across multiple modules. It should be emphasized that this co-occurrence network, inferred from correlation-based eDNA data, describes statistical association patterns among species rather than direct ecological interactions or trophic relationships. This modular organization reflects non-random patterns of species co-occurrence at the temporal scale.

3.3. Assembly Processes of Fish Communities Shaped by Environmental Gradients

At the spatial scale, fish communities exhibited clear structural differentiation and contrasting assembly-related patterns along both vertical and horizontal dimensions (Figure 5a–d). Vertically, βNTI analysis revealed contrasting assembly patterns between surface and bottom communities (Figure 5a). Surface communities showed higher and more dispersed βNTI values, frequently extending beyond the ±2 threshold, suggesting that deterministic processes consistent with relatively stronger signals in shaping community structure, although these processes were not directly attributable to the measured environmental variables. In contrast, βNTI values for bottom communities were largely constrained within ±2, indicating a greater influence of stochastic assembly processes. Mantel test results (Figure 5b) further showed that bottom communities were significantly correlated with water temperature (r = 0.270, p = 0.009) and dissolved oxygen (r = 0.218, p = 0.025), indicating that bottom-layer community variation was statistically associated with environmental gradients, in combination with stochastic processes, contributes to the formation of bottom-layer assemblages. In contrast, surface communities did not exhibit significant correlations with the measured environmental variables, despite evidence of deterministic structuring from βNTI patterns, suggesting that the observed non-random assembly patterns in surface waters may be influenced by unmeasured or external factors rather than in situ physicochemical gradients. This discrepancy indicates that surface communities may be more strongly influenced by external disturbances, such as hydrodynamic forcing and anthropogenic activities, rather than by local environmental gradients alone.

Horizontally, fish communities also displayed pronounced differences among water masses. βNTI values were highest in EHSW, indicating strong environmental filtering and high spatial heterogeneity. In contrast, MTW and OWSW showed lower and more stable βNTI values, consistent with assembly processes driven by neutral factors (Figure 5c). Mantel correlations (Figure 5d) revealed that community composition in MTW was significantly associated with suspended solids (r = 0.584, p = 0.028) and total phosphorus (r = 0.404, p = 0.045), implying that localized nutrient enrichment plays an important role in structuring communities within this transitional water mass. In the other two water masses (EHSW and OWSW), correlations with the measured environmental variables were weaker. In EHSW, strong hydrological variability and estuarine mixing may decouple community structure from static local environmental gradients. In contrast, in OWSW, characterized by relatively stable hydrographic conditions, community structure appeared consistent with patterns expected under dispersal and ecological drift, consistent with the lower βNTI values observed.

3.4. The Underlying Influence of Seasonal Changes in Planktonic Basal Resources on Fish Community Structure

Prey community composition exhibited pronounced seasonal variation between April and November (Figure 6a). In April, the assemblage was strongly dominated by Arthropoda (57.45%), followed by Bacillariophyta (12.47%) and Chlorophyta (10.86%), whereas all other phyla individually contributed less than 5% (p < 0.05). This structure indicates a prey community composition dominated by zooplankton and phototrophic algae. In contrast, November communities displayed a markedly different composition, with Chlorophyta (34.04%) becoming the most abundant group, accompanied by substantial increases in Cnidaria (18.24%) and Pseudomonadota (7.56%), while Arthropoda declined sharply to 18.28% (p < 0.05). Bacillariophyta also decreased from 12.47% in April to 8.49% in November (p > 0.05). Minor taxa, including Rhodophyta, Chordata, Haptophyta, Oomycota, and Annelida, showed relatively low and stable contributions across seasons, each accounting for approximately 1–4% of total ASVs (p > 0.05).

Consistent with these phylum-level patterns, genus-level analysis further highlighted seasonal restructuring of representative prey taxa (Figure 6b). Zooplankton- and phytoplankton-associated genera, including Bathycoccus, Paracalanus, and Euphausia, were significantly more abundant in April than in November (Wilcoxon rank-sum test, p < 0.01). In addition, Parempheriella, a representative planktonic ciliate taxon, also showed significantly higher relative abundance in April than in November (Wilcoxon rank-sum test, p < 0.01). In contrast, Micromonas and CandidatusPelagibacter exhibited markedly higher proportions in November. Together, the phylum- and genus-level results demonstrate a coherent seasonal shift in prey community composition, from a spring assemblage dominated by zooplankton and relatively larger phytoplankton taxa to an autumn assemblage increasingly characterized by small-sized phototrophic algae and microbial-associated groups.

Corresponding seasonal contrasts were also evident in trophic associations between fish orders and prey phyla (Figure 6c,d). In April, although several prey–fish correlations were nominally significant based on uncorrected p-values (p < 0.05), none remained significant after controlling for multiple testing using the Benjamini–Hochberg FDR procedure (q > 0.05). Nevertheless, consistent prey-specific trends were apparent, particularly involving benthic-associated fish orders. Chlorophyta showed a positive association with Pleuronectiformes (r = 0.49, p = 0.007, q = 0.461) but a negative association with Gobiiformes (r = −0.46, p = 0.012, q = 0.465). Similarly, Rhodophyta and Haptophyta both exhibited negative correlations with Gobiiformes, with a stronger trend for Haptophyta (r = −0.55, p = 0.002, q = 0.288) and a weaker trend for Rhodophyta (r = −0.41, p = 0.027, q = 0.805). Beyond these relationships, several prey phyla displayed moderate but non-significant correlations with multiple fish orders. For example, Arthropoda tended to correlate positively with demersal and mid-water fish orders such as Gobiiformes and Centrarchiformes (r ≈ 0.20–0.33), while Bacillariophyta showed generally positive associations with Perciformes and Scombriformes (r ≈ 0.23–0.25), suggesting broader but relatively diffuse association patterns during spring. In November, prey-fish associations appeared more concentrated in a smaller number of prey–fish pairs, although none remained significant after FDR correction (q > 0.05). Chlorophyta showed a positive association with Pempheriformes (r = 0.51, p = 0.015, q = 0.388), while Arthropoda was positively associated with Pleuronectiformes (r = 0.56, p = 0.007, q = 0.271) and negatively correlated with Myctophiformes (r = −0.50, p = 0.018, q = 0.388). In addition, Cnidaria displayed the strongest positive association with Myctophiformes (r = 0.69, p < 0.001, q = 0.060) and a positive association with Tetraodontiformes (r = 0.56, p = 0.007, q = 0.271). Compared with April, these results suggest a shift toward more concentrated prey–fish association patterns under the autumn prey regime, but the inferred relationships should be interpreted cautiously given the lack of FDR-supported significance. However, given that prey communities were inferred from COI metabarcoding and summarized at relatively high taxonomic levels, these associations should be interpreted as statistical co-occurrence relationships rather than direct evidence of trophic coupling or energy transfer.

4. Discussion

Based on eDNA sequencing and co-occurrence network analysis, our study systematically revealed the multidimensional heterogeneity of fish communities across vertical, hydrographic, and seasonal environmental gradients. We further emphasized a multi-level pattern-based framework integrating environmental forcing, biological responses, and interaction feedbacks. Compared with traditional studies that rely on capture-based species inventories or linear diversity-environment correlations, this integrative approach provides a more comprehensive lens for decoding marine community complexity.

The vertical differentiation of fish communities in the Yangtze River Estuary and Adjacent East China Sea (YRE-ECS) arises from the combined effects of environmental stratification and biotic interactions operating across the water column. In this study, temperature, salinity, dissolved oxygen, and nutrient concentrations formed a tightly coupled environmental axis, establishing a pronounced stratified water column that delineates distinct surface and bottom habitats (Figure 5b). Consistent with this physical separation, fish assemblages exhibited clear vertical contrasts in taxonomic dominance (Figure 3a), with surface waters characterized by higher contributions of pelagic-associated orders such as Perciformes and Clupeiformes. In contrast, bottom waters were dominated by benthic-associated Gobiiformes and Pleuronectiformes. Correspondingly, surface communities exhibited higher and more dispersed absolute βNTI values (Figure 5a), indicating more substantial deviations from neutral expectations, consistent with a greater influence of deterministic assembly tendencies, likely in combination with stochastic influences under conditions of elevated temperature and intensified hydrodynamic disturbance [29,30]. Specifically, surface-water βNTI values were more dispersed and tended to be positive, with a higher frequency of βNTI > +2, suggesting stronger signals consistent with variable selection relative to null expectations, rather than uniform community turnover. This pattern is consistent with spatially heterogeneous surface conditions shaped by mixing fronts and dynamic hydrography, which can create strong environmental contrasts and patchy resources over short distances [31]. Across both surface and bottom waters, a substantial proportion of βNTI values fell within the neutral range (−2 < βNTI < +2), indicating that community turnover was dominated by weak phylogenetic deviations from null expectations at the spatial scale considered. In contrast, bottom communities showed βNTI values largely clustered within the neutral range, suggesting weaker detectable phylogenetic selection at the community level, environmental filtering associated with colder, more stable, and relatively oxygen-depleted conditions may still contribute to shaping species composition in bottom waters [32,33]. Importantly, βNTI values within the neutral range do not preclude deterministic assembly, but rather indicate that strong selection signals are not detectable relative to stochasticity under the current set of measured variables; unmeasured factors such as sediment characteristics, benthic habitat structure, or prey availability may still be associated with additional environmental filtering effects on bottom communities, while methodological limitations may constrain the detectability of such effects [34,35]. Beyond compositional differences, patterns of species co-occurrence were also structured along the vertical gradient. Co-occurrence network analysis (Figure 4) further revealed that several demersal species belonging to Gobiiformes, including Amblychaeturichthys hexanema and Parachaeturichthys polynema, occupied highly connected positions within network modules. This suggests their central positions in occupying highly connected positions within co-occurrence patterns, potentially reflecting shared habitat use or common responses to environmental gradients. While eDNA-based co-occurrence networks provide valuable insights into species co-occurrence, they primarily reflect shared occurrence rather than direct ecological interactions or trophic relationships. The nature of eDNA sequencing, including DNA degradation, shedding rates, and PCR amplification biases, can introduce methodological biases, affecting the relative abundance and accuracy of ecological interactions. Together, these results suggest that vertical stratification not only filters species according to habitat preferences but also reorganizes interaction networks by promoting distinct sets of dominant orders and key species across layers, a pattern consistent with vertical differentiation reported in other stratified coastal and shelf ecosystems [36,37].

Horizontally, inter-water-mass variation revealed that geographic gradients play an important role in shaping fish community differentiation (Figure 2b). In the MTW zone, fish community structure showed significant positive correlations with suspended solids and total phosphorus (Figure 5d), suggesting that nutrient enrichment and hydrodynamic disturbance in this estuarine-offshore transition zone enhance habitat heterogeneity and promote species coexistence [38,39]. This pattern is consistent with previous studies highlighting ecotonal characteristics of estuarine-coastal transition zones, where sharp environmental gradients may enhance habitat heterogeneity and species turnover [40]. OWSW communities exhibited lower and more stable βNTI values that were largely clustered within the neutral range (−2 < βNTI < +2) (Figure 5c), together with weak correlations with measured environmental variables (Figure 5d), indicating that community structure is maintained mainly through dispersal stochasticity and ecological drift under relatively stable hydrographic conditions. Similar assembly mechanisms have been reported in the North Sea [41] and the southwestern Atlantic [42], where low environmental gradients favor community stability driven by neutral stochasticity and geographic isolation. By comparison, EHSW displayed higher βNTI values and stronger signatures of environmental filtering (Figure 5c), reflecting its location at the freshwater-marine interface where steep nutrient gradients and frequent disturbances generate high spatial heterogeneity, as also evidenced by greater compositional variability among stations (Figure 2b). Importantly, across all water masses, the predominance of βNTI values within the neutral range indicates that differences among water masses primarily reflect relative shifts in the strength and direction of phylogenetic deviation, rather than discrete transitions between deterministic and stochastic assembly regimes. Collectively, these results indicate that horizontal structuring of fish communities in the Yangtze River Estuary reflects a balance between environmentally associated assembly tendencies and stochastic dispersal-related processes operating across contrasting hydrographic regimes. Consistent assembly mechanisms have been reported in the adjacent East China Sea [43,44], suggesting that the Yangtze River Estuary and Adjacent East China Sea functions as an important estuarine-coastal transition zone.

Specifically, prey communities transition from a spring assemblage enriched in zooplankton-related and phytoplankton-associated taxa to an autumn assemblage increasingly dominated by small-sized phototrophic algae and microbial components, accompanied by changes in the structure of prey-fish associations. From a trophic perspective, spring communities were supported by a diverse and zooplankton-rich prey base dominated by Arthropoda, together with substantial contributions from Bacillariophyta and Chlorophyta (Figure 6a). This pattern was further supported by the prevalence of zooplankton- and phytoplankton-associated genera such as Paracalanus and Euphausia (Figure 6b). Within this context, the increased spring abundance of planktonic ciliates may reflect enhanced trophic coupling by mediating the transfer of carbon and energy from pico- and nano-sized primary producers to higher trophic levels. This interpretation is consistent with a conceptual microbial loop framework, in which dissolved organic carbon is incorporated into bacterial biomass and subsequently transferred to higher trophic levels via protistan consumers such as ciliates [45]. Planktonic ciliates are recognized as important trophic intermediates that consume pico- and nanoplankton and serve as prey for larger zooplankton and fish, consistent with previous observations of their role in coastal food webs [46]. Such prey conditions provide multiple parallel energy pathways and promote relatively diffuse fish-prey associations (Figure 6b), which are characteristic of periods with elevated primary productivity and enhanced trophic transfer efficiency [44,45,47,48]. However, direct empirical evidence for microbial loop–mediated coupling in our system is still lacking, and future studies integrating measurements of pico- and nano-phytoplankton, bacterial production, and ciliate grazing rates are needed to validate this mechanism. Under these conditions, higher fish diversity and broader species overlap were observed in spring assemblages (Figure 2c and Figure 3c), indicating increased niche availability and functional redundancy. In addition to prey availability, seasonal shifts in species co-occurrence structure were observed, indicating changes in structuring seasonal dynamics. Co-occurrence network analysis revealed that several mid- to high-trophic species occupied central positions (Figure 4), including demersal Gobiiformes such as Amblychaeturichthys hexanema and Parachaeturichthys polynema (Gobiiformes), as well as more mobile taxa such as Scomber japonicus (Scombriformes) and Nibea albiflora (Perciformes). These species are predominantly carnivorous and occupy intermediate to high trophic levels, enabling them to be associated with a wide range of prey groups and co-occurring species [49,50,51]. Moreover, A. hexanema and P. polynema are demersal or burrowing species with strong sediment associations, whereas S. japonicus and N. albiflora exhibit high mobility and flexible habitat use across seasons, placing them in central positions within the co-occurrence network across habitats, suggesting broad habitat associations rather than demonstrated functional connectivity between benthic and pelagic compartments [52]. During spring, the influence of these key species is likely buffered by high prey diversity and functional redundancy, reducing competitive constraints and stabilizing community organization.

In contrast, autumn communities exhibited a pronounced restructuring of prey composition toward dominance by Chlorophyta and increased contributions of Cnidaria and microbial-associated taxa (Figure 6a), accompanied by enrichment of small-sized phototrophic and microbial genera such as Micromonas and CandidatusPelagibacter (Figure 6b). This shift is consistent with a contraction of effective energy pathways and reduced trophic breadth, consistent with the fewer and more concentrated fish-prey correlations observed in November (Figure 6c). This more simplified co-occurrence structure may reflect a concentration of statistical associations, rather than direct evidence of reduced energy flow or trophic resilience. Similar patterns have been reported in the South China Seas [53] and NE Atlantic Ocean [54], where warming and phenological shifts have led to shorter periods of primary productivity and reduced fishery yields [55]. These results indicate that seasonal fish community dynamics emerge from the interaction between prey-driven bottom–up forcing and the concentration of interspecific interactions around ecologically important taxa. Spring conditions promote trophic diversification and interaction buffering, whereas autumn conditions favor trophic compression and heightened dependence on structurally important species, potentially reducing functional redundancy and increasing ecological sensitivity [56,57].

These findings have direct implications for fisheries management and conservation in the Yangtze River Estuary and Adjacent East China Sea, particularly during the early implementation of the Yangtze River Basin Ten-Year Fishing Ban. The pronounced vertical and horizontal structuring of fish communities indicates that effective protection strategies should explicitly account for both water-column stratification and hydrographic heterogeneity. Notably, consistent with observations from fish communities surrounding oyster aquaculture areas [58], seasonal variation emerged as an important structuring factor of ecological dynamics in the Yangtze River Estuary. Against the background of the early implementation of the Yangtze River Basin Ten-Year Fishing Ban, seasonal signals may become more apparent as fishing disturbance is reduced, allowing natural environmental variability to play a clearer role in shaping community organization. Spring represents a critical ecological window, characterized by higher biodiversity and stronger trophic connectivity across surface and bottom layers as well as within estuarine and mixed transitional waters, underscoring the need to prioritize the protection of spawning and nursery habitats during this period. In autumn, the more simplified and vulnerable community structure observed, particularly in offshore and bottom habitats, suggests that more conservative harvesting practices are necessary to avoid disrupting recovering fish populations. Integrating long-term monitoring of planktonic basal resources with fish community indicators may further enhance ecosystem-based and adaptive management. These recommendations align with emerging conservation frameworks that emphasize functional diversity and network resilience [59,60].

5. Conclusions

This study provides a comprehensive assessment of spatial and seasonal variability in fish community structure across the East China Sea, using eDNA-based diversity profiling, co-occurrence networks, and environmental association analyses. Our results show that fish community organization varies systematically with vertical stratification, hydrographic gradients, and seasonal shifts in planktonic basal resources, reflecting the combined influence of environmental filtering, dispersal-related processes, and changes in species association patterns. Vertically, water-column stratification was associated with contrasting community assembly tendencies between surface and bottom layers. Bottom communities showed patterns more consistent with stochastic influences under relatively stable conditions, whereas surface communities exhibited stronger deviations from neutral expectations. Horizontally, estuarine-to-offshore transitions foster ecological heterogeneity, especially in the MTW zone, which exhibited higher compositional heterogeneity under enhanced nutrient availability and hydrodynamic disturbance. Seasonally, spring communities are characterized by high diversity and more complex association patterns, coinciding with enriched and diverse phytoplankton assemblages. In autumn, co-occurrence networks appeared more simplified and concentrated around fewer basal taxa, suggesting a more concentrated co-occurrence structure and potentially reduced functional buffering capacity, rather than direct evidence of altered trophic redundancy. These findings highlight the sensitivity of fish community dynamics to environmental gradients and seasonal fluctuations, particularly through associations consistent with potential bottom–up influences. Given the observed seasonal shift toward fewer significant prey-fish associations and a more concentrated interaction structure in autumn, maintaining the diversity and stability of basal producers may be critical for sustaining marine ecosystem resilience. From a conservation perspective, our study supports prioritizing the protection of spring spawning grounds and transitional water masses such as the MTW zone. It advocates adaptive fisheries management in autumn to avoid overexploitation during ecologically vulnerable periods. Incorporating long-term monitoring of plankton-fish coupling into ecosystem-based management may enhance predictive capacity and support sustainable stewardship in rapidly changing coastal regions.