The 2023 research cruises were conducted along the 33° N latitude in the northeastern ECS, covering nine stations within the longitudinal range of 125–127.5° E (Fig. 1). The survey was performed from May 16th to 18th during the spring (R/V Saedongbaek, Chonnam National University), and from July 3rd to 7th during the summer and October 10th to 12th during the autumn (R/V Cheongyung, Chonnam National University). At each station, a conductivity-temperature-depth (CTD) profile was collected using a CTD (SBE 911plus; Sea-Bird Electronics, Bellevue, WA, USA), while chlorophyll fluorescence was assessed with a CTD equipped with a fluorescence sensor (WET Labs ECO-AFL/FL; WET Labs, Philomath, OR, USA). Odd-numbered stations followed a CTD-based depth scheme (1 m intervals). Biological and chemical samples were collected at standard water depths (i.e., 0, 10, 20, 30, 50, 75, 100, and 125 m) from all even-numbered stations. The water depths sampled were derived from the CTD profiles obtained in situ.

For 18S rRNA gene amplicon sequencing, 1.5 L samples were collected and filtered through 0.45 μm pore size cellulose nitrate membrane filters (47 mm; Whatman, Florham Park, NJ, USA) under gentle vacuum pressure (<100 mmHg) to collect biological material, and the processed filters were stored at −80°C in a deep freezer until further analysis. For cell counts, 0.5 L samples were mixed with Lugol’s iodine solution to a final concentration of 2% for preservation. The Lugol’s-fixed protist samples were stored in the dark at 4°C until further analysis.

For nutrient analysis, concentrations of nitrate + nitrite (NO₃ + NO₂), ammonium (NH₄), phosphate (PO₄), and silicic acid (Si[OH]₄) were determined using a nutrient auto-analyzer (New QuAAtro39; SEAL Analytical, Fareham, UK). Each seawater sample, totaling 15 mL, was filtered through a glass microfiber GF/F filter (nominal pore size 0.7 μm, 47 mm; Whatman) using an acid-washed filtration unit under gentle vacuum, and collected in a polypropylene Falcon tube. Reference standards for each nutrient were used for quality control. For stoichiometric ratio calculations, when PO₄ concentrations precluded reliable ratio calculation, the N:P and Si:P ratios were assigned the maximum N:P and Si:P values observed within the same season, respectively. The robustness of this approach was confirmed by a sensitivity analysis on both ratios across three detection limit (DL)-handling treatments (Supplementary Fig. S1). To assess their sensitivity to the DL-handling, the environmental vectors (envfit) were likewise re-fitted onto the non-metric multidimensional scaling (NMDS) under the same three treatments using both raw and log₁₀-transformed ratios, with largely consistent vector significance (Supplementary Table S1).

For 18S rRNA gene sequencing analyses, DNA was extracted from the collected filters using the DNeasy PowerSoil Kit (Qiagen, Hilden, Germany), following the manufacturer’s instructions. The extracted DNA was quantified with a Qubit fluorometer using Quant-iT PicoGreen (Invitrogen, Carlsbad, CA, USA) and stored at −80°C until use in polymerase chain reaction (PCR).

Sequencing libraries were prepared following the Illumina Metagenomic Sequencing Library protocols for gene amplification (San Diego, CA, USA). The V4 region of the 18S rRNA gene was PCR-amplified using the universal forward (TAReuk454FWD1; 5′-CCA GCA SCY GCG GTA ATT CC-3′) and reverse (TAReukREV3; 5′-ACT TTC GTT CTT GAT YRA-3′) primers (Stoeck et al. 2010). PCR was carried out in 1× PCR buffer with 1 μL genomic DNA template (5 ng μL⁻¹), 0.3 μM of each primer, dNTP mixture (0.2 mM of each dNTP), 0.5 U Takara Ex Taq polymerase (Takara, Osaka, Japan), and DNase-free water (Bioneer, Daejeon, Korea) to reach a final volume of 10 μL. Negative controls were included for all PCRs by replacing the DNA template with an equal volume of DNase-free water. The thermocycler protocol for the initial PCR included 3 min at 95°C, followed by 25 cycles of 30 s at 95°C for denaturation, 30 s at 55°C for annealing, and 30 s at 72°C for extension; followed by a final 5-min extension at 72°C. The resulting products were purified using ExoSAP-IT Express PCR Product Cleanup (Applied Biosystems, Foster City, CA, USA). Following purification, 2 μL of each initial PCR product underwent a second amplification for final library construction using Nextera XT indexed primers in a 20 μL reaction with the same final reagent concentrations as the initial PCR, corresponding to 0.2 U of Takara Ex Taq polymerase per reaction. The thermocycler conditions for the second PCR were 3 min at 95°C; 8 cycles of 30 s at 95°C, 30 s at 55°C, and 30 s at 72°C; followed by a final 5-min extension at 72°C. PCR products were visualized by electrophoresis using 1 μL of amplicons on 1% agarose gels at 110 V for 30 min. No visible amplification was detected in the negative controls. Amplicons were purified to a final volume of 30 μL using yesC PCR purification kit (GensGen, Busan, Korea). The purified final product was quantified using a Qubit 3 Fluorometer (Thermo Fisher Scientific Inc., Waltham, MA, USA). Paired-end (2 × 250 bp) sequencing was conducted at CJ BioScience (Seoul, Korea) using a single lane on the Illumina MiSeq platform.

The resulting paired-end reads were submitted to the NCBI Short Read Archive (PRJNA1194217), with per-sample accessions and metadata (Supplementary Table S2). Demultiplexed paired-end reads were processed in QIIME 2 v2021.4 (Bolyen et al. 2019) (Supplementary Text S1). Primer sequences were removed using q2-cutadapt with the same primer pair used for amplification, allowing an error rate of 0.4 and a minimum overlap of 3 bp (Martin 2011). Amplicon sequence variants (ASVs) were inferred using the q2-dada2 plugin with the following parameters: truncation quality score = 2, forward and reverse truncation lengths = 200 bp, maximum expected errors = 2 for both forward and reverse reads, number of reads used to learn error rates = 1,000,000, and chimera detection method = pooled (Callahan et al. 2016). This step denoised reads, merged paired-end reads, removed chimeras, and generated the ASV table and representative sequences. Taxonomic assignment was performed using a naïve Bayes classifier trained on the PR2 database v5.0.0 after extraction of the target 18S rRNA V4 region. The resulting ASV table, representative sequences (Supplementary Table S3), and taxonomy assignments were then exported from QIIME 2 for downstream analyses (Guillou et al. 2013). ASVs identified as non-target taxa (Metazoa, Fungi, and Streptophyta) were excluded, and the filtered table was normalized by proportional scaling, in which each ASV count was divided by the total read count of the sample (Supplementary Table S4).

Light microscopy was utilized to identify and count protist species. For microscopic enumeration, 500 mL subsamples preserved in acidic Lugol’s solution (final concentration of 2%) were concentrated 5–10 times using the 2-d settlement technique. After thorough mixing, either all or at least 300 cells from three to five 1-mL Sedgwick-Rafter counting chambers were identified morphologically and counted at 100× and/or 200× magnification.

To classify the spatio-vertical sampling stations into distinct clusters according to oceanographic characteristics, K-means clustering was conducted using R version 4.3.1 (R Foundation for Statistical Computing, Vienna, Austria). The input variables were temperature, salinity, and nutrient concentrations (i.e., NO₃ + NO₂ + NH₄, PO₄, and Si[OH]₄). Prior to clustering, each variable was standardized to account for variability due to differences in units. To determine the optimal number of clusters, the elbow method (Supplementary Fig. S2) was employed in conjunction with the NbClust package (Supplementary Table S5). Principal component analysis was subsequently applied to the standardized data to visualize the clustering results. Then, sampling stations were assigned to a specific water mass and grouped in each season by integrating the clustering outcomes with information on seasonal water masses within the region (Li et al. 2006, Kim et al. 2009, Lie and Cho 2016). Cluster stability was evaluated using non-parametric bootstrap resampling (B = 1,000) with the clusterboot function in the fpc R package, based on the Jaccard similarity between original and resampled clusters (Supplementary Table S6).

To assess the statistical significance of differences in oceanographic variables and biological factors, including alpha diversity and species richness, across different water masses, linear mixed-effects models (LMMs) were fitted using the lme4 package in R version 4.3.3. LMMs were fitted separately for each cruise with water mass as a fixed effect, depth as a covariate, and station as a random intercept. Post hoc Tukey-adjusted pairwise comparisons were applied using the emmeans package to assess differences between specific water masses when applicable.

Additionally, NMDS based on Bray-Curtis dissimilarities was used to visualize pairwise dissimilarities among communities at different 18S rRNA gene sequencing sampling stations (stress values are shown in the respective figures). Environmental factors were overlaid on the NMDS ordination using the “envfit” function (9,999 permutations) to illustrate their influence on protistan community distribution patterns. Both the NMDS and permutational multivariate analysis of variance (PERMANOVA) analyses were conducted using the vegan package in R version 4.3.3. The ASV data were analyzed to evaluate differences in community composition among water masses using PERMANOVA based on Bray-Curtis dissimilarities (9,999 permutations, adonis2 function), where R² represents the proportion of community variation explained by the grouping variable. Group dispersion homogeneity was tested using betadisper (9,999 permutations). p-values from both analyses across seasons were adjusted using the Benjamini-Hochberg procedure. To assess the sensitivity of biological interpretations to the choice of grouping variable, PERMANOVA was additionally performed using both K-means cluster assignments and final water mass classifications as grouping variables. Finally, the agreement between the microscopy- and 18S rRNA gene amplicon sequencing-based community structures was assessed using a Procrustes test on the corresponding Bray-Curtis principal coordinates (9,999 permutations), which indicated significant concordance (r = 0.622, p < 0.001).

The water masses of the northeastern ECS during the spring (May), summer (July), and autumn (October) of 2023 were grouped using K-means clustering based on oceanographic variables and information on seasonal water masses within the region (Fig. 2).

In the spring, the samples were clearly divided by K-means clustering into three distinct clusters corresponding to three water masses: YSW, KSW, and KSSW (Fig. 2A & B). At the westernmost sampling location, all sampling points, from the surface to the bottom layer, were categorized as YSW. In contrast, the surface sampling points in the warmer and more saline eastern region were identified as KSW, while the points at the lower depths of the eastern region were classified as KSSW.

In the summer, the clustering analysis initially divided water masses into two groups (Fig. 2C & D), one representing the sampling points at greater depths, corresponding to the KSSW, and the other representing the sampling points nearer the surface. However, based on regional water mass information, the surface cluster was further subdivided into CDW and ECS Shelf water. This subdivision aligned well with the distribution patterns of the sampling points observed in the clustering analysis.

In the autumn, the clustering analysis separated the sampling points into three groups, clearly distinguishing Cluster 2, corresponding to the KSW, from the lower depths of the eastern region (Fig. 2E & F). While Cluster 1 consisted of only surface sampling points from the western region, hydrographic knowledge led to the inclusion of the middle and bottom-layer sampling points from the same region into a single water mass corresponding to the YSW. The remaining sampling points in Cluster 3 were categorized as KSSW.

During the spring season, YSW displayed relatively cooler and less saline characteristics, while both KSW and KSSW exhibited significantly higher salinities, exceeding 34.0 (LMM, Tukey-adjusted p < 0.001) (Table 1, Fig. 3, Supplementary Table S7). In terms of nutrient concentrations, dissolved nitrogen levels in the surface waters were below 1 μM across all stations along the 33° N line (Fig. 4A). In contrast, concentrations of Si[OH]₄ (p < 0.01) and NH₄ (p < 0.001) were relatively higher in the surface of the western water mass (YSW) than the eastern water mass (KSW) (Table 1, Fig. 4B, Supplementary Table S7). However, PO₄ levels in these western surface waters were low, near the detection limit, indicating an inverse trend (Fig. 4C). Consequently, high N:P and Si:P ratios were observed in the surface waters during spring, particularly in YSW.

During the summer, the surface water mass CDW, located in the western region, exhibited a significantly lower average salinity, 30.3, and higher temperatures than the ECS Shelf water in the eastern region (both p < 0.001) (Table 1, Fig. 3). Regarding nutrient concentrations, NH₄ levels exhibited seasonal specificity, with average concentrations ranging from 0.1 to 0.2 μM across all water masses in the summer (Table 1). In the surface waters, total dissolved nitrogen and Si[OH]₄ were slightly higher in CDW than in ECS Shelf water, but PO₄ levels did not differ between the two water masses (Fig. 4).

During the autumn season, KSW in the east exhibited a notably high average temperature of 25.5°C and a slightly higher salinity compared to YSW in the west (both p < 0.01) (Table 1, Fig. 3). Notably, in autumn, water masses showed mixing near the surface and mid-layer throughout the region, which was attributed to intense wind-driven mixing. In terms of nutrient concentrations, dissolved nitrogen (p < 0.05) and silicate (p < 0.01) levels were elevated in YSW relative to KSW, remaining consistent with the summer observations (Fig. 4A & B). Ammonium concentrations were not detected at significant levels across any of the sampling points (Table 1). Conversely, PO₄ levels in the surface waters were relatively high in the surface of the eastern water mass during autumn (Fig. 4C).

When comparing the oceanographic characteristics of water masses across different seasons, KSSW, situated below the mixed layer, consistently displayed lower temperatures and higher salinities, which remained near 34.0 throughout the study period, relative to the warmer surface water masses (Fig. 3). This water mass also had higher dissolved nutrient concentrations (Fig. 4). Overall, nutrient levels in the water masses were elevated in autumn compared to spring and summer. When comparing the same water masses across different seasons, the N:P and Si:P ratios in autumn were generally lower than those in other seasons, which is attributed to a modest increase in PO₄ levels (Table 1).

To understand how the protist communities responded to the changing oceanographic factors, an NMDS analysis of community dissimilarity was conducted (Fig. 5, Supplementary Table S8). During the spring survey, the communities formed clusters generally corresponding to their respective water masses (Fig. 5A). The YSW communities were more distinctly separated from the other two. Specifically, YSW communities were more closely associated with chlorophyll fluorescence and NH₄, with similar trends for N:P and Si:P ratios, while temperature and salinity were more strongly correlated with KSW. The eastern bottom water mass, KSSW, showed a strong positive correlation with nutrient concentrations.

During the summer survey, pronounced structural differences were observed between communities in the CDW and KSSW water masses (Fig. 5B). The CDW communities exhibited strong positive correlations with Si:P and temperature, with a similar trend for N:P. Conversely, the KSSW communities showed correlations with nutrient concentrations and salinity. Additionally, NH₄ concentrations were more strongly associated with surface sampling points from the ECS Shelf water and CDW than the deeper sampling points from the KSSW water mass.

In the autumn sampling, protist communities were again distinctly separated according to the water mass they inhabited, each of which showed unique correlations with oceanographic factors (Fig. 5C). The YSW communities exhibited positive correlations with chlorophyll fluorescence and N:P and Si:P ratios, whereas the KSW communities were associated with relatively warmer and more saline conditions. Meanwhile, communities in the KSSW showed correlations with salinity and nutrient concentrations, as was observed in other seasons.

PERMANOVA detected significant differences in protist community composition among water masses across all three seasons (R² = 0.215–0.272, p < 0.001) (Supplementary Table S9). Unlike summer, group dispersions were heterogeneous in spring and autumn (betadisper, p < 0.01); however, the PERMANOVA effect sizes (R² = 0.258 and 0.272, both p < 0.001) suggest that water mass differences are also evident in these seasons. These compositional differences were consistent regardless of whether K-means cluster assignments or final water mass classifications were used as grouping variables, indicating that the observed differences are robust to the choice of grouping variable (ΔR² ranged from −0.019 to 0.051).

As assessed using the Shannon diversity index and species richness based on ASVs, protist communities in the northeastern ECS in 2023 displayed seasonal variations in diversity that could reflect the influence of distinct water masses (Fig. 6A & B). In the spring, the Shannon diversity index for the protist communities in YSW was similar to that observed in the KSW communities (LMM, Tukey-adjusted p = 0.611), but numerically lower than that of the KSSW communities (Fig. 6A, Supplementary Table S10). In the summer, compared with the CDW, protist diversity was significantly higher in the KSSW (p < 0.05) and tended to be higher in the ECS Shelf water (p = 0.07). In the autumn, all three communities within the region showed high diversity values, all exceeding 4. Similarly, species richness, inferred from ASV numbers, was higher in the KSSW than in the YSW and KSW in the spring, although these differences were not statistically significant (p = 0.55–0.90) (Fig. 6B, Supplementary Table S11). During the summer, protist communities in the CDW tended to exhibit lower species richness than those of the other water masses (p = 0.06–0.08), representing the lowest richness observed among the water masses across all seasons. Diatoms showed the highest richness in the ECS Shelf water compared to other water masses across all observed seasons. In the autumn, species richness was highest in KSSW, followed closely by the KSW and YSW communities, and these differences were not statistically significant (p = 0.82–0.83). Across the three sampled seasons, the contribution of Rhizaria and Syndiniales was higher in deeper KSSW than in surface water masses.

The community structure within each water mass by season is summarized in Fig. 6C. Lower relative abundances of dinoflagellates were seen in autumn compared to spring and summer, while diatoms were particularly prevalent in summer, especially in the ECS Shelf water. In the YSW and KSW water masses during spring, Kathablepharidacea exhibited notably higher relative abundances, whereas Rhizaria were more abundant in all water masses during autumn. Regardless of season, the relative abundance of dinoflagellate reads was higher in the surface water masses compared to the bottom water mass, KSSW. In spring, dinoflagellates accounted for 67.1% of the reads from the surface water at station E42 in YSW, and in summer, they represented 72.3% of the reads at a depth of 20 m at station E50 in the ECS Shelf water (Fig. 7). However, in the KSSW water mass during autumn, they constituted only 11.9% of the community on average. Diatoms exhibited significantly high read counts in the ECS Shelf water during summer, with a particularly notable relative abundance of 49.7% observed at the surface of station E50. Although Kathablepharidacea was not a prevalent group in the overall protist community, it exhibited notable read counts, averaging 8.2% across specific stations, in the surface communities of YSW and KSW during spring (i.e., at the 0 and 10-m depths at stations E42, E44, and E46). The Rhizaria group exhibited a high average read abundance of 32.7% in the KSSW water mass during autumn. Additionally, in CDW during summer (specifically at the 0 and 10-m depths at station E42), it reached high relative abundances, up to 49.1% at 10 m. Read abundances of Chlorophyta and Stramenopiles (except diatoms) were notably high in the KSW water mass (Figs 6C & 7). Chlorophyta represented 11.8% in spring and 16.5% in autumn within the KSW, but only 5.9% in the same region (i.e., E46–E50) when it was occupied by ECS Shelf water during summer. Similarly, Stramenopiles constituted 7.9% in spring and 7.8% in autumn within the KSW, but only 3.9% in the ECS Shelf water during summer. Finally, the relative abundances of the parasitic Syndiniales were consistently high across all communities, with particularly high levels in the KSSW communities regardless of the season (spring, 59.8%; summer, 49.1%; and autumn, 45.4%) (Figs 6C & 7).

Community compositions at the species level in the water masses were examined using cell counts (Fig. 8, Supplementary Table S12). During the study period, total protistan cell abundances were highest in the autumn, regardless of water mass, followed by the YSW communities in the spring. Cell abundances were greater in YSW during the spring and autumn seasons than in CDW, which occupied the same stations in the summer. Moreover, in spring and autumn, the cell abundances of protists were generally about twice as high in YSW compared to KSW.

Overall, large cryptophyte (>10 μm) cells dominated most communities in the northeastern ECS, with the highest abundances observed in YSW and KSW in autumn, and in YSW in spring, at 77.6, 38.5, and 22.9 cells mL⁻¹, respectively. Small cryptophytes (<10 μm) were particularly abundant in spring and summer, contributing 17.2 cells mL⁻¹ in YSW and 6.1 cells mL⁻¹ in CDW, respectively. Gyrodinium spp. had the highest abundance among dinoflagellates, representing 41.7, 11.9, and 7.7 cells mL⁻¹ in the YSW and KSW communities in autumn, and in the ECS Shelf water community in summer, respectively. Gymnodinium sp. (>30 μm) exhibited high abundances in all three water masses in autumn and in the YSW community in spring, while Prorocentrum sp. showed a particularly high abundance, 10.2 cells mL⁻¹, in YSW in spring but was rare otherwise. The abundance of diatoms was highest in YSW communities, with Skeletonema costatum showing the highest abundance, 11.5 cells mL⁻¹, in spring. In autumn, Chaetoceros socialis was the most commonly encountered diatom at 31.4 cells mL⁻¹, followed by Nitzschia sp. and Cylindrotheca closterium. Finally, the abundance of oligotrich ciliates (20–40 and 40–60 μm) was also highest in autumn, ranging from 8.2 to 14.1 cells mL⁻¹.

The locations and identities of seasonal water masses, classified using K-means clustering in conjunction with established knowledge of the regional water masses, closely correspond to those reported in previous studies conducted over multiple years (Li et al. 2006, Kim et al. 2009, Lie and Cho 2016). In 2023, the spring and autumn water masses in the northeastern ECS were influenced by KSW (i.e., Tsushima Warm Current and Jeju Warm Current) in the east and by YSW in the west (Pang et al. 1992, Kim et al. 2009, Lee et al. 2014b). During autumn, strong winds at the time of the survey caused mixing down to the mid-layer at shallow western stations, pushing stratification to a deeper depth. Although K-means clustering classified the western lower layer as a different water mass, we grouped it with YSW based on established hydrographic knowledge of the region (Kim et al. 2009, Choi et al. 2022). Nevertheless, considering the high nutrient concentrations in some samples classified as being in the YSW region, it appears that the YSW water mass was partially mixed with the KSSW water mass. In contrast, during summer, as in other studies, different water masses were found to influence the northeastern ECS, with CDW forming most of the western region’s surface water (Pang et al. 1992, Kim et al. 2009, Lie and Cho 2016, Liu et al. 2021). The eastern region in this study appears to be strongly influenced by ECS Shelf water, which flows northeast along the topography of the ECS Shelf (Ichikawa and Beardsley 2002, Yanao and Matsuno 2013).

The Changjiang discharge water expands northward starting in spring, reaching its maximum extent in summer, and as it progresses from the estuary to the northeastern ECS, nutrient depletion occurs, resulting particularly in phosphorus depletion by the time it reaches the western region of our study area (Wang et al. 2003, Liu et al. 2016, Kim et al. 2023). The underlying cause of the persistent elevated N:P ratio observed in the northeastern ECS remains unclear (Wang et al. 2003, Moon et al. 2021). However, it is plausible that PO₄ depletion in the YSW water mass is driven by a relatively limited supply of additional PO₄ in the region, while bioavailable nitrogen may be supplemented by various potential sources, such as atmospheric deposition of nitrogen from pollution or nitrogen fixation by cyanobacteria (Kim et al. 2011, Zhang et al. 2012). In spring, the YSW water mass exhibited significantly elevated NH₄ concentrations, despite a high N:P ratio, along with increased chlorophyll fluorescence values. These observations suggest that elevated NH₄ levels are a characteristic feature that may shape the protist community in the YSW water mass during this season, and chlorophyll fluorescence likely reflects a phytoplankton response to the nutrient conditions.

During summer, intense solar radiation elevates surface water temperatures, while increased terrestrial runoff severely lowers salinity in the western region, resulting in the development of strong stratification—a pattern commonly observed in the ECS (Kim et al. 2009, Son and Choi 2022). Distinctly higher NH₄ concentrations were detected throughout the water column in summer, while NO₃ + NO₂ and chlorophyll levels remained similar across seasons, suggesting that this is a seasonal characteristic.

The elevated surface nutrient concentrations we observed in autumn are likely attributed to wind-driven vertical mixing, facilitating the upward transport of nutrients from deeper layers to the surface (Kim et al. 2007, 2009). Furthermore, in autumn, the relatively high nutrient concentrations in the bottom waters are presumed to result from an increased influx to the west due to an intensification of the Kuroshio Current (Guo et al. 2006, Isobe 2008).

In the western surface waters, nutrient concentrations exhibited notable seasonal fluctuations compared to those observed in the east. In the classified water masses, the eastern region consistently displayed a boundary between the upper and lower layers, while in the west, a single homogeneous water mass extended from the top to the bottom during spring and autumn. This suggests that vertical mixing likely occurred more frequently in the western region across the observed seasons, driven by the shallower depth in this area, which facilitated enhanced nutrient transport to the upper layers (Kim et al. 2009, Wang et al. 2019). These seasonal changes in the physico-chemical parameters are important factors influencing the distribution of protist communities.

The oceanographic conditions observed in 2023 are not unique to that particular year but reflect a recurring pattern consistently observed in similar seasons and locations in previous studies (Lie et al. 2000, Wang et al. 2003, Kim et al. 2009, Lee et al. 2014b). Notably, the persistence of phosphate concentrations below 0.1 μM in the western region during spring appears to be a distinctive seasonal feature of this area, indicative of a region-specific pattern of nutrient limitation (Wang et al. 2003, Kim et al. 2020, 2023).

The NMDS analysis highlights a prominent feature of the study region: a persistent phosphorus deficiency in the western surface waters, specifically associated with the YSW and CDW water masses, throughout the study period. Previous studies have documented protist communities’ responses to phosphorus limitation during spring in this region (Kim et al. 2020, 2023). This research expands on those findings by demonstrating that P depletion patterns persist and that this, in turn, could influence the structure of protist communities in the region. In contrast, the NMDS analysis revealed that the typical oceanic-origin KSW community, characterized by high temperature and salinity during the spring, exhibits correlations opposite those of the YSW community. This suggests that the protist communities in the northeastern ECS may diverge significantly between the two upper water masses during the spring season. Similarly, Kim et al. (2023) analyzed that the western region of the northeastern ECS was strongly associated with the Si:P ratio, ammonium, and chlorophyll-a concentration, whereas the eastern region, particularly influenced by Jeju Warm Current Water (JWCW), was associated with temperature and salinity. In this study, YSW and KSW corresponded to western and eastern regions, respectively, suggesting that similar oceanographic characteristics may have repeatedly contributed to the formation of comparable protist community structures.

In spring and autumn, protist communities in the eastern surface waters are expected to be shaped by the high temperatures and salinity typically associated with the Kuroshio branch. However, in summer, the NMDS analysis revealed a positive relationship between ECS Shelf water communities and temperature, suggesting that this community is primarily driven by Shelf water dynamics, distinguishing it from those in other seasons. Additionally, during spring and autumn, surface water communities in the western and eastern regions occupied contrasting positions in the NMDS analysis. However, in summer, the eastern water mass appears to be less distinct from the western water mass, likely because several oceanographic variables (i.e., NH₄, temperature, and nutrient ratio) become less different between the two regions. This indicates that, in summer, the eastern surface water community is more strongly influenced by Shelf water than by the Kuroshio Current, distinguishing it from other seasonal patterns (Zhang et al. 2019).

The observed seasonal variations in protist community structure in the northeastern ECS reflect the dynamic interactions between water masses and environmental factors, such as temperature, nutrient availability, and mixing regimes. The data show a clear seasonal shift in community diversity and composition, highlighting the ecological influences of distinct water masses.

In spring, the Shannon diversity index and species richness were greater in KSSW, revealing that it hosts a more diverse protist community than YSW and KSW. This corresponds to previous findings that, despite relatively low biomasses in bottom water masses, high biodiversity is supported by the diverse survival strategies employed by protists, such as cyst formation and phagotrophy (Kremp et al. 2009, Jang et al. 2019). Additionally, lower diversities in YSW and KSW may be the result of more dynamic and variable conditions near the surface, potentially leading to fewer niche opportunities for diverse protist taxa.

During summer, a marked decline in both diversity and species richness was observed in the CDW, likely due to the influx of freshwater, which decreases the salinity through dilution, thereby impacting protist communities. This contrasts with the higher diversities seen in the ECS Shelf water and KSSW, indicating that these water masses play a stabilizing role during the summer season. High diversities in autumn across all water masses suggest a convergence of optimal conditions in terms of temperature and nutrient cycling, supporting rich protist communities.

The community composition data showed seasonal influences on specific taxa. Dinoflagellates were prevalent during the spring and summer seasons in both surface and sub-surface waters, particularly in YSW and ECS Shelf water. Despite the highly diverse, species-specific traits of the dinoflagellate group, making broad generalizations challenging, their high relative abundance in water masses with high N:P ratios, such as YSW in spring and ECS Shelf water in summer, where phosphorus is depleted, may be explained by their mixotrophic strategy. This adaptive mechanism could provide a competitive advantage in nutrient-depleted environments by enabling the acquisition of nutrients from organic matter through phagotrophy (Li et al. 2000, Jeong et al. 2021). Their decline in autumn, however, which coincided with an increase in Chlorophyta and Rhizaria, points to a seasonal succession likely driven by changing nutrient regimes and mixing dynamics in the autumn season. Furthermore, Rhizaria prevalence could serve as a bioindicator of KSSW; however, this hypothesis requires validation through advanced methodological approaches (Dufrêne and Legendre 1997, Pawlowski et al. 2018).

The significant increases in the relative abundance of diatoms in summer, particularly in the ECS Shelf water, suggest that these primary producers thrive in warmer temperatures and higher nutrient conditions. Considering the relatively low rDNA copy numbers of diatoms among protists (Yarimizu et al. 2021), the ecological significance of diatoms within this water mass may have been underestimated. Moreover, the macronutrient concentrations in the ECS Shelf water, including silicate at over 9 μM, along with its average temperature of 22.2°C, both surpassing the general half-saturation constants for nutrient uptake in diatoms, suggest that conditions may have been conducive to diatom proliferation (Eppley et al. 1969, Paasche 1973, Suzuki and Takahashi 1995).

The contrasting patterns of Kathablepharidacea, i.e., its notable presence in spring, and Rhizaria, which was particularly prevalent in autumn, also illustrate niche partitioning among protists in response to environmental changes. Kathablepharidacea, though a minor component of the community, appears to flourish in the surface waters during spring, potentially functioning as generalist grazers (Kwon et al. 2017, Cui et al. 2021).

The consistent prevalence of parasitic Syndiniales across all water masses and sampled seasons, particularly in KSSW, suggests an important ecological role in regulating protist population dynamics, likely through the parasitism of major host taxa such as dinoflagellates and diatoms (Anderson and Harvey 2020, Yang et al. 2024). Their higher relative abundance in KSSW may be attributable to the relatively stable and nutrient-rich sub-surface conditions, which could support their hosts and, thus, their parasitic lifecycle.

Although comparable spatiotemporal studies on protist communities in this region remain limited in terms of both cell counts and environmental DNA approaches, hindering direct interannual comparisons with previous years, available data from May 2022 indicate a trend of increasing species richness toward the eastern region, where the influence of the Kuroshio-branch Current (i.e., JWCW) is more pronounced (Kim et al. 2023). Overall, this interannual comparison further supports that water mass dynamics in the region likely play a key role in shaping protist community structure. Although read proportions do not directly reflect cell abundance owing to rDNA copy-number variation, a Procrustes test confirmed that the community structures inferred from microscopy and 18S rRNA gene amplicon sequencing were significantly concordant. Consistent with this, previous studies have demonstrated that the 18S rRNA gene amplicon sequencing approach provides reliable ecological interpretations of protist communities (Pawlowski et al. 2016, Piwosz et al. 2020), supporting the validity of the community-level inferences drawn here.

The distribution of protist species across the different water masses, as determined through cell counts, revealed clear seasonal and spatial variability. The YSW consistently exhibited higher cell abundances in both spring and autumn compared to the CDW and KSW, indicating that YSW provides a more stable environment conducive to protist proliferation. Moreover, the stratification patterns seen in the summer, especially with the CDW, likely inhibit nutrient mixing, creating oligotrophic conditions that suppress overall protistan cell abundance. By contrast, autumn mixing events due to wind-induced turbulence bring nutrients to the surface, explaining the peaks in cell abundances during this period, particularly in the YSW and KSW communities.

The seasonal variation in community composition further highlights the adaptability of protist species to changing environmental conditions. For example, cryptophytes exhibited significant seasonal shifts, with larger cryptophytes (>10 μm) dominating in autumn, particularly in YSW and KSW, while smaller cryptophytes (<10 μm) thrived in spring and summer. This shift could be explained by the nutrient environment: larger cryptophytes may be more competitive under nutrient-rich conditions following autumn mixing, whereas smaller cryptophytes could exploit more stable, oligotrophic conditions in spring and summer (Raven 1998).

Similarly, dinoflagellates—particularly Gyrodinium and Gymnodinium species—exhibited higher abundances in several water masses, including those where phosphorus was limited, despite being larger than other protists. The phagotrophic nature of many dinof-lagellates could enable them to thrive in environments where nutrients fluctuate, offering a competitive advantage when inorganic nutrients are scarce but organic sources of carbon are abundant (Arenovski et al. 1995, Jeong et al. 2010). Correspondingly, when the cells of Gyrodinium and Gymnodinium were abundant, our data also showed high abundances of cryptophytes, diatoms, and dinoflagellate Prorocentrum species, which are generally recognized as prey species of these dinoflagellates (Kim and Jeong 2004, Yoo et al. 2009, Lee et al. 2014a, Kang et al. 2020). This trait may explain their ability to persist across multiple seasons and water masses. Moreover, the presence of Prorocentrum spp. exclusively in the YSW during spring further highlights the role of specific nutrient regimes in supporting certain dinoflagellate species, possibly driven by shifts in the N:P ratio, as PO₄-limited conditions have been shown to favor Prorocentrum growth (Ou et al. 2008).

Diatoms dominated in certain water masses during different seasons. For example, Skeletonema costatum dominated in the spring YSW community, while Chaetoceros socialis dominated in the autumn YSW community. These diatom blooms likely contribute significantly to primary production, underpinning the seasonal food web dynamics in the region (Zhou et al. 2019). Furthermore, the high abundance of oligotrich ciliates in autumn across all water masses suggests their role in controlling protist community structure by preying on smaller protists and bacteria. Their peak in autumn could also reflect increased grazing pressure on phytoplankton populations, contributing to the observed shift in community compositions (Haraguchi et al. 2018).

The dominance of certain taxa, such as cryptophytes and diatoms, during nutrient-rich periods indicates their key roles in primary production, while the higher abundance of heterotrophic groups, like dinoflagellates and ciliates, in autumn reflects complex trophic interactions. Understanding these patterns helps illuminate how environmental changes, such as shifts in nutrient regimes, may affect future protist communities in the northeastern ECS.

This study reveals that spatiotemporal shifts in water masses within the northeastern ECS lead to intricate physico-chemical changes, which are key drivers of protist community structure. In particular, the potential phosphorus limitation in the western waters, along with the elevated N:P ratios that accompany it, may have influenced community structure, especially in YSW. Protist groups, such as dinoflagellates and diatoms, exhibited clear seasonal prevalence patterns, with dinoflagellates thriving in the nutrient-fluctuating environments of YSW in spring and diatoms being particularly prevalent in the ECS Shelf water during summer. Additionally, nutrient mixing in autumn enhanced protist diversity and abundance, particularly in YSW, where Chaetoceros socialis and other diatoms flourished. These findings underscore the importance of water mass dynamics in shaping protist community structure and provide insights into how these communities might respond to future environmental changes. A deeper understanding of these interactions will be essential for predicting shifts in marine ecosystems, particularly under ongoing environmental changes driven by climate variability and anthropogenic impacts. Given the scarcity of comparable interannual studies in this region, the present results also contribute valuable context for future comparisons and long-term ecological assessments of protist communities in the northeastern ECS.