Protist Interactions and Seasonal Dynamics in the Coast of Yantai, Northern Yellow Sea of China as Revealed by Metabarcoding

FU Yingjun, ZHENG Pengfei, ZHANG Xiaoli, ZHANG Qianqian, *,and JI Daode, *

Protist Interactions and Seasonal Dynamics in the Coast of Yantai, Northern Yellow Sea of China as Revealed by Metabarcoding

FU Yingjun1), ZHENG Pengfei2), 3), ZHANG Xiaoli2), 3), ZHANG Qianqian2), 3), *,and JI Daode1), *

1),,264003,2),,,264003,3),,266071,

Facilitated by the high-throughput sequencing (HTS) technique, the importance of protists to aquatic systems has been widely acknowledged in the last decade. However, information of protistan biotic interactions and seasonal dynamics is much less known in the coast ecosystem with intensive anthropic disturbance. In this study, year-round changes of protist community composition and diversity in the coastal water of Yantai, a city along the northern Yellow Sea in China, were investigated using HTS for the V4 region of 18S rDNA. The interactions among protist groups were also analyzed using the co-occurrence network. Data analyses showed that Alveolata, Chlorophyta, and Stramenopiles are the most dominant phytoplanktonic protists in the investigated coastal area. The community composition displayed strong seasonal variation. The abundant families Dino-Group-I-Clade-1 and Ulotrichales_X had higher proportions in spring and summer, while Bathycoccaceae exhibited higher ratios in autumn and winter. Alpha diversities (Shannon and Simpson) were the highest in autumn and the lowest in spring (ANOVA test,0.05). Nutrients (SiO42−, PO43−), total organic carbon (TOC), and pH seemed to drive the variation of alpha diversity, while temperature, PO43−and TON were the most significant factors influencing the whole protist community. Co-variance network analyses reveal frequent co-occurrence events among ciliates, chlorophytes and dinoflagellate, suggesting biotic interactions have been induced by predation, parasitism and mixotrophy.

coastal zone; biotic interaction; high-throughput sequencing; northern Yellow Sea; protist diversity; seasonal dynamic

1 Introduction

Protists have long been recognized as either photoautotroph (typically diatoms and coccolithophores) or heterotrophic predators, which are major contributors to the ocean productivity (Falkowski, 2007) and carbon consump- tion (Sherr and Sherr, 2002). In addition to the photosynthetic and heterotrophic growth, protists exhibit a range of other trophic modes, including mixotrophy, parasitism, sym-biosis, osmotrophy, and saprotrophy (Worden, 2015). These complex behavioral strategies and organismal interactions have led to challenges for the comprehensive model of the marine carbon cycle (Worden, 2015). Protists are morphologically and genetically diverse and widespread in all types of marine habitats (Massana, 2015). In temperate seas, protist community undergoes strong seasonal changes in composition and abundance, as a re- sult of alterations in both abiotic factors,., irradiance, temperature and nutrient levels, and biotic factors, such as grazing, competition and interaction with pathogens and symbionts (Gran-Stadniczeñko, 2018).

Traditional studies of diversity, distribution and dyna- mics of marine protists were mainly based on microsco- pical identification and enumeration (Song., 2009, 2019; Liu, 2017, 2019; Pan., 2020). Due to the poor efficiency of the traditional methods, the protist diversity might have been underestimated (Keeling and Del Campo, 2017). The PCR-based approaches on the diversity of 18S rRNA gene have been widely applied in the last decade and revealed massive new lineages of aquatic protists (Moon-van der Staay, 2001; Lefranc, 2005). Recently, the high throughput sequencing techniques (HTS) has been applied to analyze the hyper-variable regions of the 18S rRNA gene, which has highlighted remarkable seasonal and spatial variations in marine protists (Gong, 2015; Sun, 2017; Li, 2019) and indi- cated that a large portion of taxa were still uncovered (De Vargas., 2015; Gran-Stadniczeñko, 2018). Com- pared to the microbiota in pelagic ocean, it is less known about the information of protist community in coastal areas (., estuary and intertidal zone), the involved biotic interactions, and how the seasonal andchemical- physical factors determine protist community structure and diversity. Coastal areas are potential hotspots of studies on microbial communities interacting with biogeoche- mical processes and adapting to environmental changes, such as salinity shifts and impacts of human activities (Balzano, 2015).

During the seasonal variation, autotrophic phytoplanktons can be affected by fluctuations of light, temperature, and/or nutrient concentration in the coastal area (Montagnes and Franklin, 2001; De la Vega., 2011; Vannier, 2016). As predators, heterotrophic protozoas are also indirectly infected by the change of chemical- physical factors (Paffenhöfer, 1998). In addition, the organismal interactions (., parasitism, symbiosis, and mixo- trophy) have recently been recognized as major biotic factors contributing to the protist community structure (Ca- ron, 2017; Stoecker and Lavrentyev, 2018). To our knowledge, neither abiotic nor biotic factors affecting pro- tist communities have been addressed in the coast of the northern Yellow Sea of China, although implications from the protist community dynamics may provide a glimpse to the inducements of many endemic ecological events,., red tide, jellyfish blooming, hypoxia, which are in-creasingly occurring in this coastal area (Wei., 2008; Wang, 2009; Xu, 2013).

In this study, using high-throughput sequencing (HTS) for the V4 region of 18S rDNA gene, year-round and monthly-taken samples were collected in two coastal sites around Yantai, a city along the northern Yellow Sea of China. Seasonal changes of planktonic protist community and the chemical-physical factors were investigated. Here we addressed the following questions: i) What is the seasonal distribution pattern of the protist community in the studied area? ii) How do the biotic interactions happen among protist groups? iii) What are the main drivers attributing to the seasonal pattern? By answering the above questions, we aim to better understand the dynamic and ecological function of planktonic protist communities in coastal ecosystems.

2 Materials and Methods

2.1 Sample Collection and Processing

Water samples were collected monthly from the surface layer (about 20cm deep) during the high tide at two sites (X and Y, Fig.1) in the coast of Yantai from March 2018 to February 2019. As shown in Fig.1, station X is located at the estuary of Xinan River (37˚25΄N, 121˚37΄E), while station Y is located at an intertidal zone near the Yantai University, about 8km away from station X (37˚28΄N, 121˚27΄E).

Fig.1 Locations of the sampling sites. X, estuary of Xinan River; Y, Yantai University coast.

The water temperature, salinity, DO and pH were mea- suredusing an electronic sensor (YSI, USA). Water samples for nutrient analysis were collected with Nis- kin bottles (15mL) and triplicates were set for each sample. The ammonium (NH4+), nitrate (NO3−), nitrite (NO2−), phosphate (PO43−) and silicate (SiO42−) were analyzed using a nutrient Autoanalyzer (Seal, Germany). The total organic carbon (TOC) and total organic nitrogen (TON) of water samples were determined with a Vario Micro Cube Elemental Analyzer (Elementar, Germany).

Protist communities were collected by filtrations of themarine water (500mL). To remove large plankton, the water samples were prefiltered through a 200μm mesh. Protist cells were then collected with a peristaltic pump (JINTENG, China) through polycarbonate filters with a diameter of 45mm a pore size of 0.22μm (Millipore, Ireland). The filter membranes with protistan cells were im- mediately put into cryopreservation tubes and stored at −80℃ before DNA extraction.

2.2 DNA Extraction and High-Throughput Sequencing

Extraction and purification of DNA from the filter mem- branes were carried out using the Fast DNA SPIN Kit for Soil (Q-BIO gene, USA). The DNA concentrations were determined using a Nanodrop 2000c spectrophotometer (Thermo-Fisher, USA).

The V4 region of 18S rDNA has been successfully used to deduce the seasonal dynamics of the Skagerrak plankton communities (Gran-Stadniczenko, 2019). Accordingly, in the present research the V4 region was applied to study the protist community. PCR amplification was performed using the primer set 528F (5’-GCGGTA ATTCCAGCTCCAA-3’) and 706R (5’-AATCCRAGAA TTTCACCTCT-3’). Sequencing was performed bidirectionally on the Ion S5TMXL platform (Thermofisher,USA) in Novogene sequencing company (Tianjin, China).

2.3 Processing Analyses of High-Throughput Sequencing (HTS) Data

Raw data were processed and analyzed using QIIME v.1.8.0 (Caporaso, 2010) and Mothur v.1.34.4 (Schloss, 2009). Sequences with high quality were obtained after processes of quality controlling, denoising and discarding chimera and singleton.

The clean sequences were clustered to operational taxo- nomic units (OTU) at a sequence similarity of 97%. The OTUs were blasted, aligned, and classified against theProtist Ribosomal Reference (PR2 database; Guillou., 2012). Sequences of non protists (metazoan, macroalgae) and unclassified taxa were omitted (2% of OTUs). Alpha diversity of protists inferred by OTU richness, Shannon, Simpson, and Chao1 were calculated based on the lowest number of the high-quality sequences among all samples. Beta diversity was calculated with Bray-Curtis distances and was visualized using nonmetric multidimensional sca- ling (nMDS) in PRIMER v.6.0 (Primer-E, UK).

2.4 Network Construction

The co-occurrence network was inferred based on the spearman correlation matrix constructed following the methods of Ma(2016), using the WGCNA package (Langfelder and Horvath, 2012) appended in R v3.5.1 (R Core, 2016) and packages provided by Ma(2016). Only OTUs with relative abundance higher than 1% were subjected to the analyses, and the positive correlations be- tween OTUs were considered. The-values for multiple testing were adjusted using the Benjamini and Hochberg false discovery rate (FDR) controlling procedure (Benjamini, 2006). The direct correlation dependencies were distinguished using the network deconvolution method (Feizi, 2013). The co-occurrence networks were con- structed based on correlation coefficients and FDR adjusted-values for correlation. The cutoff of FDR-adjust- ed-values was set as 0.05. The cutoff of correlation coefficients was 1.25 for station X and 1.4 for station Y through the random matrix theory-based methods (Luo, 2006).

2.5 Statistical Analyses

One-way ANOVA analysis was performed to examine the differences in environmental factors, α-diversity estimators and proportions of abundant (>1%) groups at family level among four seasons. Spearman’s correlation ana- lysis was conducted to examine the relationships between α-diversity estimators or relative abundance of protistan families and environmental factors. The normality of data was tested and log-transformed when it was necessary. All these analyses were executed using SPSS v. 20.0 (SPSS,Chicago, IL, USA, Kirkpatrick and Feeney, 2012). To sta- tistically test the variations in community structure of planktonic protists in different seasons, ANOSIM (analysis of similarity) was performed using PRIMER V.6 (Clarke and Gorley, 2006). The major environmental factors shap- ing the community structure of planktonic protists were determined by conducting Redundancy analysis (RDA) using packages in R v3.5.1 (R Core, 2016).

3 Result

3.1 Seasonal Variations of Environmental Factors

Generally, physicochemical properties of marine water (temperature, DO and nutrients) showed clear seasonal patterns (Table 1). The water temperature increased up to 25.4℃ in summer and decreased to a minimum of 2.5℃ in winter (ANOVA,<0.01). Ammonium, phosphate, silica and nitrate exhibited the contrast trend to temperature, with the highest values in winter and lowest values in sum- mer (<0.02). The concentration of dissolved oxygen (DO) also changed with season significantly. The highest DO concentration was found in February (14.20mgL−1) and the lowest value was recorded in August and September (5.52mgL−1) (=0.03). Other chemical factors (., the total organic carbon and nitrogen, C:N ratio, and N:P ratio) showed no significant difference with seasonal variation (>0.05). No significant regional difference be- tween X and Y stations in environmental factors were de- tected (Table 1).

3.2 Overview of Sequencing Data and Community Composition

After quality filtering, a total of 1745616 reads were obtained in the 24 samples. Further on, 1231976 protistan sequences were retained after removing reads classified as metazoan and multicellular organisms. Sequence num- bers of the 24 samples varied from 8680 to 76258 (Table 2).

By comparing the relative abundance of super-phylum groups between two sites among different seasons, we found that there was no statistically significant difference between two sites (T-test,>0.05; Table 3); however, most groups showed significant difference among seasons (ANOVA,<0.05). Therefore, we combined the data of two sampling sites and focused on the seasonal variation of the protozoan community.

Table 1 Seasonal and regional differences of the environmental factors

Notes: C:N means the ratio of the total organic carbon to the total organic nitrogen; N:P means the ratio of dissolved inorganic nitrogen (the sum of NH4+, NO3−and NO2−) to PO43−. Significant differences (<0.05) are highlighted in bold.

Table 2 Summary of the pyrotags and alpha diversity of planktonic protists in all samples

Table 3 Statistical analyses of the difference of protist communities at super-phylum level between two sites and among four seasons

Note: Significant differences (<0.05) are highlighted in bold.

Generally, infrakingdom Alveolata dominated in the protist communities, accounting for 54% of the total protistan sequences, and was observed in all samples (Fig.2A). Most of the Alveolata sequences were clustered into phylum Dinoflagellata (74%) and Ciliophora (20%) (Fig.2B). Higher abundance of Dinoflagellata (40%–80%) was ob- served from April to July (spring and summer), while the lowest abundances appeared in autumn (mean 22%). Ci- liophora was more abundant in winter than in other seasons (Fig.2B). Within Dinoflagellata, the most abundant classes were Dinophyceae (50%) and Syndiniales (47%) (Fig.2C); the most abundant orders were Dino-Group-I (attributing to 34% of Dinoflagellata), Peridiniales (20%) and Dinophyceae_X (16%) (Fig.2D). In Ciliophora, Spirotrichea was the most representative class (Fig.2C), in which Tintinnida (31%), Strombidiida (32%) and Choreo- trichida (11%) were dominant (Fig.2D).

Super-group Archaeplastida (22.85%±3.81% of the to- tal sequences) was the second abundant group in the protist community (Fig.2A), which was primarily represented by the phylum Chlorophyta (95%; Fig.2B). Mamiellophy- ceae (45%), Trebouxiophyceae (25%) and Ulvophyceae (16%) were the most abundant classes in Chlorophyta(Fig.2C). At the order level, Mamiellales (97%) was do- minant in class Mamiellophyceae, while class Trebouxiophyceae was dominated by order Chlorellales (99%; Fig.2D). Notably, the order Ulotrichales of Ulvophyceae was markedly dominant in March, which accounted for almost 45% of the sequences of the month (Fig.2D).

Subkingdom Stramenopiles accounted for 12.57% of the total reads, with Ochrophyta as the major phylum (65%), followed by Sagenista (16%) and Pseudofungi (14%; Figs. 2A, B). Class Bacillariophyta was dominant in phylum Ochrophyta (70%). Labyrinthulomycetes (86%) was the most abundant class in Sagenista, followed by MAST-6 (6.2%), MAST-9 (3.9%) and MAST-7 (2.1%) (Fig.2C). Oomycota (50%) and MAST-1 (40%) were the most abundant classes in Pseudofungi (Fig.2C). At the order level, the most abundant groups were Bacillariophyta (46%), Chrysophyceae_X (8.5%) and Thraustochytriales (2.9%; Fig.2D).

In addition, super-phylum Opisthokonta (4.9%±0.89%), Rhizaria (2.6%±0.32%), and Hacrobia (2.3%±0.56%) exhibited relatively low abundance. Other super-phyla (., Ameobozoa, Apusozoa and Excavata) appeared to be mi- nor ones (less than 1%; Fig.2A).

Fig.2 The relative abundance of taxonomic groups of protists at superphylum (A), phylum (B), class (C), and order (D) levels.

3.3 Alpha Diversity of Planktonic Protists

Alpha-diversity estimators of protists varied greatly among the monthly samples, with OTU richness ranging from 548 to 1499, Simpson indices from 0.296 to 0.981, Shannon indices from 1.90 to 7.35, and Chao1 indices from 1348 to 3487 (Table 2). Simpson and Shannon in- dicators were significantly different among four seasons (<0.05, Fig.3), while the highest value was found in autumn and the lowest value was recorded in spring.Otherwise, Chao1 and Richness indices showed no significant difference among seasons (>0.05, Fig.3).

Fig.3 Comparisons of alpha diversity estimators of planktonic protists. P values are given for the comparisons among seasons using one-way ANOVA.

3.4 Seasonal Distribution of Planktonic Protists

Seasonality of protist community was indicated by the ordination analyses based on the Bray-Curtis dissimilarities. In the nMDS plot (Fig.4), summer, autumn and win- ter samples were clearly separated from each other, while spring samples were evenly distributed between the spaces occupied by other season samples. The result of ANOSIM confirmed that samples among four seasons were significantly distinct (<0.05, Table 4).

Fig.4 Nonmetric multidimensional scaling (nMDS) plot showing the difference in protist community structures collected from different seasons. The stress value of 0.018 indicates a good ordination.

Table 4 ANOSIM statistical analyses of the difference of protist communities among seasons

Note: Significant differences (<0.05) are highlighted in bold.

The heterogeneity of protist community among four sea- sons was further determined by a comparison of relative abundance of the major (relative abundance>1%) groups at family level (Fig.5). Significantly, the relative abun- dances of Chlorellales_X and Suessiaceae were the high- est in summer and the lowest in winter (<0.05). In con- trary, Bathycoccaceae showed higher ratios in autumn(11.6%) and winter (11.4%) than in spring (0.9%) and sum- mer (0.7%) (<0.01). Similarly, Mamiellaceae, Ceratia- ceae, Cryptomonadales_X, Strombidiidae_R, Tintinnida_ TIN_9, Raphid-pennate, and Dino-group-II-Clade-1 had the highest proportions in autumn or winter (<0.05). In spiteof no statistical significance, the abundant families Dino- Group-I-Clade-1 and Ulotrichales_X had obviously high- er proportions in spring and summer.

3.5 Co-Occurrence Network Among Protistan Groups

At both sites, significantly high betweenness centrality scores were observed for OTUs affiliated to two prevailing classes, the ciliate Spirotrichea and the dinoflagellate Syndiniale (Figs.6, 7). The involved Spirotrichea comprised mainly tintinnids (,,, and) and a few species from Hypotrichia () and Choreotrichida (). The top frequent co-occurred Syndiniales belonged to Di- no-Group-I, with few exceptions from Dino-Group-II. Co- occurrences were also observed among classes Spirotri- chea (mainly Strombidiida), Mamiellophycease (genera,, and), and/or Di- nophyceae (.,) (Figs.6, 7). Ad- ditionally, OTUs belonging to the same class Dinophyceae (blue circles) significantly co-occurred at station X (Fig.6). These OTUs affiliated to genera,and one unidentified Dinophyceae species. Class Bacil- lariophyta (represented by,and one unidentified Polarcentric-Mediophyceae), one speciesfrom Dictyochophycea, and seve- ral OTUs from Mast-groups and Chrysophyceae were po- sitively co-varied (Fig.6). Additionally, OTUs affiliated to genusin class Trebouxiophycease, showed po- sitive co-variations to two Dinophyceae OTUs, affiliated to generaand(Fig.6). At station Y(Fig.7), we specifically observed significantly high be- tweenness centrality scores among groups Labyrinthulomy- cetes (green circles) and two OTUs from Bacillariophyta () and Ulotrichales () (Fig.7).

Fig.5 Comparisons of relative abundance of the major taxa (relative abundance>1%) of microeukaryotes at family level among four seasons. Asterisks and double asterisks indicate P values<0.05 and <0.01, respectively.

Fig.6 The co-occurrence network interactions of main groups (at OTU level, relative abundance>1% in each sample) in station X. The nodes represented unique OTU in the data sets. OTUs of the top eight abundant protistan classes were labeled in color, while those of the rest classes were in gray. The connection stands for a strong (Spearman’s ρ>1.25) and significant correlation (P<0.05).

3.6 Environmental Drivers of Variation in Protist Diversity and Community Structure

Spearman’s correlation coefficients between alpha di- versity estimators and environmental factors were given in Table 5. Briefly, Shannon index had the strongest cor- relation with SiO42−(=0.56,=0.01), followed by PO43−(=0.54,=0.01), TOC (=0.52,=0.01), NO2−(=0.49,=0.02) and pH (=0.41,=0.05), respectively. Simpson index had the strongest correlation with the C:N ratio (=0.56,=0.01), followed by pH (=0.50,=0.01), PO43−(=0.50,=0.01), TOC (=0.46,=0.03), and SiO42−(=0.46,=0.01), respectively. No signifi- cant correlation between OTU richness and Chao1 index with environmental variables was observed, as shown in Table 5.

Table 5 Spearman’s correlation coefficients (R) between alpha diversity estimators and environmental factors

Notes: **≤0.01; *≤0.05.

Fig.7 The co-occurrence network interactions of main groups (at OTU level, relative abundance>1% in each sample) in station Y. The nodes represented unique OTU in the data sets. OTUs of the top eight abundant protistan classes were labeled in color, while those of the rest classes were in gray. The connection stands for a strong (Spearman’s ρ>1.4) and significant correlation (P<0.05).

The redundancy analyses (RDA) plots (Fig.8) suggest- ed that temperature, PO43−and TON were the most sig- nificant factors influencing the whole protist community structure among seasons (<0.05).

Furthermore, the associations between relative propor- tion of major families (relative abundance>1%) and en- vironmental factors were also explored (Table 6). The abun- dant families Dino-Group-I (Clade-1 and X) and Dinophy- cease_XX did not show significant relationship with any en- vironmental factors (>0.05). Heterocapsaceae only had significant and positive correlation with temperature (=0.43,<0.05). In contrast, Chlorellales_X was positively related to temperature (=0.88,<0.01), but negatively related to DO (=−0.56,<0.01) and nutrients including NO3−(=−0.78,<0.01), NH4+(=−0.71,<0.01), PO43−(=−0.46,<0.05) and SiO42−(=−0.42,<0.05). Si- milarly, Suessiaceae correlated positively with tempera- ture (=0.61,<0.01), and negatively with DO (=−0.60,<0.01), NO3−(=−0.42,<0.05) and NH4+(=−0.44,<0.05). The proportion of Bathycoccaceae primarily ex- hibited significant positive relationship with nutrientsPO43−, SiO42−and NO3−(<0.05). Most of the other groups,., Ulotrichales_X, Polar-centric-Mediophyceae, Cryptomo- nadales_X, Tintinnida_TIN_09, and Dino-Group-II-Clade- 1, presented a negative correlation with temperature and a positive correlation with nutrients (<0.05).

Fig.8 The redundancy analyses (RDA) plots showing relationships between seasonal communities and environmental factors.

Table 6 Spearman’s correlation coefficients (r) between the major taxa (relative abundance>1%) of microeukaryotes at family level and environmental factors

Notes:**≤0.01;*≤0.05, and insignificant correlations were left blank.

4 Discussion

Previously, the protist diversity and community com- position in coastal water and sediments in the Yellow Sea were mainly studied using morphology-based methods (Song, 2009; Zhang and Xu, 2015). In this study, we firstly investigated the seasonal dynamic of protist com-munity in the coastal Yellow Sea based on the HTS tech- nology, which will broaden our insight to the dynamic and environmental response of protists in coastal ecosys- tems. Using HTS technology, a total of 26038 OTUs from different protistan groups were revealed, which far ex- ceeded the number of taxa detected through morphologi- cal observation.

4.1 Dominant Groups in the Surficial Water of the Yantai Coast, Northern Yellow Sea

Alveolata was revealed as the most dominant infrak- ingdoms in this study, which is comparable with the clone library analyses and HTS studies on microeukaryotes in other water columns (Pernice, 2013; Gran-Stadnic- zeñko, 2018) and the surficial sediment in the north Yellow Sea (Gong, 2015).

Dinoflagellata was the most abundant phylum in Al- veolata, dominated by classes Syndiniales Dino-Group-I and Peridiniales. Syndiniales is a group comprised exclu- sively of parasitic species that infect tintinnid ciliates, crustaceans, dinoflagellates and fish (Guillou, 2008). The high abundance of Syndiniales detected in this study corresponds well to the previous findings in several me-tabarcoding studies on planktonic microeukaryotes (Geor- ges, 2014; Christaki, 2017), and supports thehypothesis that parasitic interactions could dominate plank- ton networks (De Vargas, 2015; Lima-Mendez, 2015). Otherwise, it should be aware that potentially high rRNA gene copy numbers and inactive spores could lead to overestimation on the proportion of Syndiniales (Not, 2009; Sassenhagen, 2019). Another dominant order Peridiniales is mainly composed of small mixotro- phic species (<200μm), harboring a number of red tide- forming taxa(Iwataki, 2008). In our samples, one OTU related to speciesdominated the Peridiniales (98%). Consistently,was also recorded as one significantly abundant OTU in the surfi- cial sediments of Yellow Sea coast (Gong, 2015). These findings suggest a common presence of red tide- forming species in the coast of the Yellow Sea.

Archaeplastida was the second dominant super-group, mainly comprised of phylum Chlorophyta. The propor- tion of Chlorophyta detected in the Yellow Sea samples is much higher than that revealed in other water samples by HTS analyses (Pernice, 2013; Gran-Stadniczeñko, 2018),., in the surface water of estuarine in Skagerrak coast of Europe (22%. 7.6%) (Gran-Stadniczeñko, 2018).

Lower proportion of subkingdom Stramenopiles was detected in our Yellow Sea coastal samples, compared to those reported previously in the estuarine samples in the Skagerrak coast (12%. 20%) (Gran-Stadniczeñko, 2018), and in the surficial sediments (12%. 15%) of the Yellow Sea (Gong, 2015). Moreover, our study re- vealed less proportion of Oomycota (0.8%. 3.6%) and Labyrinthulea (0.6%. 3.4%) in the surface water com- pared to that recorded in the sediment (Gong, 2015). This difference between surface water and sediment sam- ples might attribute to the distinct ecological niches. La- byrinthulomycetes was reported playing a key role on the initial stage of the microbial chain food as an organicmatter degrader (Raghukumar, 2002), while Oomycota was predicted being saprobic or parasitic on terrestrial or aqua- tic plants and animals (Nigrelli and Thines, 2013).

4.2 Co-Occurrence of Protist Groups

Significant correlations between OTUs in dinoflagel- late Syndiniales and ciliate Spirotrichea were indicated by the co-occurrence analysis (Figs.6, 7). The co-occurrence might reflect a parasite-host relationship. A few genera of Syndiniales Dino-Group-I were reported as common pa- rasites in tintinnids of Spirotrichea (Coats, 2010;Bachvaroff., 2012). An increase of tintinnid hosts’ population might provide breeding ground for the active parasite or dinospores of syndiniales (Coats, 2010). Other spirotrichs (., Hypotrichia and Choreotrichida) involved in this relationship may indicate a wider range of host-parasite relationships than previously known, which need furtherobservations to validate theinferences (Figueroa, 2010). Other biotic interac- tions, such as the parasites favoring the growth of one species through the infection of its grazers or competitors, may also contribute to this ciliate-parasite co-occurrence (Coats and Heisler, 1989).

The co-occurrence of ciliates and Syndiniales in this study is different from the result of Sassenhagen.(2019) on Kerguelen Islands, which showed negative re- lationship between two tintinnid genera (and) and Syndiniales Group-II-clade-12. The contradiction might attribute to the variable efficiencies that Syndiniales infect ciliates. The life cycles of the Syn- diniales parasitic on ciliates encompass a bi-flagellated in- fective stage, and the non-infectious dinospore stage that Syndiniales may only penetrate the host cell membrane or be ingested by the host as preys (Coats, 2010; Je- phcott, 2016). Indeed, Sassenhagen(2019) observed high concentrations of free Syndiniales spores but relatively few amounts of Syndiniales associated with microplankton.

Remarkable co-occurrence was observed among ciliate Spirotrichea (mainly Strombidiida), Mamiellophyceae (and), and Dinophyceae (.,). This could represent a predation relationship considering the latter two small phytoplank- ton are well known food resources for phagocytic ciliates (Lynn, 2008). On another hand, the co-variation might associate with the mixotrophic nature of organisms in the three groups. Strombidiids are typical plastidic ciliates ca- pable of mixotrophic by acquiring phototrophy (Stoecker and Lavrentyev, 2018). Other ciliate species involved in this network,.,sp. andsp. (Figs.6, 7), are also known mixotrophs with plastid-re- tention (Stoecker and Lavrentyev, 2018). Meanwhile, the photosynthetic algaeandare well known mixotrophic cases which sometimes are more bac- terivory than heterotrophic flagellates (Jeong, 2005; Stoecker and Lavrentyev, 2018). In this sense, the co- occurrence of Spirotrichea, Mamiellophyceae, and Dino- phyceae in this study could also be a result of homolo- gous ecological niches with similar requirement of nutri- tion and environmental conditions (., light, bacteria food).

4.3 Seasonal Variation of the Protist Community and the Environmental Drivers

Our analyses revealed a significant seasonal pattern of protist community structure (Figs.2, 5). The Shannon and Simpson index also showed significant seasonal fluctua- tion (Fig.3). The highest micro-eukaryotic diversity oc- curred in autumn, while the lowest one appeared in spring (Fig.3), which is different to the commonly accepted no- tion that phytoplankton generally produce maxima in spring and autumn, minima in summer and winter (Wet- zel, 1983). This divergence could be interpreted by thespecific features in community structures and nutrientfluctuations in the studied coast of the Yellow Sea. The lowest protistan diversity in spring might attribute to the significant blooms of several families in this season, in- cluding two parasitic families of Syndiniales Dino-group- I and another family Ulotrichales_X (probably in forms of spore with diameter<200μm, or retentions from large foliate bodies) (Fig.5). On another hand, the protist diver- sity strongly positively correlated with nutrients (Table 6), which may explain the highest diversity of these groups in autumn. Indeed, we observed higher nutrient concen- trations in autumn of the year 2018 in this area (Table 1).

The family Ulotrichales_X blooming in spring was do- minantly represented by one green algain our data (Fig.5), which has been reported usually occurring in cold-water sublittoral in winter and spring (Sumner and Osburn, 1913; Levring., 2019). This is consistent to its negative correlation to temperature in the Spearman’s analysis (Table 6). Moreover, the positive Spearman’s correlation of Ulotrichales_X with NH4+in- dicates a stimulation by ammonia nitrogen to its spring blooming.

Bathycoccaceae in order Mamiellales of phylum Chlo- rophyta was significantly higher in autumn and winter than in spring and summer (Fig.5). This finding agrees with several other studies showing that members of Ba- thycoccaceae have advantages under low-light conditions or persist throughout winter darkness in marine water (Joli, 2017). In this study, the distribution of this family could be mainly explained by the changes of SiO42−and NO3−in water, but had no correlation with tempe- rature (Table 6). The family was represented by genusin our data (42%). The rich nutrients in the winter water may favor the growth of, which has shown high abundance and wide distribution in nu- trient-rich waters but relies less on the temperature and light time (Vannier, 2016).

The family Strombidiidae_R of phylum Ciliophora show- ed significantly higher abundance in winter than in other seasons (Fig.5). The prevailing of Strombidiidae over other heterotrophic microzooplankton in winter may at- tribute to their advantages in mixotrophy (Stoecker and Lavrentyev, 2018). These mixotrophic ciliates are less re- stricted by the decreased food resources of small phyto- plankton in winter, compared to the strictly heterotrophic microzooplankton (Mitra, 2016).

5 Conclusions

In this study, we obtained a total of 26038 protistan OTUs in the northern Yellow Sea coastal samples using HTS. Alveolata, Archaeplastida, and Stramenopiles are the dominant super-phylums in this coastal area. The protist community composition and relative abundance showed dramatic seasonal variation through the year. However,we did not observe significant difference between the estu- ary and intertidal zones, which may be caused by the close distance of the two sampling sites, and should be considered in the future sampling strategies. Blooms of the chlorophyte family Ulotrichales_X and the red tide- forming taxawere recorded in some seasons, which may indicate a sign of eutrophica- tion in the coast of Yantai. The co-occurrence analysis method applied in this work turned out to be efficiencyand revealed frequent and strong biotic interactions, which could be widely used in the future ecological studies on protists.

Acknowledgements

This work was supported by the National Natural Sci- ence Foundation of China (Nos. 31672251, 31772413), the Youth Innovation Promotion Association, CAS (No. 2019216), the Strategic Priority Research Program of the Chinese Academy of Sciences (No. XDA23050303), and the Key Research Project of Frontier Science, CAS (No. QYZDBSSW-DQC013-1). We specially thank Prof. Bin Ma from Zhejiang University for his kindly help in the co-occurrence network analyses.

Bachvaroff, T. R., Kim, S., Guillou, L., Delwiche, C. F., and Coats, D. W., 2012. Molecular diversity of the syndinean ge- nusbased on single-cell PCR analysis., 78: 334-345.

Balzano, S., Abs, E. C., and Leterme, S., 2015. Protist diversity along a salinity gradient in a coastal lagoon., 74: 263-277.

Benjamini, Y., Krieger, A. M., and Yekutieli, D., 2006. Adaptive linear step-up procedures that control the false discovery rate., 93: 491-507.

Caporaso, J. G., Kuczynski, J., Stombaugh, J., Bittinger, K., Bushman, F. D., Costello, E. K., and Fierer, N., 2010. QIIME allows analysis of high-throughput community sequencing data., 7: 335-336.

Caron, D. A., Alexander, H., Allen, A. E., Archibald, J. M., Armbrust, E. V., Bachy, C., and Heidelberg, K. B., 2017. Pro- bing the evolution, ecology and physiology of marine protists using transcriptomics., 15 (1): 6.

Christaki, U., Genitsaris, S., Monchy, S., Li, L. L., Rachik, S., and Breton, E., 2017. Parasitic eukaryotes in a meso-eutro- phic coastal system with markedblooms., 4: 416.

Clarke, K. R., and Gorley, R. N., 2006.. Version 6.0. Plymouth Marine Laboratory, Ply- mouth, UK, 93pp.

Coats, D. W., and Heisler, J. J., 1989. Spatial and temporal occurrence of the parasitic dinoflagellateand its tintinnine hostin Chesapeake Bay., 101 (3): 401-409.

Coats, D. W., Kim, S., Bachvaroff, T. R., Handy, S., and Del- wiche, C. F., 2010.n. gen., n. sp. (Phy- lum Dinoflagellata), an ectoparasite of the ciliateDaday 1887, and its relationship toGrassé 1952., 57: 468-482.

De la Vega, M., Diaz, E., Vila, M., and León, R., 2011. Isola- tion of a new strain ofspand characterization of its potential biotechnological applications., 27 (6): 1535-1543.

De Vargas, C., Audic, S., Henry, N., Decelle, J., Mahé, F., Lo- gares, R., and Carmichael, M., 2015. Eukaryotic plankton di- versity in the sunlit ocean., 348 (6237): 1261605.

Falkowski, P. G., and Raven, J. A., 2007.. Princeton University Press, Princeton, 484pp.

Feizi, S., Marbach, D., Médard, M., and Kellis, M., 2013. Net- work deconvolution as a general method to distinguish direct dependencies in networks., 31: 726-733.

Figueroa, R. I., Garcés, E., and Camp, J., 2010. Reproductive plasticity and local adaptation in the host-parasite system form- ed by the toxicand the dinoflagellate parasite., 10 (1): 56-63.

Georges, C., Monchy, S., Genitsaris, S., and Christaki, U., 2014. Protist community composition during early phytoplankton blooms in the naturally iron-fertilized Kerguelen area (South- ern Ocean)., 11: 5847-5863.

Gong, J., Shi, F., Ma, B., Dong, J., Pachiadaki P., and Zhang, X. L., 2015. Depth shapes α- and β-diversities of microbial eu- karyotes in surficial sediments of coastal ecosystems., 17 (10): 3722-3737.

Gran-Stadniczeñko, S., Egge, E., Hostyeva, V., Logares, R., Wen- che, E., and Edvardsen, B., 2018. Protist diversity and sea- sonal dynamics in Skagerrak plankton communities as reveal- ed by metabarcoding and microscopy., 66: 494-513, DOI: 10.1111/jeu.12700.

Guillou, L., Bachar, D., Audic, S., Bass, D., Berney, C., Bittner, L., and Del Campo, J., 2012. The Protist Ribosomal Refer- ence database (PR2): A catalog of unicellular eukaryote small subunit rRNA sequences with curated taxonomy., 41 (D1): D597-D604.

Guillou, L., Viprey, M., Chambouvet, A., Welsh, R. M., Kirk- ham, A. R., Massana, R., and Worden, A. Z., 2008. Widespread occurrence and genetic diversity of marine parasitoids be- longing to Syndiniales (Alveolata)., 10 (12): 3349-3365.

Iwataki, M., 2008. Taxonomy and identification of the armored dinoflagellate genus(Peridiniales, Dinophyceae), 3 (3): 135-142.

Jeong, H. J., Doo, D., Park, J. Y., Song, J. Y., Kim, S. T., and Lee, S. H., 2005. Feeding by phototrophic red-tide dinoflag- ellates: Five species newly revealed and six species previ- ously known to be mixotrophic., 40 (2): 133-150, DOI: 10.3354/ame040133.

Jephcott, T. G., Alves-de-Souza, C., Gleason, F. H., Van Ogtrop, F. F., Sime-Ngando, T., Karpov, S. A., and Guillou, L., 2016. Ecological impacts of parasitic chytrids, syndiniales and per- kinsids on populations of marine photosynthetic dinoflagel- lates., 19: 47-58.

Joli, N., Monier, A., Logares, R., and Lovejoy, C., 2017. Sea- sonal patterns in Arctic prasinophytes and inferred ecology ofunveiled in an Arctic winter metagenome., 11 (6): 1372-1385.

Keeling, P. J., and Del Campo, J., 2017. Marine protists are not just big bacteria., 27 (11): R541-R549.

Kirkpatrick, L., and Feeney, B., 2012.Nelson education, Scarborough, 128pp.

Langfelder, P., and Horvath, S., 2012. Fast R functions for robust correlations and hierarchical clustering.,46 (11): i11.

Lefranc, M., Thenot, A., Lepère, C., and Debroas, D., 2005. Ge- netic diversity of small eukaryotes in lakes differing by their trophic status., 71: 5935- 5942.

Levring, T., Hoppe, H. A., and Schmid, O. J., 2019.. Walter De Gruyter, Ber- lin, 423pp.

Li, X. H., Huang, J., Filker, S., Stoeck, T., Bi, Y. H., and Song, W. B., 2019. Spatio-temporal patterns of zooplankton in a mainstem dam affected tributary: A case study in the Xiangxi River of the Three Gorges Reservoir., 62: 1058-1069, DOI: 10.1007/s11427-018-9523-0.

Lima-Mendez, G., Faust, K., Henry, N., Decelle, J., Colin, S., Carcillo, F., and Bittner, L., 2015. Determinants of commu- nity structure in the global plankton interactome., 348 (6237): 1262073.

Liu, W. W., Jiang, J., Xu, Y., Pan, X., Qu, Z., Luo, X., El- Serehy, H. A., Warren, A., Ma, H., and Pan, H., 2017. Di- versity of free-living marine ciliates (Alveolata, Ciliophora): Faunal studies in coastal waters of China during the years 2011–2016.,61: 424-438.

Liu, W. W., Zhang, K. X., Chen, C., Li, J. Q., Tan, Y. H., War- ren, A., Lin, X. F., and Song, W. B., 2019. Overview of the biodiversity and geographic distribution of aloricate oligo- trich ciliates (Protozoa, Ciliophora, Spirotrichea) in coastal waters of southern China., 17 (8): 787-800, DOI: 10.1080/14772000.2019.1691081.

Luo, F., Zhong, J. X., Yang, Y. F., Scheuermann, R. H., and Zhou, J. Z., 2006. Application of random matrix theory to biological networks., 357: 420-423.

Lynn, D. H., 2008.. Springer Verlag, Germany, 604pp.

Ma, B., Wang, H. Z., Dsouza, M., Lou, J., He, Y., Dai, Z. M., Brookes, P. C., Xu, J. M., and Gilbert, J. A., 2016. Geographic patterns of co-occurrence network topological features for soil microbiota at continental scale in eastern China., 10 (8): 1891-1901.

Massana, R., 2015. Protistan diversity in environmental mole- cular surveys. In:Ohtsuka, S.,., eds., Springer, Tokyo, 3-21.

Mitra, A., Flynn, K. J., Tillmann, U., Raven, J., and Caron, D., 2016. Defining planktonic protist functional groups on mecha- nisms for energy and nutrient acquisition: Incorporation of diverse mixotrophic strategies., 167: 106-120.

Montagnes, D. J., and Franklin, M., 2001. Effect of temperature on diatom volume, growth rate, and carbon and nitrogen content: Reconsidering some paradigms., 46 (8): 2008-2018.

Moon-van der Staay, S. Y., De Wachter, R., and Vaulot, D., 2001. Oceanic 18S rDNA sequences from picoplankton reveal un- suspected eukaryotic diversity., 409 (6820): 607-610.

Nigrelli, L., and Thines, M., 2013. Tropical oomycetes in the German Bight–Climate warming or overlooked diversity?, 6 (2): 152-160.

Not, F., Del Campo, J., Balagué, V., De Vargas, C., and Mas- sana, R., 2009. New insights into the diversity of marine pi- coeukaryotes., 4 (9): e7143.

Paffenhöfer, G. A., 1998. Heterotrophic protozoa and small me- tazoa: Feeding rates and prey-consumer interactions., 20 (1): 121-133.

Pan, H. B., Zhang, Q. Q., Dong, J. Y., and Jiang, J. M., 2020. Morphology and phylogeny of two novel Pleurostomatids (Ciliophora, Litostomatea), with establishing a new genus., 67 (2): 252-262, DOI: 10.1111/ jeu.12779.

Pernice, M. C., Logares, R., Guillou, L., and Massana, R., 2013. General patterns of diversity in major marine microeukaryote lineages., 8 (2): e57170.

R Core, T., 2016.. R Foundation for Statistical Computing, Vi- enna, 3165pp.

Raghukumar, S., 2002. Ecology of the marine protists, the La- byrinthulomycetes (Thraustochytrids and Labyrinthulids)., 38 (2): 127-145.

Sassenhagen, I., Irion, S., Jardillier, L., Moreira, D., and Christaki, U., 2019. Protist interactions and community structure during early autumn in the Kerguelen region (Southern Ocean)., 171 (1): 125709.

Schloss, P. D., Westcott, S. L., Ryabin, T., Hall, J. R., Hartmann, M., and Hollister, E. B., 2009. Introducing mothur: Open- source, platform-independent, community supported software for describing and comparing microbial communities., 75 (23): 7537-7541.

Sherr, E. B., and Sherr, B. F., 2002. Significance of predation by protists in aquatic microbial food webs., 81 (1-4): 293-308.

Song, W. B., Warren, A., and Hu, X. Z., 2009.. Science Press, Bei- jing, 518pp.

Song, W., Xu, D. P., Zhang, Q. Q., Liu, W. W., Warren, A., and Song, W. B., 2019. Taxonomy and phylogeny of two poorly studied genera of marine oligotrich ciliates including descrip- tions of two new species:sp. n. andsp. n. (Ciliophora: Spirotri- chea).,70: 1-16.

Stoecker, D. K., and Laventyev, P. J., 2018. Mixotrophic plank- ton in the polar seas: A pan-arctic review., 5: 292.

Sumner, F. B., and Osburn, R. C., 1913.. Bulletin of the Bureau of Fisheries, United States, 545pp.

Sun, P., Huang, L. Y., Xu, D. P., Huang, B. Q., Chen, N. W., and Warren, A., 2017. Marked seasonality and high spatial variation in estuarine ciliates are driven by exchanges be- tween the ‘abundant’ and ‘intermediate’ biospheres., 7 (1): 1-12.

Vannier, T., Leconte, J., Seeleuthner, Y., Mondy, S., Pelletier, E., Aury, J. M., and Wincker, P., 2016. Survey of the green picoalgagenomes in the global ocean., 6: 37900.

Wang, B. D., 2009. Hydromorphological mechanisms leading to hypoxia off the Changjiang Estuary., 67 (1): 53-58.

Wei, G. F., Tang, D. L., and Wang, S. F., 2008. Distribution of chlorophyll and harmful algal blooms (HABs): A review on space-based studies in the coastal environments of Chinese marginal seas., 41 (1): 12-19.

Wetzel, R. G., 1983. Attached algal-substrata interactions: Fact or myth, and when and how? In:. Springer, Dordrecht, 207-215.

Worden, A. Z., Follows, M. J., Giovannoni, S. J., Wilken, S., Zim- merman, A. E., and Keeling, P. J., 2015. Rethinking the marine carbon cycle: Factoring in the multifarious lifestyles of microbes., 347 (6223): 1257594, DOI: 10.1126/sci- ence.1257594.

Xu, Y. J., Ishizaka, J., Yamaguchi, H., Siswanto, E., and Wang, S. Q., 2013. Relationships of interannual variability in SST and phytoplankton blooms with giant jellyfish () outbreaks in the Yellow Sea and East China Sea., 69 (5): 511-526.

Zhang, W., and Xu, H. L., 2015. Seasonal shift in community pattern of periphytic ciliates and its environmental drivers in coastal waters of the Yellow Sea, northern China., 95 (2): 277-288.

. E-mail: qqzhang@yic.ac.cn

E-mail: daodeji@126.com

January 18, 2020;

March 17, 2020;

May 8, 2020

(Edited by Qiu Yantao)