Spatiotemporal heterogeneity in diazotrophic communities reveals novel niche zonation on the continental shelf of the East China Sea

The East China Sea (ECS) continental shelf is a hotspot for studying nitrogen fixation in the marginal seas of the western Pacific, where this microbially mediated process is profoundly influenced by both the coastal and oceanic current systems. Yet, how physical forcing controls the biogeography of diazotrophs and regional nitrogen budget on the ECS shelf remains poorly characterized. Here, we carried out a cross-season survey and demonstrated dynamics in diazotrophic communities that is tightly linked to distinct water masses on the ECS shelf. An overall spatial heterogeneity among some of the major diazotrophic phylotypes was unveiled, with the filamentous cyanobacteria Trichodesmium, diatom-diazotroph symbioses (Het-1 and Het-2), the unicellular cyanobacterial diazotrophs (UCYN-B) and Haptophyta-associated nitroplasts (UCYN-A) dominating the upper 30 m of the warm, nitrogen-limited offshore region intruded by the Kuroshio and Taiwan Strait water, whereas diatom-associated putative nitroplasts (UCYN-C and γ-24774A11) were abundant both at the surface and 50 m depth. The nitrogen fixation rates were generally higher in autumn than in spring, particularly in Kuroshio-affected waters dominated by Trichodesmium, Het-1 and Het-2, suggesting a seasonal variability primarily regulated by hydrographic conditions (mainly temperature and salinity) associated with water mass movement. Modeling the distribution of diazotrophs in the water masses identified three taxon-specific niches occupied by eight distinct diazotrophic phylotypes. Taken together, our analyses provide mechanistic insights into the role of dominant forms of physical forcing in driving the spatiotemporal variability in diazotrophic distribution and activity on the ECS shelf, which is of important reference in assessing diazotrophs adaptation in a changing marine ecosystem.

1 Introduction

Nitrogen is an essential nutrient for sustaining life on Earth. Although nitrogen is abundant in the atmosphere as dinitrogen gas (N₂), it is largely unusable in this form to most organisms on the planet. For nitrogen to be able to make biomolecules such as proteins, nucleic acids, and chlorophyll, it must first be converted into biologically available ammonia (NH₃) through N₂ fixation, a process mediated by a select group of prokaryotes termed diazotrophs (Zehr and Capone, 2020). This biologically fixed nitrogen fuels primary production and carbon export in the ocean, thus ultimately affecting global carbon and nitrogen cycles (Falkowski, 1997; Karl et al., 1997; Wang et al., 2019). While much research has been devoted to N₂ fixation in the open ocean (e.g. Chen et al., 2019; Shao et al., 2023; Wen et al., 2022), accumulating evidence suggests the potential of coastal aquatic systems as a significant yet underappreciated source of global N budget (Fulweiler et al., 2025; Tang et al., 2019). As the vital ecotones bridging terrestrial and marine ecosystems, coastal waters exhibit pronounced environmental gradients shaped by land-ocean interactions. These complex conditions challenge our understanding of the physicochemical processes involved in regulating the distribution and activity of marine diazotrophs (Fontela et al., 2025; Jiang et al., 2023a, b; Tang et al., 2019). Delineating the dynamics of diazotrophic biogeography in the nearshore environment is, therefore, key to modeling regional nitrogen cycle and resolving disparate accounts of fixed nitrogen budget globally (Fulweiler et al., 2025; Tang et al., 2019; Tang and Cassar, 2019; Zehr and Capone, 2020).

The East China Sea (ECS) continental shelf is a temperate marginal sea in the western North Pacific with complex hydrodynamic conditions governed by the East Asian monsoon as well as the intrusion of Kuroshio, a vigorous western boundary current as it enters the ECS shelf off the northeast Taiwan (Fig. 1) (Cui et al., 2021; Qiu and Imasato, 1990). The intruded Kuroshio is typically divided into two layers, the Kuroshio surface water (KSW, upper 120 m) and the Kuroshio subsurface water (KSSW, 120–250 m) (Chen et al., 1995). At approximately $\mathrm{27.5}\mathit{°}\mathrm{N}/\mathrm{122}\mathit{°}\mathrm{E}$, the KSSW bifurcates into the offshore and nearshore Kuroshio branch currents, with the latter transporting rich nutrients into the shelf region, particularly in summer (Yang et al., 2012, 2018). Meanwhile, Taiwan warm current (TWC), a mixture of waters from the Taiwan Strait and the Kuroshio, flows northward on the ECS shelf and occasionally reaches the Changjiang Estuary (Zhou et al., 2015). It has been shown that the southwest monsoon can enhance the transport of Kuroshio and TWC during warm seasons (Yang et al., 2018; Zhou et al., 2015), whereas the inshore ECS region is strongly influenced by the Changjiang diluted water (CDW) and the coastal water (CW), both of which are notably expanded seaward under the influence of the northeast monsoon (Yue et al., 2021). Thus, the dynamic interplay between these water masses creates steep biogeochemical gradients across the ECS shelf (Chen, 2009; Zhang et al., 2007). The upwelling of phosphorus-rich KSSW, coupled with terrestrial inputs from CDW and CW, may fuel phytoplankton growth in the coastal regions (e.g. Gao et al., 2025; Sun et al., 2025). In contrast, the warm, saline, and N-limited KSW and TWC establish favorable habitats for diazotroph communities and N₂ fixation (Jiang et al., 2018, 2019, 2023a, b). Previous studies have demonstrated that changes in hydrographic conditions can control nutrient fluxes and planktonic community shifts throughout the ECS ecosystem (Jiang et al., 2018, 2019, 2023a, b; Sun et al., 2025), however, the detailed patterns of diazotrophs distribution on the ECS shelf in relation to different water masses remain poorly characterized.

Figure 1 Sampling on the East China Sea (ECS) shelf during the 2023 autumn and 2024 spring cruises. A total of 42 stations along 5 transects were selected for the collection of biological samples and environmental parameters. The station numbers are positioned adjacent to corresponding points. Transects A (yellow line, stations 22–30) and B (yellow line, stations 36–42) were chosen to investigate variations in biological and environmental factors along the vertical gradient extending from inshore to offshore. Major circulations are indicated, including the Changjiang diluted water (CDW), Coastal water (CW), Taiwan warm current (TWC), Kuroshio current, Kuroshio surface water (KSW) and Kuroshio subsurface water (KSSW) (Yang et al., 2012, 2018). Arrow sizes denote specific discharge rates (Liu et al., 2021). Land topography and ocean bathymetry data were obtained from the General Bathymetric Chart of the Oceans (GEBCO, https://www.gebco.net/, last access: 24 January 2025).

Recent studies have elucidated Kuroshio intrusion as a key driver in transporting filamentous cyanobacterial diazotrophs, including the free-living Trichodesmium and diatom-associated symbionts (Richelia and Calothrix) (Jiang et al., 2018, 2019, 2023a, b), to the Changjiang Estuary during warm seasons, significantly elevating the rates of N₂ fixation in the Kuroshio waters expanding toward the shelf region compared to the CDW and CW (Jiang et al., 2023a, b). However, the influence of water mass-driven hydrographic processes on other globally distributed diazotrophs on the ECS shelf, such as the unicellular cyanobacterial diazotroph (e.g. UCYN-B) and Haptophyta-associated nitroplasts (early-stage N₂-fixing organelles) (Coale et al., 2024; Cornejo-Castillo et al., 2024), have not been adequately examined. Given the dominance of these diazotrophs in Kuroshio (Chen et al., 2014; Cheung et al., 2017, 2019; Shiozaki et al., 2018; Wu et al., 2018) and a currently disproportionate research focus on the filamentous Trichodesmium populations in individual seasons (Jiang et al., 2017, 2019, 2023a, b; Yue et al., 2021), systematic investigation into the compositional dynamics of diazotrophs across the ECS shelf is urgently needed. Furthermore, water mass movements on the ECS shelf are seasonally modulated by the East Asian monsoon (Yin et al., 2018), and the seasonal variability in current patterns is projected to be strengthened in climate models (Vélez-Belchí et al., 2013; Yang et al., 2024), thus necessitating a comprehensive understanding of how increasingly intensified water mass movements may restructure the distribution of distinct diazotrophic phylotypes on the ECS shelf (Tang and Cassar, 2019).

Realized niche modeling provides a trait-based framework for analyzing biogeographic patterns by examining species fitness and their ecological responses to environmental changes (Edwards et al., 2013; Irwin et al., 2012; Litchman et al., 2012), which has been integrated into biogeochemical models to predict species distribution under future climate change scenarios (Irwin et al., 2012, 2015). The Maximum Entropy model has proven particularly effective for characterizing species distribution in realized niches, even from sparse field observations (Irwin et al., 2012; Phillips et al., 2006). This model is capable of identifying multidimensional niche spaces along key environmental gradients (e.g. temperature and nutrients) through regularization that maximizes environmental dependency while minimizing observational bias (Irwin et al., 2012; Phillips et al., 2006). A major advantage of this method is its ability to operate without absence data or uniform sampling effort in space or time, making it suitable for oceanographic studies that rely on opportunistic sampling (Brun et al., 2015; Irwin et al., 2012). This approach has been successfully applied across diverse oceanic regimes, from the western Pacific marginal seas (Xiao et al., 2018; Zhong et al., 2020) and the North Atlantic (Irwin et al., 2012, 2015) to the global open ocean (Brun et al., 2015), consistently revealing that functionally analogous phytoplankton taxa (e.g. coccolithophores) occupy more similar niches than taxonomically distinct groups (e.g. diatoms). Despite the progress in modeling phytoplankton ecology in general, our understanding of the realized niches of diazotrophic communities remains fragmented, with only limited documentation on the modeling of realized niches of Trichodesmium and Richelia in the open ocean (Brun et al., 2015). However, these niche traits have yet to be extended to highly variable temperate marginal seas, thereby constraining our ability to predict N₂ fixation across the global ocean under changing climatic conditions.

To illuminate the patterns and drivers of N₂ fixation on the ECS shelf, we conducted a cross-season survey analyzing the distribution, abundance and activity of the major diazotrophic phylotypes known to date, alongside hydrographic analysis and integrated niche modeling approaches. Our analyses reveal spatiotemporal heterogeneity in diazotrophic distribution and N₂ fixation rate, which is highly correlated with distinct water masses. Overall, the realized niche partitioning presented here emphasizes the role of physical forcing (e.g. Kuroshio intrusion) in shaping the diversity and biogeography of diazotrophs in a marginal sea with complex land-ocean interactions.

2 Materials and methods

2.1 Field sampling and environmental parameter collection

We conducted a cross-season survey at 42 stations on the ECS shelf aboard the RV Xiang Yang Hong 18 during 2023 autumn (13–30 October) and 2024 spring (9–24 April) (Fig. 1). At each station, temperature and salinity were acquired from a Seabird 911plus conductivity–temperature–depth sampler (Sea-Bird Electronics, USA) and binned over 1 m depth intervals. For biological and environmental sample collection, seawater was collected at 3 or 7 discrete depths throughout the upper 100 m using Niskin bottles. Approximately 2 L of seawater was prefiltered through a 200 µm pore-size nylon mesh to remove large zooplankton, followed by filtration onto a 0.2 µm pore-size PES membrane filter (Pall, USA). The membranes were flash-frozen in liquid N₂ for subsequent DNA extraction. Seawater samples for nutrient analyses were filtered through 0.2 µm pore-size polycarbonate membranes (Millipore, USA) and then dispensed into 100 mL acid-washed polyethylene bottles and stored at −20 °C. Nitrate plus nitrite (hereafter abbreviated as NO_x), soluble reactive phosphorus (SRP) and dissolved silicon (DSi) concentrations were analyzed using a Technicon AA3 AutoAnalyzer (Bran Luebbe, Germany). The N:P ratios were calculated based on the molar ratio of NO_x and SRP. Seawater samples (prefiltered with 200 µm pore size mesh) for RNA analysis were collected at 7 stations in the upper 50 m using acid-rinsed Nalgene polycarbonate bottles during the day (13:00–15:30 LT) or at night (21:40–03:40 LT) to mirror diel variations in nitrogenase reductase gene (nifH) expression (Church et al., 2005b; Moisander et al., 2014). Each sample was taken by filtering 4.5 L of seawater onto a 0.2 µm pore-size PES membrane (Pall, USA), which was soaked in 1 mL of RNAlater (ThermoFisher, USA) and then instantly refrigerated in liquid N₂ for later extraction. Surface seawater was also obtained from 24 stations of the study area to measure the rates of N₂ fixation, with details described below. Additionally, two transects (A and B) spanning from nearshore to offshore waters were chosen to analyze vertical gradients of biological and environmental factors (Fig. 1).

2.2 DNA and RNA extraction

DNA was extracted using the cetyltrimethylammonium bromide method as described by Zhang et al. (2024). The concentration and purity of DNA were determined using a NanoPhotometer N50 spectrophotometer (Implen, Germany). RNA was extracted with TRI Regent (Molecular Research Center, USA), and the extracted RNA was subjected to DNA digestion and purification using a Direct-zol RNA Miniprep kit (Zymo, USA). DNA contamination in purified RNA was checked by amplifying the bacterial 16S rRNA gene with universal primers 341F and 806R (Wasimuddin et al., 2020). The nifH mRNA was reverse transcribed to synthesize complementary DNA (cDNA) using a RevertAid RT kit (Thermo Fisher Scientific, USA) with the gene-specific reverse transcription (RT)-primers nifH2 and nifH3 (Moisander et al., 2014), following the manufacturer's protocol.

2.3 TaqMan qPCR assay of targeted nifH genes and transcripts

Recent studies have demonstrated UCYN-A2 as a nitroplast in the haptophyte Braarudosphaera bigelowii (Coale et al., 2024; Cornejo-Castillo et al., 2024). Similarly, UCYN-C and γ-24774A11 have been proposed as putative nitroplasts in diatoms (Schvarcz et al., 2022, 2024; Tschitschko et al., 2024). In contrast, other UCYN-A sublineages such as UCYN-A1 have not yet been confirmed to possess the same defining characteristics in their association with haptophytes (Coale et al., 2024; Kantor et al., 2024). Therefore, for clarity and consistency with prior literature, we classified UCYN-A2 as a haptophyte nitroplast, UCYN-C and γ-24774A11 as putative diatom nitroplasts, while designating UCYN-A1 as other UCYN-A sublineages. We then performed TaqMan quantitative PCR (qPCR) to quantify nifH gene and transcript abundances of nine diazotrophic phylotypes, including Trichodesmium, Richelia (Het-1 and Het-2 associated with diatoms in the genera Rhizosolenia and Hemiaulus, respectively), Calothrix (Het-3 associated with Chaetoceros), UCYN-B as well as the N₂-fixing haptophytes and diatoms mentioned above, using established primer and probe sets (Table S1 in the Supplement) (Church et al., 2005a, b; Foster et al., 2007; Langlois et al., 2008; Moisander et al., 2008, 2010; Thompson et al., 2014). The reference nifH sequences used for making standard curves were synthesized at Sangon Biotech (Shanghai, China) (Table S1). Duplicate qPCR assays were run for each DNA/cDNA sample and the standards on a CFX96 real-time system (Bio-Rad Laboratories, USA). Each reaction mixture (10 µL) contained 5 µL of Premix Ex Taq (Takara Bio, Japan), 0.4 µM each of the forward and reverse primers, 0.4 µM of TaqMan probe, 2 µL of template DNA or cDNA, and nuclease-free water (Thermo Fisher Scientific). The thermal cycling conditions were 50 °C for 2 min, 95 °C for 2 min, and 45 cycles of 95 °C for 15 s, followed by 60 °C (64 °C for UCYN-A2) for 1 min. A standard curve was generated by serial dilutions (between 10–10⁷ copies per well) of a known concentration of the nifH references, and the negative controls were also included in each run to check contamination. The limit of detection (LOD) for the qPCR reactions was 10 nifH gene copies per reaction, corresponding to approximately 56–250 copies per liter of seawater. Considering that this qPCR assay could amplify UCYN-A2 together with relatively high abundances of UCYN-A3 and UCYN-A4 (Cheung et al., 2017, 2019; Farnelid et al., 2016), we referred to the target collectively as UCYN-A2/A3/A4.

2.4 Measurement of nitrogen fixation rate

Nitrogen fixation rates (NFRs) were determined using the ¹⁵N₂ gas dissolution method (Mohr et al., 2010; Großkopf et al., 2012; Montoya et al., 1996). Briefly, 0.22 µm-filtered seawater was degassed using a Sterapore membrane unit (20M1500A) (Mitsubishi Rayon, Japan) and then stored in a Tedlar bag without headspaces (Shiozaki et al., 2015). After that, the ¹⁵N₂ gas stock (99 %) (Cambridge Isotope Laboratories, USA) was injected into the bag at a ratio of 10 mL of ¹⁵N₂ gas per 1 L of seawater. The bag was gently tapped until the gas was fully equilibrated. At designated stations, duplicate acid-cleaned 2.3 L polycarbonate bottles were filled with bubble-free surface seawater, each spiked with 200 mL of ¹⁵N₂-enriched seawater and then incubated on deck for 24 h. Incubated samples were filtered (vacuum pressure <200 mm Hg) onto pre-combusted (450 °C, 4 h) GF75 membranes (Advantec, Japan), which were immediately stored in liquid N₂ until onshore processing for particulate organic N (PON) and ¹⁵N abundance measurements. Unamended natural seawater samples for PON were also collected for comparison as a blank control. The PON and ¹⁵N abundances were analyzed using an elemental analyzer interfaced with an isotope ratio mass spectrometry (Thermo Fisher Scientific). The NFRs were determined following Montoya et al. (1996). The LOD for the NFRs was calculated according to White et al. (2020), defined as a minimum difference of 0.00146 atom % in ¹⁵N enrichment of particulate nitrogen between the initial and the final measurements. This LOD corresponds to a range of 0.11–0.76 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ across the stations (Table S2).

2.5 Statistical analyses

2.5.1 Water mass analyses

To identify the hydrographic dynamics induced by distinct water masses, an optimum multiparameter analysis (OMPA) (Tomczak and Large, 1989) was applied to calculate their contribution at each sampling depth. The method quantitatively evaluates the mixing ratios of source water types (SWTs) based on quasi-conservative properties (e.g. temperature and salinity) by solving a series of linear equations. Based on previous studies (Chen, 2009; Yang et al., 2012) and the temperature-salinity diagram (Fig. S1 in the Supplement), the ECS shelf is primarily influenced by the CDW/CW, TWC, KSW and KSSW, and we determined the properties (temperature and salinity) of CDW/CW by the observational data from stations SF1 ($\mathrm{31.16}\mathit{°}\mathrm{N}/\mathrm{122.56}\mathit{°}\mathrm{E}$) and SF2 ($\mathrm{30.49}\mathit{°}\mathrm{N}/\mathrm{122.67}\mathit{°}\mathrm{E}$) (Table S3). Given that the TWC represents a mixture of Kuroshio and Taiwan Strait water (TSW), we used TSW, rather than TWC, as an end-member to better evaluate the influence of Kuroshio intrusion on diazotrophs. Accordingly, the properties of TSW, KSW and KSSW were obtained from the World Ocean Atlas 2023 dataset (Fig. S2 and Table S3) (Locarnini et al., 2024; Reagan et al., 2024). The OMPA was carried out with Python package “pyompa” (Shrikumar et al., 2022). The OMPA results were validated by the model residuals, which lay within the acceptable thresholds for temperature (−0.80–0.86) and salinity (−4.51–0.07) as defined by Lawrence et al. (2026).

2.5.2 Significance test, correlation analysis and multivariate model

A Wilcoxon signed-rank test was employed to assess the significant differences (p<0.05) in biological and environmental variables across regions and sampling seasons because of the violation of assumptions related to normality and homogeneity of variance. Pearson correlation analysis was adopted to explore relationship among nifH gene and transcript abundances (log ₁₀ transformation), NFRs, environmental factors (z score scaling), and distinct water masses (centered log ratio transformation) obtained from OMPA. Furthermore, a multivariate generalized additive model (GAM) with Tweedie distribution was used to estimate the relative contribution of environmental variables to variations in diazotroph abundances and NFRs, with explained deviance hierarchically partitioned (Lai et al., 2024). Concurvity among variables was assessed, and all pairwise indices below 0.8 were kept to enhance model robustness. To identify key environmental drivers, we applied an automatic smoothness selection method, which effectively shrinks non-significant terms to zero degrees of freedom and thus enables automatic variable selection (Marra and Wood, 2011).

2.5.3 Univariate realized niche parametrization

Realized niches of diazotrophs in relation to environmental variables on the ECS shelf were determined in a combined framework implementing the univariate Maximum Entropy (MaxEnt) (Phillips et al., 2006) and GAM, as described by Irwin et al. (2012) and Xiao et al. (2018). This method was chosen over the standard multivariate generalized linear model (GLM) or GAM for several reasons to overcome the challenges in analyzing marine ecological data. While the multivariate models can capture species-environment relationships and nonlinear effects, they require reliable presence/absence records and treat non-detections as true absences, an assumption often invalid in undersampled, patchy marine systems where absences frequently indicate missed detection rather than actual exclusion. The MaxEnt component of our framework circumvents this issue by modeling presence data against background environmental conditions, making it more robust for sparse and patchy biological observations. Additionally, environmental variables in oceanic systems are often highly collinear (e.g. NO_x and SRP), which in a multivariate model can obscure individual predictor effects, destabilize parameter estimates, and produce ecologically ambiguous response curves. Although collinearity may be reduced by variable selection techniques (e.g. variance inflation factors), it comes at the cost of excluding some key drivers that affect species niche traits. In contrast, the univariate approach examines the response along each environmental gradient, providing clear and interpretable niche estimate for robust, cross-phylotype comparisons (Irwin et al., 2012). Finally, the combined framework used here can independently model two ecologically distinct niche components: the probability of occurrence (via MaxEnt) and the expected abundance when present (via GAM), thus avoiding the conflation and bias caused by the zero values in abundance data inherent to multivariate models. This framework has proven robust for sparse and heterogeneous biological data within environmental matrices, as demonstrated in prior studies of phytoplankton niches in the western Pacific marginal seas (Xiao et al., 2018; Zhong et al., 2020).

The MaxEnt-GAM framework is formulated as follows:

$$\begin{array}{}\text{(1)}& {\displaystyle }& {\displaystyle }f\left(x\right)=P\left(y=\mathrm{1}\mathrm{|}x\right)\times C\left(x\right)\text{(2)}& {\displaystyle }& {\displaystyle }P\left(y=\mathrm{1}\mathrm{|}x\right)=P\left(y=\mathrm{1}\right)\times g_{\mathrm{1}}\left(x\right)/g\left(x\right)\text{(3)}& {\displaystyle }& {\displaystyle }C\left(x\right)=\mathit{\alpha }+s\left(x\right)+\mathit{\epsilon }\end{array}$$

where f(x) denotes the diazotroph abundance determined by a specific environment, x. The conditional probability of detecting diazotrophic phylotype in the environment, P(y=1|x), is evaluated using Bayes' theorem of the MaxEnt. P(y=1) denotes the probability that diazotrophic phylotype would be found in a random sample. The probability distribution functions g(x) and g₁(x) are estimated from the environmental condition of all available background observations and from the conditions where the phylotype is present, respectively. C(x) represents the estimated abundance of the diazotrophic phylotype as derived from GAM, using the same observational data as g₁(x), but with abundance data (log ₁₀ transformation) instead of presence-only data. The term s(x) refers to a one-dimensional nonlinear function based on cubic regression splines, α is the grand mean, and ε denotes the error term.

Following Irwin et al. (2012), we characterized the univariate response curves, f(x), with two parameters: niche mean (μ) and breadth (σ), calculated as:

$$\begin{array}{}\text{(4)}& {\displaystyle }& {\displaystyle }\mathit{\mu }={\displaystyle \frac{\int xf\left(x\right)\mathrm{d}x}{\int f\left(x\right)\mathrm{d}x}}\text{(5)}& {\displaystyle }& {\displaystyle }{\mathit{\sigma }}^{\mathrm{2}}={\displaystyle \frac{\int {\left(x-\mathit{\mu }\right)}^{\mathrm{2}}f\left(x\right)\mathrm{d}x}{\int f\left(x\right)\mathrm{d}x}}\end{array}$$

where μ represents the abundance-weighted mean environmental condition (i.e. the realized niche occupied by the phylotype in the study area), while σ indicates the phylotype's tolerance range for a given variable.

We employed MaxEnt software (version 3.4.3) for model fitting and related statistics. To minimize model overfitting and complexity, threshold features were disabled (Irwin et al., 2012), while linear, quadratic, product and hinge features were enabled. All other parameters remained at default values. Model performance was assessed using the receiver operating characteristic (ROC) metric, which compares the true positive rate to the false positive rate across multiple probability thresholds applied to continuous probabilistic output from the MaxEnt. The area under the ROC curve (AUC) served as a summary metric, where 1 indicates perfect discrimination and 0.5 represents random prediction. For each diazotrophic phylotype, the univariate model with a mean cross-validation AUC significantly greater than 0.5 at the 95 % significance level (one-sample Student's t test, five replicates) were retained. Subsequently, 75 % of presence records were randomly selected with replacement for model training, and the remainder were used for testing. We performed 1000 bootstrap iterations for every phylotype-environmental variable combination, recording the logistic probability for each resampled set. Preliminary sensitivity analysis confirmed that 1000 iterations yielded stable niche estimates and confidence intervals (data not shown). The GAMs with default settings were similarly applied to the same bootstrapped samples using abundance data (log ₁₀ transformation).

Owing to limited data for some diazotrophic phylotypes, we integrated published data (Cheung et al., 2019; Sato et al., 2025b; Shiozaki et al., 2018) with our own field observations on the ECS shelf (this study plus unpublished surface data from October–November 2022 and April–May 2023). The total number of data points utilized in the MaxEnt-GAM framework for each diazotrophic phylotype was presented in Table S4. Principal component analysis (PCA) was utilized to discern the relationships between diazotrophic phylotypes and the realized niches (z score scaling). Ward's minimum variance clustering was conducted based on Euclidean distances calculated from standardized scores of the first two PCA components. The optimal cluster number was determined by analyzing within-cluster sum of squared errors through an elbow plot across candidate cluster numbers.

The Wilcoxon signed-rank test and one-sample Student's t test were performed with R (version 4.1.3, http://www.r-project.org, last access: 30 December 2025) package “stats”, Pearson correlation analysis with R package “psych” (Revelle, 2025), PCA and clustering analysis with R package “vegan” (Oksanen et al., 2024), GAM with R package “mgcv” (Wood, 2011), hierarchical variance partitioning for GAM with R package “gam.hp” (Lai et al., 2024), and pie chart with R package “scatterpie” (Yu, 2025).

3 Results

3.1 Environmental conditions

The physicochemical parameters on the ECS shelf showed apparent spatial variations during autumn and spring (Fig. 2). Sea surface temperature was significantly higher in autumn (23.1–27.6 °C) than that in spring (14.0–25.8 °C; p<0.05) (Fig. 2A and F), while no significant seasonal difference was detected in salinity, though autumn values had lower range (29.7–34.4 PSU) compared to spring (26.6–34.6 PSU) (Fig. 2B and G). Horizontally, the highest temperature (>25 °C) and salinity (>34 PSU) were observed on the outer ECS shelf, accompanied by high proportions (∼80 %) of KSW and TSW (Fig. S3B, C, F, and G), indicating strong intrusions from these water masses. In contrast, the lowest temperature and salinity values were found mainly in the coastal regions influenced by CDW/CW (>60 % in proportion; Fig. S3A and E). Particularly in spring, the seaward extension of CDW/CW (temperature <18 °C and salinity <31 PSU) was strongly enhanced under the East Asian monsoon (Fig. S3E). Vertically, the cold (<20 °C), saline (>34 PSU) waters dominated the bottom layer along the transects A and B during the sampling seasons (Fig. S4), indicating significant KSSW intrusions (>80 % in proportion) on the ECS shelf (Fig. S3D and H). The mixed-layer depth was high (>40 m) on the eastern and southeastern ECS shelf (Fig. S5), corresponding approximately to Kuroshio and TSW intrusion (Fig. S3B, C, F, and G). Surface NO_x concentration in autumn (0.44–17.71 µM) was in a narrower range than that in spring (0.16–19.11 µM), but the values were 6–12 times higher in nearshore waters (stations 1, 12, 21, 22, 31, and 36) relative to offshore regions (stations 5–11, 15–20, 25–30, 33–35, and 38–42) in both seasons (p<0.05; Fig. 2C and H). However, SRP concentration was found to be comparable between these regions during autumn (0.14–0.91 µM; Fig. 2D) and spring (0.14–0.84 µM; Fig. 2I), suggesting that CDW/CW and KSW may be important sources for SRP in nearshore and offshore environments, respectively. A trend of increasing SRP with depth was observed, with the concentrations ranging from ∼0.3 µM at the surface to 0.7–1.2 µM at the bottom (Fig. S3I–L). The N:P ratio exhibited a lower range in autumn (1.3–19.5) compared to spring (0.4–37.3), but exceeded the Redfield ratio (16) in coastal waters in both seasons (autumn: 18.2±1.2 at stations 1, 22, and 31; spring: 25.7±7.0 at stations 1, 4, 12, 13, 14, 21, and 31), mirroring the spatial distribution pattern of NO_x (Fig. 2E and J).

Figure 2 Variations in surface temperature (A, F), salinity (B, G), NO_x (C, H), SRP (D, I) and N:P ratio (E, J) during autumn (A–E) and spring (F–J) on the East China Sea shelf. The map was created using Ocean Data View 5.7.2 (Schlitzer, Reiner, Ocean Data View, https://odv.awi.de, last access: 24 January 2025).

3.2 Variations in the distribution of major planktonic diazotrophs

Significant spatial heterogeneity in nifH gene abundance between autumn and spring was observed in almost all the diazotrophic phylotypes examined except Het-3, which remained undetectable in qPCR assay (Fig. 3). In autumn, Trichodesmium emerged as the dominant diazotroph, reaching up to 6.6×10⁵ copies L⁻¹ at surface and 9.2×10⁹ copies m⁻² throughout the water column (Fig. 3A and C). Trichodesmium accounted for approximately 90 % of the nifH gene pool of the targeted diazotrophs on the ECS shelf, with both the surface and depth-integrated abundances increasing along the coast-to-offshore transects. Additionally, Hets (Het-1: 11 533±11 039 copies L⁻¹; Het-2: 5439±7317 copies L⁻¹), UCYN-B (9107±13 009 copies L⁻¹) and γ-24774A11 (1961±2316 copies L⁻¹) were found primarily on the southeastern ECS shelf, with 1–2 orders of magnitude lower in abundance than Trichodesmium. The abundances of the UCYN-A (UCYN-A1 and UCYN-A2/A3/A4) were mostly below 1.0×10³ copies L⁻¹ in the surveyed areas. In spring, most of the targeted diazotrophs cannot be detected, except at some stations on the southeastern ECS shelf, where Trichodesmium was the most abundant (up to 6.0×10⁵ copies L⁻¹ at surface and 1.1×10¹⁰ copies m⁻² throughout the water column), followed by UCYN-A1 (up to 1.7×10⁵ copies L⁻¹ and 8.4×10⁹ copies m⁻²) (Fig. 3B and D). The abundances of the other targeted diazotrophs (Hets, UCYN-A2/A3/A4, UCYN-B, UCYN-C and γ-24774A11) varied between 1.0×10³–3.0×10⁴ copies L⁻¹ on the ECS shelf.

Figure 3 Surface (A, B), and depth-integrated (C, D) abundances of the eight major diazotrophic phylotypes on the East China Sea shelf during autumn and spring as determined based on the quantification of the nifH gene (i.e. DNA-based) with TaqMan qPCR. The station numbers are positioned adjacent to corresponding pie charts and crossover points. Land topography was obtained from the General Bathymetric Chart of the Oceans (GEBCO, https://www.gebco.net/, last access: 24 January 2025).

The targeted diazotrophs also exhibited seasonal variations in vertical distribution (Figs. 4 and S6). The average nifH gene abundances of Trichodesmium and Hets in the upper 50 m were considerably higher in autumn than those in spring, but this seasonal difference diminished gradually with increasing depths (Fig. 4). Remarkably, UCYN-B displayed the seasonal abundance difference continuously throughout the water column (Fig. 4). In contrast, the UCYN-A and UCYN-C in the upper 50 m were markedly less abundant in autumn than that in spring, while γ-24774A11 showed no depth-dependent seasonal variation (Fig. 4). Profiling along the transects A (stations 22–30) and B (stations 36–42) further revealed that Trichodesmium occurred more widely in autumn (recorded at 57 stations) compared to spring (20 stations), but its abundance consistently peaked (>10⁵ copies L⁻¹) in the upper 30 m during both seasons (Fig. S6A–D). The Hets were mainly found in the upper 30 m, and most abundant (>10⁴ copies L⁻¹) on the southeastern ECS shelf during spring (Fig. S6E–L). Among them, Het-2 showed a broader distribution (detected at 47 stations in autumn and 38 stations in spring) than Het-1 (12 stations per season) and peaked in the subsurface (30 m depth) (Figs. 4 and S6E–L). The putative diatom nitroplasts (i.e. UCYN-C and γ-24774A11) were most abundant (>10³ copies L⁻¹) at depths of 0–50 m (Figs. 4 and S6U–Z), demonstrating moderately deeper distributions than Hets. UCYN-A2/A3/A4 was similar to Trichodesmium in vertical distribution but with lower abundances (∼10⁴ copies L⁻¹) in spring, while UCYN-A1 and UCYN-B reached their highest abundance (>10⁴ copies L⁻¹) at depths of 0–30 m (Fig. S6M–T). Overall, diazotroph abundances decreased markedly below the 50 m layer, except for UCYN-A1 and UCYN-B which exhibited a moderate decline, possibly due to low temperature caused by KSSW intrusion (Figs. S4 and S6).

Figure 4 Vertical distributions in the average abundances of the eight major diazotrophic phylotypes on the East China Sea shelf during autumn and spring as determined based on the TaqMan qPCR assay of the nifH gene (i.e., DNA-based). The grey solid and dashed lines represent ±1 standard deviation for autumn and spring, respectively.

3.3 Pattern of nifH transcription in the targeted diazotrophs

Transcription of the nifH gene in the targeted diazotrophs generally mirrored the observed gene abundance pattern but displayed greater seasonal variability across the stations (Fig. 5). In autumn, the nifH transcripts of Trichodesmium were detected as the most abundant among 6 of the 7 stations sampled, ranging from 6.8×10³ to 1.7×10⁵ nifH transcripts L⁻¹ at surface and from 5.1×10⁴ to 2.4× 10⁶ nifH transcripts m⁻² in the upper 50 m (Fig. 5A and C). The nifH transcripts of Het-2 and UCYN-B were similarly high on the southeastern ECS shelf (8.0×10³ nifH transcripts L⁻¹ at surface and 1.7×10⁵ nifH transcripts m⁻² in the upper 50 m water column). The nifH transcripts of Het-1 and γ-24774A11 were 1–2 orders of magnitude lower relative to the other targeted diazotrophs. Overall, Trichodesmium dominated the nifH transcript pool (>70 %) of the targeted diazotrophs, except at station 30 where Het-2 (40 %) and UCYN-B (40 %) exhibited relatively higher abundance. In spring, the nifH transcripts were only detected at 2 of the 7 stations sampled (Fig. 5B and D). At station 30, Trichodesmium, Het-1, UCYN-A1 and UCYN-B equally dominated the nifH transcript pools of the targeted diazotrophs at surface and in the upper 50 m water column. At station 40, however, Trichodesmium, Het-1 and UCYN-A1 transcripts were abundant at surface, whereas UCYN-B prevailed at depths of 10–50 m. Despite the station-to-station variations, some diazotrophs exhibited distinct diel variations in nifH transcription (Fig. S7). The nifH transcription levels of Trichodesmium, Het-1 and UCYN-A were significantly higher at daytime (Fig. S7A, B, D, and E), whereas Het-2 and UCYN-B mainly expressed the nifH gene at night (Fig. S7C and F). The γ-24774A11 displayed no detectable diel pattern in nifH transcription (Fig. S7G).

Figure 5 Surface (A, B), and depth-integrated (C, D) abundances of the seven major diazotrophic phylotypes on the East China Sea shelf during autumn and spring as determined based on TaqMan qPCR assay of the nifH transcript (i.e. RNA-based). Note that the transcript abundance of UCYN-B was derived from nighttime samples, whereas that of the other groups was from daytime samples. The depth-integrated transcript abundance was calculated based on the integration of the upper 50 m of the water columns. The station numbers are positioned adjacent to corresponding pie charts and crossover points. Land topography was obtained from the General Bathymetric Chart of the Oceans (GEBCO, https://www.gebco.net/, last access: 24 January 2025).

3.4 Seasonal variations in the rates of nitrogen fixation

The NFRs in surface waters exhibited distinct seasonal variations on the ECS shelf (Fig. 6 and Table S2). In autumn, the NFRs ranged from undetectable to 6.28 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ (1.35±1.45 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$, on average). In spring, the NFRs ranged from undetectable to 1.67 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ (0.74±0.42 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ on average) and were about half of the autumn rates on average. In both seasons, the rates were higher at Kuroshio-influenced stations (e.g. 20, 30, and 42), but lower in the CDW/CW-affected regions (>60 % of water mass; e.g. stations 1 and 31; Fig. S3A, C, E, and G).

Figure 6 N₂ fixation rates in surface waters at designated stations on the East China Sea shelf during autumn and spring as determined with in situ isotope tracing. The station numbers are positioned adjacent to corresponding bar charts and crossover points. Sampling excluded stations 7, 37, and 42 in autumn, and stations 8 and 36 in spring. Land topography was obtained from the General Bathymetric Chart of the Oceans (GEBCO, https://www.gebco.net/, last access: 24 January 2025).

3.5 Relative contribution of environmental variables to diazotroph abundances and N₂ fixation rates across distinct water masses

In both seasons, diazotroph abundances were positively correlated with temperature and salinity, and negatively with sampling depth, nutrient concentrations (NO_x, SRP and DSi), and N:P ratio, although UCYN-B deviated from these trends in autumn (Fig. 7A). Similarly, the NFRs were also positively correlated with temperature and salinity and negatively with nutrient levels (except for SRP in spring) and N:P ratio. Regarding water masses, diazotroph abundances and activity (i.e. NFR) were positively associated with KSW in both seasons and TSW in spring, but negatively associated with CDW/CW and KSSW in spring. Additionally, the NFRs correlated positively with diazotroph abundances but negatively with their nifH gene transcripts across both seasons (Table S5), indicating that nifH transcript levels are not a reliable proxy for N₂ fixation. Hierarchical variance partitioning in GAMs further indicated that the measured factors explained 0.6 %–97.5 % of the spatial variation in diazotroph abundances and NFRs (Fig. 7B). Temperature and salinity accounted for the majority (59.2 % on average) of the explained variation, followed by nutrients and N:P ratio (37.5 % on average). Collectively, these findings demonstrate that the distribution and activity of targeted diazotrophs were primarily constrained by KSW, characterized by high temperature and salinity, shallow water depth, and limited N.

Figure 7 Relative contribution of environmental variables to diazotroph abundances and N₂ fixation rates (NFRs) across distinct water masses on the East China Sea shelf during autumn and spring. (A) Correlation of diazotrophic phylotypes with environmental variables, water masses and NFRs, as indicated with Pearson coefficients. (B) Relative contribution of environmental variables to diazotroph abundances and NFRs, as indicated with the total deviance explained. Total, combined nifH gene abundances of the eight diazotrophs detected; CDW, Changjiang diluted water; CW, Coastal water; TSW, Taiwan Strait water; KSW, Kuroshio surface water; KSSW, Kuroshio subsurface water; *p<0.05; **p<0.01; ***$p<\mathrm{0.001}.$

3.6 Model performance and diazotrophic realized niches

Prior to the MaxEnt-GAM framework, model performance for the realized niches of each diazotrophic phylotype in relation to each individual environmental factor was assessed using AUC scores from MaxEnt. The AUC values ranged from 0.61 to 0.94 across diazotrophic phylotypes and environmental variables, with an average of 0.77 (Table S4), suggesting better predictive accuracy of the univariate models than the random predictions. Notably, UCYN-C achieved the highest model performance with an average AUC of 0.87 across all the variables, likely due to the relatively small dataset (N=30; Table S4). In contrast, Trichodesmium yielded the lowest average AUC (0.66), which may be attributed to its larger sample size (N=242).

We characterized the realized niche mean (μ) and breadth (σ) of diazotrophs on the ECS shelf with combined hydrographic, environmental and biological data, identifying three distinct clusters among the diazotrophs (Figs. 8 and S8). Trichodesmium exhibited low, broad temperature and salinity niches (μ_Temperature=23.5 °C, σ_Temperature=3.5 °C, μ_Salinity=33.2 PSU, and σ_Salinity=1.2 PSU) as well as high nutrient niches (${\mathit{\mu }}_{{\mathrm{NO}}_x}=\mathrm{6.93}$ µM, μ_SRP=0.61 µM, μ_DSi=11.12 µM, and ${\mathit{\mu }}_{\mathrm{N}:\mathrm{P}}=\mathrm{7.64}$). UCYN-A1, UCYN-A2/A3/A4 and UCYN-B generally cohabited in high-salinity (34.3 PSU on average) and moderate-temperature (24.1 °C on average) environments. UCYN-B, however, was slightly separated from UCYN-A due to its relatively high nutrient niches (${\mathit{\mu }}_{{\mathrm{NO}}_x}=\mathrm{2.91}$ µM and μ_SRP=0.48 µM; Fig. 8G). The diatom-diazotroph symbioses (i.e. Hets, UCYN-C and γ-24774A11) were associated with high-temperature (25.7 °C on average), low-NO_x (0.98 µM on average) and moderate-SRP conditions (0.39 µM on average), demonstrating substantial ecological overlap among them.

Figure 8 Scatter plots (A–F) and principal component analysis (G) depicting niche mean (μ) and breadth (σ) of each diazotrophic phylotype relative to environmental variables on the East China Sea shelf. Lines (A–F) across the data points indicate 95 % confidence intervals for each variable derived from 1000 bootstrap resampling. Clustering among the diazotrophs (G) is shown in grey lines. Elbow plot (H) shows the determination of the number of optimal clusters based on within-cluster sum of squares.

4 Discussion

4.1 Spatial distribution of key planktonic diazotrophs on the ECS shelf

Temperature and nutrient availability are widely recognized as key environmental factors governing the biogeographical distribution and physiological constraints of filamentous diazotrophs (Chen et al., 2019; Tang and Cassar, 2019; Tuo et al., 2021). As expected, our quantification of the nifH gene copies showed that Trichodesmium and Hets were more abundant on the warm, NO_x-limited outer ECS shelf (Fig. 3) frequently intruded by the Kuroshio and TSW (Fig. S3B, C, F, and G). In contrast, their abundances were considerably lower in the cold, NO_x-rich coastal waters influenced by the CDW/CW (Fig. S3A and E), showing positive correlation with temperatures but negative correlation with nutrient concentrations and N:P ratios (Fig. 7). The distribution of Trichodesmium and Hets in the CDW/CW of the ECS shelf has been demonstrated via microscopic observations in previous studies (Jiang et al., 2018, 2019, 2023a, b), however, this pattern is different from the observations of filamentous diazotrophs that are abundant (10⁴–10⁶ copies L⁻¹) in the Amazon (Foster et al., 2007) and Mekong River plumes (Bombar et al., 2011), where the concentration of NO_x is low but that of SRP, DSi and iron is high. The ratios of N:P (autumn: 18.2; spring: 25.7) and NO_x:DSi (0.8 in both seasons) in the CDW of the ECS shelf are much higher than those in the Amazon (N:P: 1.3; NO_x:DSi: 0.01) (Foster et al., 2007) and Mekong River plumes (N:P: 3.7; NO_x:DSi: 0.1) (Nguyen-Ngoc et al., 2023), implying that nutrient stoichiometry plays a key role in structuring diazotroph communities in these regimes (Gomes et al., 2018). Notably, while our modeling result indicated Trichodesmium can survive in relatively cool (23.5 °C), NO_x-rich environments (6.93 µM; Fig. 8), it was detected as the most abundant on the warm, NO_x-depleted outer ECS shelf (Fig. 3), a pattern that is consistent with the basal physiology of Trichodesmium (Brun et al., 2015; Jiang et al., 2025; Tang and Cassar, 2019).

Additionally, non-filamentous diazotrophs (UCYN-A, UCYN-B, UCYN-C and γ-24774A11) were abundant ($>\mathrm{5}\times {\mathrm{10}}^{\mathrm{4}}$ nifH gene copies L⁻¹) in Kuroshio-influenced waters, but not detectable in the CDW/CW-dominated regions (Figs. 3 and S3). Highly abundant UCYN-A1, UCYN-B and γ-24774A11 (>10⁴ nifH gene copies L⁻¹) were consistently found in Kuroshio, particularly in its upstream region (Chen et al., 2019; Wen et al., 2022) and surface waters spanning from northeastern Taiwan to south coast of mainland Japan (Cheung et al., 2019; Shiozaki et al., 2018). It is, therefore, reasonable to attribute the prevalence of these diazotrophs on the outer ECS shelf to the lateral transport by Kuroshio mainstream, as proved by the positive correlation between their abundance and Kuroshio intrusion (Fig. 7). Overall, filamentous diazotrophs were the primary N₂ fixers on the ECS shelf (Jiang et al., 2023a, b), as illustrated by the prevalence of their nifH gene (∼84 % of all the nifH copies detected; Fig. 3) that was positively correlated with NFRs (Fig. 7A) across most of the surveyed areas. On the other hand, UCYN-A1, UCYN-B and γ-24774A11 were only occasionally found in the Kuroshio-influenced regions (Figs. 3 and S3), and the correlation of their nifH gene abundances with NFRs (Fig. 7A) suggests their potential contribution to regional N₂ fixation (Chen et al., 2014; Wu et al., 2018). It should be noted that the abundance of colonial diazotrophs like Trichodesmium and Hets may be potentially underestimated due to the use of <200 µm size fraction in our study (e.g. Jiang et al., 2023a).

Despite a gradual decrease in abundance with increasing depths, the diazotrophs exhibited taxon-specific vertical distributions on the ECS shelf (Figs. 4 and S6). Specifically, Trichodesmium showed maximum abundance in the upper 30 m (Fig. S6A–D), consistent with earlier reports of surface or subsurface maxima (10–30 m) on the ECS shelf (Jiang et al., 2018, 2019) and other regions (Carpenter et al., 2004; Chen et al., 2003; Lu et al., 2019). This depth distribution may be attributed to favorable nutrient conditions (e.g. SRP and iron) and light intensity (Turk-Kubo et al., 2018; Wen et al., 2022), which support Trichodesmium growth on the site of the ongoing surveys (Fig. 2D and I). Additionally, the absence of anticyclonic eddies during our surveys is likely to eliminate the downward entrainment of Trichodesmium in deeper water column (e.g. 50–60 m) (Fig. S9) (Jiang et al., 2018). The filamentous diazotrophs Het-1 and Het-2 showed distinct maximum depths (surface versus 30 m; Fig. 4), contrasting with previously documented subsurface maxima (10–40 m) for both phylotypes on the ECS shelf (Jiang et al., 2019), the tropical North Atlantic (Foster et al., 2007; Goebel et al., 2010) and the western tropical South Pacific (Stenegren et al., 2018). This may be due to different habitat preferences between the Het hosts, as the Het-1 host was found mainly under NO_x-depleted conditions but not in the deeper, NO_x-repleted water column (Tuo et al., 2014). Alternatively, the differences in thermal sensitivity between Het-1 and Het-2 may explain their vertical distribution, as Het-1 lives within a relatively narrower temperature range (20.4–27.5 °C) than Het-2 (17.6–27.6 °C; Fig. S4). In contrast, the other two diatom-diazotroph symbioses (UCYN-C and γ-24774A11) showed both surface and subsurface maxima (∼50 m; Fig. 4). The deep distributions may result from the low light tolerance of their diatom hosts, as both symbionts lack photosynthetic pigments and rely entirely on carbon fixed by their hosts (Schvarcz et al., 2022, 2024; Tschitschko et al., 2024). Similar subsurface distributions have been observed in UCYN-C at Station ALOHA (Schvarcz et al., 2022) and γ-24774A11 in the northern South China Sea (SCS; Chen et al., 2019; Lu et al., 2019; Shao and Luo, 2022).

Similar to the filamentous diazotrophs, UCYN-A and UCYN-B were mostly confined to the upper 30 m, matching observations from the northern SCS (Lu et al., 2019) and global oceans (Tang and Cassar, 2019). However, UCYN-A1 and UCYN-B were moderately abundant (∼10³ copies L⁻¹) in deep layers (Fig. S6M, N, and Q–T), potentially due to their low light saturation coefficients (Garcia et al., 2013; Gradoville et al., 2021; Shen et al., 2024). Moreover, the KSSW intrusion-induced low temperature (∼20.5 °C) at depths of 30–100 m may also favor the prevalence of UCYN-A1 (Figs. S4 and S6), the optimal growth temperature of which (∼22–24 °C) (Jiang et al., 2025; Moisander et al., 2010; Tang and Cassar, 2019) is moderately lower than that of the filamentous groups (∼26 °C for Trichodesmium and Hets) (Breitbarth et al., 2007; Fu et al., 2014; Tang and Cassar, 2019). However, such low temperatures may inhibit UCYN-B, which has an optimal growth temperature around 25–30 °C (Fu et al., 2014; Jiang et al., 2025; Tang and Cassar, 2019; Webb et al., 2009). Nevertheless, studies have indicated that UCYN-B can persist in cold environments, albeit at low abundances (Jiang et al., 2025; Tang and Cassar, 2019). Low-temperature tolerance of UCYN-B has also been demonstrated in the laboratory (Deng et al., 2022), however, the lowest threshold temperature used for the bottle incubation (25 °C) is substantially higher than the in situ temperatures (∼10–30 °C) in the field (Tang and Cassar, 2019). Therefore, the distribution of UCYN-B under such low-temperature conditions warrants further investigation.

4.2 Seasonal variations in diazotrophic distribution and nitrogen fixation rates on the ECS shelf

In addition to the spatial heterogeneity, our results also revealed pronounced seasonal dynamics in diazotrophic distribution on the ECS shelf (Figs. 3 and S6), primarily driven by the synergistic interplay between monsoonal hydrographic forcing and diazotrophic taxon-specific characteristics. In autumn, the southwest monsoon intensified the northward and shoreward penetration of Kuroshio and TSW onto the ECS shelf (∼80 %; Fig. S3B and C) (Oh et al., 2024), establishing warm (>24 °C), NO_x-deplete (<1 µM) water mass (Figs. 2 and S3) that favored the proliferation of Trichodesmium (Fig. 3A and C), while increasing the coastal shoaling of CDW/CW. It has been suggested that iron deficiency limits Trichodesmium growth (Berman-Frank et al., 2001). Given that SRP (0.2–0.4 µM) did not appear to be a limiting factor (Fig. 2D), the dominance of Trichodesmium on the ECS shelf is likely due to the high level of dissolved iron (0.76–30 nM) transported via the Kuroshio and TSW (Shiozaki et al., 2015; Su et al., 2015) or delivered via aerial dust deposition (Guo et al., 2014). This pattern is consistent with the extensive observations in Kuroshio-affected, island-adjacent waters where low NO_x but sufficient SRP and iron collectively support Trichodesmium growth (Chen et al., 2003; Shiozaki et al., 2015). In spring, however, the northeast monsoon facilitated the offshore expansion of the CDW/CW while weakening Kuroshio and TSW intrusion into the ECS shelf (Fig. S3A–C and E–G), resulting in decreased Trichodesmium abundance on the cold (<20 °C), NO_x-replete (>1 µM) northern ECS shelf (Fig. 3B and D). The high Trichodesmium abundance on the southern ECS shelf (Fig. 3B and D) may be related to its physiological preferences for warm, NO_x-deplete waters brought by Kuroshio and TSW intrusion (Fig. 2F and H) (Jiang et al., 2018, 2019). In addition to Trichodesmium, UCYN-A and UCYN-C were 4–5 orders of magnitude more abundant in spring than in autumn, but the seasonal variation was not significant for UCYN-B (Figs. 3 and S6M–V). The distribution of UCYN-A and UCYN-C was positively correlated to higher salinity in spring (>34.5 PSU; Figs. 7A and S4), a pattern also observed in previous studies (Li et al., 2021; Shiozaki et al., 2018), but the underlying mechanism is unknown.

Our data also showed seasonal dynamics in N₂ fixation on the ECS shelf, with higher NFRs in autumn relative to spring, particularly in Kuroshio-affected waters dominated by filamentous diazotrophs (i.e. stations 5, 20, and 30; Figs. 3, 6, and 7). The averaged autumn NFR (1.35 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$) was similar to what has been reported during summer on the ECS shelf (1.54 $\mathrm{nmol}\mathrm{N}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$) (Jiang et al., 2023a), despite the differences in NFR measurement between the two studies (i.e. dissolution versus bubble methods). The decreased NFRs in spring may be related to a substantial decrease in abundance and activity of filamentous diazotrophs (Figs. 3A and B and 7), which has been suggested to account for approximately 60 % of the bulk NFRs on the outer ECS shelf (Jiang et al., 2023a, b). On the other hand, the low temperature (14–20 °C) in spring could otherwise promote the prevalence of non-cyanobacterial diazotrophs (NCDs) and their contribution to N₂ fixation, particularly in CDW/CW-affected waters (Jiang et al., 2023b), but this pattern was not comprehensively analyzed in this study since only one phylotype of NCD (γ-24774A11) was quantified with qPCR (Figs. 3 and 6). Notably, the seasonal dynamics of NFRs showed negative, but not significant, correlation with transcriptional abundances of diazotrophs (Figs. 5 and 6 and Table S5), suggesting that the expression of nifH gene in diazotrophs may not be synchronized with actual nitrogenase activity (Turk-Kubo et al., 2012).

To assess the contribution of key diazotrophs to the N budget, we extrapolated water-column NFRs using published, cell-specific rates for the targeted diazotrophs (Table S6). We based our depth-integrated assessment on two assumptions: (1) cell-specific NFRs were considered constant throughout the water column, neglecting the regulation by environmental variables such as light, temperature and nutrient availability (e.g. Lu et al., 2018; Jiang et al., 2023a, b); (2) nifH-specific NFRs were derived from cell-specific rates using average polyploidy factors, despite the environmentally and/or phylogenetically driven variations in the cellular nifH gene copies (Gradoville et al., 2022; Sargent et al., 2016; Shao et al., 2023). Under these assumptions, the estimated and measured surface NFRs were highly correlated (Slope =1.04, R²=0.39, p<0.001; Fig. S10). Despite the uncertainties and low level of explained variance (Fig. S10), the significantly positive correlation between the two approaches validates the use of cell-specific NFRs in assessing regional N budget. Given the depth-integrated diazotroph abundances, the averaged NFRs on the ECS shelf ranged from 29.51 $\mathrm{\mu }\mathrm{mol}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$ in spring to 44.14 $\mathrm{\mu }\mathrm{mol}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$ in autumn (Table S6). The marked seasonal disparity may result from the below-detection abundances of filamentous diazotrophs at most of stations surveyed in spring (Fig. 3). Notably, the estimated NFRs in KSW-affected waters were considerably higher in spring (averaged at 344 $\mathrm{\mu }\mathrm{mol}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$ at stations 30, 41, and 42) than in autumn (averaged at 110 $\mathrm{\mu }\mathrm{mol}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$ at stations 18, 20, 28, 30, 40, and 41), yet both values fell within previously reported ranges on the southeastern ECS shelf (100–428 $\mathrm{\mu }\mathrm{mol}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$) (Jiang et al., 2023a; Sato et al., 2025b; Zhang et al., 2012). In terms of community composition, filamentous diazotrophs contributed 75 % and 37 % of the bulk NFRs in autumn and spring, respectively (Table S6), consistent with earlier observations (26 %–80 %) on the ECS shelf (Jiang et al., 2023b), the northern SCS and upstream Kuroshio (Chen et al., 2014), and the western Pacific (Bonnet et al., 2009; Kitajima et al., 2009). However, caution is warranted in interpreting these contributions, as the values (25 % and 63 %) for non-filamentous diazotrophs may be overestimated due to the low cellular nifH gene copy number used for UCYN-C and γ-24774A11. Compounding this uncertainty, the lack of cell-specific data for other NCDs means their full ecological role, particularly in cold environments where they sustain active N₂ fixation (Figs. 2 and 6) (Lin et al., 2013), is likely underrepresented. It should also be noted that the calculated NFRs reflect an estimate of the maximum diazotrophic potential, not the actual in situ activity, as the method cannot distinguish between metabolically active (N₂-fixing) and inactive (non-N₂-fixing) cells. To better elucidate the distinction between filamentous and non-filamentous diazotrophs in their contributions to N₂ fixation on the ECS shelf, highly-resolved NFR measurements completely spanning the four seasons and across fine scales of size-fractionated sampling should be emphasized in future studies (Jiang et al., 2023b).

4.3 Spatiotemporal dynamics in diazotrophic distribution reveals niche zonation

Our analysis revealed taxon-specific diazotrophic distribution that is closely linked to adjacent water masses (Fig. 7), identifying three distinct clusters among them in the multidimensional niche space of the ECS shelf (Figs. 8 and S8). Notably, several diazotrophic phylotypes were not restricted to narrow environmental windows but occurred across a broad spectrum of temperature and nutrient conditions. For instance, while Trichodesmium typically dominated in high-temperature, NO_x-deficient waters (e.g. Jiang et al., 2025; Nguyen et al., 2025; Tang and Cassar, 2019), our modeling results indicate its realized niches characterized by cool temperatures (μ_temperature=23.5 °C) and elevated NO_x levels (${\mathit{\mu }}_{{\mathrm{NO}}_x}=\mathrm{6.93}$ µM). This change in niche means reflects a displacement of Trichodesmium population on the shelf during the physical transport of water masses along the temperature and nutrient gradients. Meanwhile, the broad niche breadths (σ_temperature=3.5 °C; ${\mathit{\sigma }}_{{\mathrm{NO}}_x}=\mathrm{5.11}$ µM) also confirm the capacity of Trichodesmium to persist across these sub-optimal temperatures and variable nutrient conditions. This finding aligned with earlier studies reporting that Trichodesmium can survive, albeit at lower abundances, at temperature around 20 °C (Breitbarth et al., 2007; Rivero-Calle et al., 2016) and/or at NO_x levels of 10–20 µM, provided that SRP is sufficient to sustain growth (Knapp et al., 2012; Knapp, 2012). However, the realized niches here were slightly cooler in temperature but substantially higher and broader in NO_x relative to those of the open ocean (μ_temperature=24.2 °C, σ_temperature=3.2 °C; ${\mathit{\mu }}_{{\mathrm{NO}}_x}=\mathrm{2.5}$ µM, ${\mathit{\sigma }}_{{\mathrm{NO}}_x}=\mathrm{2.54}$ µM) (Brun et al., 2015). The expanded NO_x niche breadth of Trichodesmium on the ECS shelf appears to arise from intra- and/or inter-specific physiological variations. Specifically, intra-specific plasticity may enable endemic Trichodesmium populations to modulate N₂ fixation and growth by exploiting ambient SRP under fluctuating NO_x regimes (Knapp et al., 2012; Knapp, 2012), or utilizing NO_x directly (Boatman et al., 2018). Meanwhile, inter-specific distributions of Trichodesmium were evident, with T. erythraeum dominating coastal waters whereas T. thiebautii prevailing on the ECS shelf (Jiang et al., 2018; Zhang et al., 2019), possibly owing to species-specific physiological optima and NO_x tolerance thresholds (Carpenter et al., 1993; Confesor et al., 2022; Rodier and Le Borgne, 2008, 2010). These adaptive strategies collectively explain the broad niche underpinning the widespread distribution of Trichodesmium on the ECS shelf (Fig. 3).

Our analysis showed that both UCYN-A1 and UCYN-A2/A3/A4 favored high-salinity (>34 PSU) and moderate-temperature niche (∼24 °C), consistent with previous findings (Jiang et al., 2025; Tang and Cassar, 2019). It has been suggested that UCYN-A1 is an open-ocean ecotype and UCYN-A2 a coastal ecotype (Henke et al., 2018; Turk-Kubo et al., 2017, 2021). However, we found a co-occurrence pattern between these two ecotypes, with UCYN-A2/A3/A4 exhibiting lower abundance and narrower vertical distribution range (Figs. 3 and S6M–P) (Cabello et al., 2016; Gérikas Ribeiro et al., 2018; Wen et al., 2022). The distinct UCYN-A2 distributions between this and previous studies may be due to the phenotypical plasticity of the UCYN-A2-haptophyte symbioses, allowing them to inhabit either in the oligotrophic open oceans with small cell appearance (4–5 µm) or large one (7–10 µm) in the eutrophic coastal regions (Cabello et al., 2016; Hagino et al., 2013). Although genes related to N₂ fixation and energy production are actively transcribed in both UCYN-A ecotypes, the mechanisms underlying their niche overlap and differentiation require further studies (Muñoz-Marín et al., 2023; Nguyen et al., 2025).

Although clustered together in niche space (Fig. 8G), UCYN-B differed slightly from UCYN-A by occupying relatively nutritious niches (Fig. 8C and D), allowing it to inhabit both the warm, saline, oligotrophic KSW, and the cool, equally saline but eutrophic KSSW (Fig. S6Q–T). The latter seems in contrast to the well-known preference of UCYN-B for warm, oligotrophic waters (Jiang et al., 2025; Webb et al., 2009). Nevertheless, global metadata analysis suggests that UCYN-B may survive under low temperatures (20–23 °C) and high NO_x (up to 20 µM) and SRP levels (1–2 µM) (Tang and Cassar, 2019). Both field and laboratory studies have shown that UCYN-B remains viable under nitrate-rich conditions (5–20 µM) when SRP is sufficient (Knapp et al., 2012; Knapp, 2012; Turk-Kubo et al., 2018). These results indicate that UCYN-B has the capacity to cope with lower temperatures and elevated nutrient concentrations. Moreover, there was a clear niche differentiation between UCYN-B and Trichodesmium (Fig. 8), which has been reported in the northern SCS (Zhang et al., 2024) and the western Pacific (Chen et al., 2019), but this niche differentiation was not observed in the south Pacific (Tang and Cassar, 2019).

The third cluster comprised diatom-diazotroph symbioses (Hets, UCYN-C and γ-24774A11), which preferred high-temperature, low-NO_x and moderate-SRP environments (Fig. 8), causing their distribution confined to KSW enriched with terrestrial and/or aerosol-derived SRP (Figs. 2D and I and 3). Such environments have been shown to enhance symbiotic associations between Hets and diatom hosts (Tuo et al., 2014) and to promote UCYN-C growth (Turk-Kubo et al., 2015) in tropical oceans. Accordingly, a significant positive correlation in distribution and abundance between Het-1 and UCYN-C has been documented in the western Pacific (Chen et al., 2019). Compared with the niche conditions suitable for open-ocean populations (temperature: 24.2 °C; NO_x: 1.38 µM) (Brun et al., 2015), Hets inhabited niches of slightly higher temperature but lower NO_x levels (Fig. 8A and C). This result is different from previous findings based on studies analyzing all the Richelia diatom hosts, including Rhizosolenia clevei (Het-1 host), Hemiaulus hauckii (Het-2 host) and H. membranaceus (Brun et al., 2015; Luo et al., 2012), which may have distinct temperature and nutrient optima (Stenegren et al., 2017; Tuo et al., 2014, 2021). In fact, the different niche preferences for temperature and SRP between Het-1 and Het-2 were verified in this study (Fig. 8A and D), as well as in other studies (Foster et al., 2007; Tuo et al., 2021; Stenegren et al., 2017). Furthermore, the modeled realized niche space based on regional environmental gradients (e.g. temperature: 17.6–27.6 °C) may be narrower than that derived from global datasets (16.5–31.6 °C) (Brun et al., 2015; Jiang et al., 2025; Tang and Cassar, 2019). Thus, the prediction of realized niche of Hets may be affected by a variety of factors, including physiological traits of the diazotrophs, their hosts, and the scale of spatial sampling. Finally, our niche modeling result suggests a cosmopolitan distribution of γ-24774A11 in oligotrophic surface waters (Fig. 8) (Cheung et al., 2020; Shiozaki et al., 2018; Tschitschko et al., 2024), despite the fact that it is frequently found in environments with elevated NO_x concentrations (Bird and Wyman, 2013; Shao and Luo, 2022; Shiozaki et al., 2014). This feature may be related to its diatom host, which acquires a substantial fraction of fixed nitrogen from γ-24774A11 (Tschitschko et al., 2024). Nevertheless, in high-NO_x environments, the host might leverage ambient NO_x for growth, thereby partially reducing its metabolic dependency on the symbiont. To illuminate the mechanisms governing the biogeographic distribution of the γ-24774A11-diatom association (Cornejo-Castillo and Zehr, 2021; Shao and Luo, 2022; Tschitschko et al., 2024), more systematic studies are needed.

5 Conclusions

In this study, we conducted a cross-season survey on the distribution and activity of the major diazotrophic phylotypes on the continental shelf of ECS, alongside hydrographic analysis and integrated niche modeling. Our results demonstrate pronounced spatiotemporal heterogeneity in composition, distribution, and activity of the diazotrophs that is closely related to distinct water masses. The patterns of taxon-specific niche zonation underscore how dominant physical forcing (e.g. Kuroshio intrusion) may shape the diversity and biogeography of diazotrophs in this marginal sea with complex land-ocean interactions. Consequently, diazotrophic abundances and activities are relatively higher on the outer ECS shelf frequently intruded by the warm, saline, N-limited Kuroshio and TSW. The comprehensive biological, hydrographic, and modeling datasets presented here are, therefore, of utmost importance in assessing the dynamics of diazotrophic communities on the ECS shelf and its contribution to the regional N budget. Moreover, our adoption of the univariate MaxEnt-GAM framework provides ecological insights that are less accessible through conventional multivariate approaches, elucidating trade-offs among co-occurring diazotrophic phylotypes (e.g. the separation of UCYN-A from diatom-diazotroph symbioses) and broadening global phytoplankton niche domains. However, due to the lack of summer and winter survey data, the current niche modeling may not be able to capture the full annual ecosystem dynamics, potentially affecting the modeling robustness and accurate estimate of realized niches to some extent. Additional uncertainties may arise from the omission of key environmental drivers such as light intensity, iron availability, and biotic interactions (e.g. microzooplankton grazing, competition from non-diazotrophic microorganisms), which strongly influence diazotrophic abundance and activity. Moving forward, to refine the niche modeling and resolve the spatiotemporal variations of marine diazotrophs in dynamic marginal seas such as ECS, a multidisciplinary framework should be prioritized in future research to integrate approaches such as microscopic cell counts, comprehensive ecosystem monitoring and environmental parameter collection, and molecular diagnostics (e.g. qPCR, amplicon sequencing, and multi-omics profiling).