1 INTRODUCTION

Microeukaryotic plankton, including protozoans, algae, fungi, and metazoans, are key components of the marine food webs (Anjusha et al., 2013). Eukaryotic plankton includes numerous species of different body sizes with different ecological functions. Small-sized organisms are supposed to have low sinking rate, less dispersal ability, and higher population growth rate (Raven, 1998; Shurin et al., 2009; Korhonen et al., 2010). Pico- and nanoplankton (0.2–2 μm and 2–20 μm in diameter, respectively) are important parts of marine ecosystem (Piwosz et al., 2015, 2018). Pico-phytoplankton could account up to 80% of the autotrophic biomass in oligotrophic areas (Rodríguez et al., 2005). Smallsized diatoms play an important role in carbon cycle and spring blooms (Leblanc et al., 2018). Additional to autotrophic phytoplankton, small-sized heterotrophic eukaryotic organisms form an important circle of the microbial food web. For instance, heterotrophic nanoflagellates are major grazers of bacteria and other picoplankton, and are ingested by larger protists and metazoans (Lim et al., 1999). It is therefore essential to investigate a wide size range of eukaryotic plankton to characterize the comprehensive community structure to assess the environmental state and ecosystem health.

The identification of plankton with light microscopy has some limits, especially for small organisms or those with less distinguishable morphologies (Nübel et al., 1999; Aguilera et al., 2006). It also needs specific staining assay to reveal the endoparasites in the hosts via microscopic examination (Chen et al., 2019). Furthermore, if the detection method can only cover a small part of the community, it will increase the difficulty of comparing the sample community structure (Fuhrman, 2009). More sensitive methods, including the molecular approach and pigment analysis, have been used to detect the small-sized plankton and proved effective in revealing a high biodiversity of small size eukaryotic plankton communities (Guillou et al., 2008; Cheung et al., 2010; Wohlrab et al., 2018). The biological meaningful operational taxonomic units (OTUs) clustered from high-throughput sequencing results are widely used in plankton community studies (Cheung et al., 2010; de Vargas et al., 2015; Liu et al., 2015; Piwosz et al., 2018; Wohlrab et al., 2018). High-throughput sequencing could provide millions of amplicons and detect rare species in plankton community (de Vargas et al., 2015; Liu et al., 2015).

The Beibu Gulf is a semi-enclosed gulf located in the northwest of the South China Sea (SCS) with an average depth of about 42 m. It is usually an important fishing ground (Huang et al., 2020). Bacillariophyta, Haptophytes, Prochlorococcus, Synechococcus, and Prasinophytes were the most dominant phytoplankton groups in the SCS. The coastal pattern of phytoplankton community in the north of SCS differed from that of the offshore, dominated by diatoms, and accompanied by dinoflagellates, Prasinophytes, and Cryptophytes (Xiao et al., 2018). Although the phytoplankton communities in the open sea and coast of SCS have been investigated intensively (Wu et al., 2014, 2015, 2017; Xiao et al., 2018), the information of phytoplankton community in the Beibu Gulf is limited. The diversity and spatio-temporal variation of phytoplankton with cell size > 20 μm have been well studied via light microscopy in the west coast of Hainan Island (Wang et al., 2012; Xu et al., 2015), but nano- and pico-plankton are hard to be recognized with light microscopy, which may result in underestimation of them in this area. On the other hand, chlorophyll-a (chl-a) concentration was higher in pico- and nano-fractions than in > 20 μm fraction in summer, and pico-phytoplankton, dominated by Prochlorococcus and eukaryotic picophytoplankton, contributed about half of the primary production (Cai et al., 2002). To fully understand the plankton community composition and spatio-temporal variations in Beibu Gulf, investigations using more sensitive methods are demanded for. The environmental factors in the Beibu Gulf, such as nutrient conditions, were sensitive to tropical monsoon climate and human activities (Yuan et al., 2019). In recent years, human activities resulted in eutrophication in the north coastal area of the Beibu Gulf (Kaiser et al., 2013; Qiao et al., 2014). The concentration and ratio of nutrients could affect the biomass or composition of phytoplankton (Hecky and Kilham, 1988). Previous studies reveal that nutrients could influence the spatio-temporal variation of plankton in the Beibu Gulf (Wang et al., 2012, 2014; Xu et al., 2015).

Changjiang Nuclear Power Plant (CNPP) was constructed in the coast of the Beibu Gulf based on the northwest of Hainan Island in 2010. The thermal discharge from CNPP could warm up the local water, increase chlorine concentration, decrease dissolved oxygen, and thus affect the ecosystems (Chen et al., 2010; Prince Prakash Jebakumar et al., 2018). Thermal discharge could also change the biomass, composition or biodiversity of phytoplankton community (Poornima et al., 2005; Li et al., 2011; Lin et al., 2018) and affect the primary productivity (Poornima et al., 2006). Nuclear power plant was considered responsible for the change of dominant species of phytoplankton from diatom to dinoflagellate in the Daya Bay (Li et al., 2011). In addition, thermal discharge could reduce the biomass, increase mortality rate, and change the community structure of copepods (Sanders et al., 1981; Jiang et al., 2009; Lee et al., 2018). Changes in phytoplankton and zooplankton could also result in changes of higher trophic organisms, and then influence the fishery production (Chew et al., 2015). Therefore, to assess the ecological consequences, it is necessary to have a better understanding of the impact of the CNPP operation on the local and nearby plankton communities.

In this study, the community structure of eukaryotic communities near CNPP was investigated via highthroughput sequencing. The community structure of different sizes was examined to better understand the small-sized fractions. The aim of this study is to understand: 1) the size-dependent composition of eukaryotic plankton communities; 2) the relationship between environmental factors and the spatiotemporal variations of eukaryotic plankton communities; and 3) the impact of thermal discharge of CNPP on the community structure of eukaryotic plankton.

2 MATERIAL AND METHOD 2.1 Study site and sample collection

CNPP is located in the southeast of Beibu Gulf and the west coast of Hainan Island, and started operation in 2014. The total length of the drainage project from the outlet to the sea is 4.5 km. The pipes are in the sea area and located on the west side of the water intake (Fig. 1). The outfall is at the depth of about 17 m, and about 0.5 km away from the water intake. The sampling area was about 6 km away from the estuaries of Zhubi River, Shanji River, and Nanluo River (Fig. 1), and represented a tropical coastal water body with seasonal temperature variation. The sampling sites covered a 9-km section perpendicular to the coastline started from the intake (CJ1) of CNPP, and the sampling sites near the outfall (CJ6, CJ2, regarded as outfall in this study), where the temperature was supposed to be influenced by thermal discharge based on the design. The maximum depth of sampling sites was about 20 m.

The study was conducted from February 2018 to January 2019 (Supplementary Table S1). Chl-a concentration, water temperature (WT), dissolved oxygen (DO), pH, and salinity were measured using an YSI EXO Multi-Parameter Water Quality Sonde (YSI, USA). Water samples from surface (< 1-m depth) and bottom (about 2 m above the bottom) were collected every two months. Two-liter water samples were pre-filtered through 200-μm sieve, filtered through 20 μm, 5 μm, and 0.45 μm (PC Membrane, 47 mm, Millipore, Ireland) membranes, and then 20– 200-μm, 5–20-μm, and 0.45–5-μm fractions of plankton were collected, respectively. The plankton assemblages collected through a 20-μm nylon net were also filtered through a 0.45-μm membrane. In November 2018 and January 2019, the plankton assemblages were also filtered directly through 0.22- μm membranes to calculate the proportion of each size-fraction in plankton communities. The membranes were submerged in DNA lysis buffer (10-mmol/L Tris pH 8.0; 100-mmol/L EDTA pH 8.0; 0.5% SDS) and stored in -20 ℃ until DNA extraction. 100-mL filtrate of 0.45 μm were stored at -20 ℃ for the measurements of nutrient concentration. The examination of nutrients was performed using a QuAAtro (SEAL, Germany): ammonium (NH_{+ 4}-N), QuAAtro method number Q-069-05 Rev. 3; nitrate (NO₃^--N), QuAAtro method number Q-035-04 Rev.4; nitrate (NO₂^--N), QuAAtro method number Q-070-05 Rev.4; silicate (SiO₄^4--Si), QuAAtro method number Q-066-05 Rev.2; phosphorus (PO₄^3--P), QuAAtro method number Q-064-05 Rev.4. Total dissolved inorganic nitrogen (DIN-N), was calculated as the sum of ammonium, nitrate, and nitrite.

2.2 DNA extraction and Illumina sequencing

Genomic DNA was extracted according to the protocol of Yuan et al. (2015) and performed as follows. The membranes were cut into pieces with axenic scissors, then incubated in 600-μL DNA lysis buffer (previously added during the sampling) and 10-μL proteinase K (TaKaRa, Code No. 9034, Japan) for 48 h at 55 ℃. The suspension was transferred into a new 2-mL tubes. 82.5 μL of 10% w/v CTAB (Cetyl/ Hexadecyl Trimethyl Ammonium Bromide, in 0.7-mol/L NaCl solution) and 82.5 μL of 5-mol/L NaCl solution were pipetted into each sample, and incubated at 55 ℃ for 10 min. 600 μL of chloroform was added and vortexed. The samples were centrifuged at 11 000×g for 10 min, then the upper layer was carefully transferred into a new tube. The following purification procedure was performed with DNA Clean & Concentration Kit (ZYMO Research ZRC000570, USA) following the protocol.

Pair-end sequencing of the 18S rRNA gene was performed via Illumina Hiseq 2500 using primers (forward 528F, 5'-GCGGTAATTCCAGCTCCAA-3'; reverse 706R, 5'-AATCCRAGAATTTCACCTCCAA- 3') targeting the V4 regions (Cheung et al., 2010; Tang et al., 2019). Reads of V4 regions of 18S rRNA gene with low quality (containing undetermined nucleotides (N), having scores less than 20, or less than 100 bp) were removed via Trimmomatic (V0.33), then pairend reads were concatenated into Unique Tags using FLASH (V1.2.11). The tags less than 200 bp were removed. The assignment of OTUs was set at a 97% cutoff. The most abundant sequence was selected as the representative sequence for each OTU. Chimeras and singleton OTUs were removed via USEARCH (V10). Representative sequences were aligned with the PR² database (http://github.com/vaulot/pr2_database) to obtain a taxonomic framework. After removing the unclassified sequences (9.4% of total tags), Illumina sequencing yielded a total of 19 216 111 tags. An average of 80 067 sequences were obtained for each sample, ranging from 30 711 to 216 443. OTUs of each sample were randomly re-sampled to an equal sample size for further analysis. The taxonomy of abundant and interesting OTU was tested by blast (Supplementary Table S2).

2.3 Statistical analysis

Non-multidimensional scaling (NMDS), analysis of similarities (ANOSIM), distance-based permutation multivariate analysis of variance (PERMANOVA) and Mantel's test were performed using the Vegan library program in R (http://www.r-project.org/). Kruskal-Wallis test was performed via coin package in R. Spearman's rank correlation coefficients were performed via ltm package in R. Venn program were drawn with Venny 2.1 at https://bioinfogp.cnb.csic.es/tools/venny/index.html.

3 RESULT 3.1 Spatio-temporal variation of environmental parameters

Significant seasonal variations of water temperature, salinity, DO, and pH were observed during the survey (Supplementary Fig.S1) (Kruskal-Wallis, P < 0.01). Temperature, DO, and salinity also showed marked differences between surface and bottom water. Water temperature ranged 20.7–29.2 ℃ (Fig. 2b), with the highest value observed in June. Salinity ranged from 32 to 34, and was higher in surface in June, but higher in bottom in August (Fig. 2a). The dissolved oxygen was above 5.8 mg/L during the sampling period and highest in January (Fig. 2c).

The concentration of PO₄^3--P and DIN-N varied seasonally (Fig. 3). PO₄^3--P concentrations were less than 0.015 mg/L in 84% samples, and below 0.030 mg/L in the rest samples. DIN-N of 97% of the samples were below 0.30 mg/L. The concentration of PO₄^3--P was higher in April and January, and lower from June to November. DIN-N and DIN-N꞉PO₄^3--P were lowest in April (Fig. 3 & Supplementary Fig.S2). The compositions of DIN-N also displayed seasonal variations. NO₃^--N was the main form of DIN-N in February and April, when the concentrations of NH₄⁺-N was lower at that time (Fig. 3d–e). Both ratios of Si꞉N and Si꞉P were lowest in April (Supplementary Fig.S2).

3.2 The diversity of eukaryotic plankton

The overlap of individual OTUs among different size fractions in the eukaryotic planktonic samples were analyzed and a substantial number of unique OTUs were revealed in nano- and pico-nanoeukaryotic plankton communities (Fig. 4d). Among the total 4 031 OTUs in four different fractions, 41.1% could be found in all fractions, 5.2% were unique in 0.45–5-μm fraction, 4.5% were unique in 5–20-μm fraction, and 11.0% were only found in 0.45–5-μm and 5–20-μm fractions, respectively. Overall, 20.7% OTUs were unique in small-size fractions (0.45–5 μm and 5–20 μm), and 7.3% OTUs were unique in largesize fractions (20–200 μm and 20 μm net collected).

The Shannon index and OTU number also revealed a higher diversity of pico- and nano-eukaryotic plankton communities in the west side of Hainan Island. The Shannon index at OTU level of 0.45–5- μm and 5–20-μm plankton communities were higher than that in 20–200-μm and 20-μm net-collected communities (Fig. 4a). Moreover, there were more OTUs in 0.45–5-μm and 5–20-μm fractions (Fig. 4b). On the other hand, the beta-diversity among eukaryotic plankton communities were higher in 20– 200-μm and 20-μm net-collected communities than in 0.45–5-μm and 5–20-μm fractions (Fig. 4c). The results of unique OTUs and alpha-diversity revealed that the small-size eukaryotic plankton composition largely differed from that of the large-size eukaryotic plankton at OTU level.

3.3 Eukaryotic plankton communities in different sizes

NMDS and ANOSIM analysis based on OTUs indicated that the community structure of different size was significantly different (Supplementary Fig. S3; Table 1, P < 0.05). The composition of eukaryotic communities is shown in Fig. 5a. Net-collected and 20–200-μm communities were dominated by Metazoa, occupying 75% and 51%, respectively. Meanwhile, Bacillariophyta (32%) was the most abundant class in 5–20-μm fractions. Eukaryotic communities of 0.45–5 μm were predominated by Chlorophyta (17%), Metazoa (17%), and Syndiniales (17%). Dinophyceae was most abundant in 5–20-μm fractions, accounting for about 14%. Syndiniales, Ellobiophyceae, Cryptophyta, and Chrysophyceae were most abundant in 0.45–5-μm fraction, occupying 17%, 3.3%, 5.3%, and 3.2%, respectively.

At genus level (Fig. 5b–e), Micromonas (Chlorophyta), Ostrecoccus, Mantoniella, and Ellobiopsis (Dinoflagellates) were most predominant in 0.45–5-μm plankton communities. Paraphysomonas (Chrysophyceae) and Minutocellus (Bacillariophyta) were more predominant in 0.45–5 and 5–20-μm fractions than in 20–200-μm fraction. Tripos (Dinophyceae), Actinocyclus (Bacillariophyta), Ditylum (Bacillariophyta) and Hemiaulus (Bacillariophyta) were most dominant in 20–200-μm fraction. Protoperidinium (Dinophyceae) and Coscinodiscus (Bacillariophyta) were more predominant in 20–200-μm and 20-μm net collected fractions.

The composition of eukaryotic plankton communities above 0.22 μm revealed a remarkable proportion of pico- and nano- organisms. Syndiniales, Micromonas, and Minutocellus, which was more predominant in 0.45–5- and 5–20-μm fractions, occupying over 15%, 7.5%, and 1.6% in eukaryotic plankton communities in November 2018 and January 2019, respectively.

3.4 Spatio-temporal dynamic of eukaryotic plankton communities

ANOSIM and NMDS results revealed the significant difference between plankton communities in surface and bottom waters. On the other hand, plankton communities were not significant different among sampling sites (Supplementary Fig.S3b, Table 1). Autotrophic Chlorophyta and Bacillariophyta were more predominant in surface water, and heterotrophic Syndiniales, Ellobiophyceae were more predominant in bottom water (Fig. 6).

ANOSIM and NMDS results also showed a seasonal variation of eukaryotic plankton community in each size-dependent fraction (Supplementary Fig. S3c–S3f, Table 1). The seasonal dynamic of taxonomic groups is shown in Fig. 5. Chlorophyta was more predominant in June and August. Specifically, Micromonas, the dominant genus of Chlorophyta, was more predominant in June and could reach 37% in the 0.45–5-μm plankton communities. Bacillariophyta was more predominant in 20–200-μm fractions in January and February, and more predominant in 5–20-μm fractions in February. Chrysophyceae was more predominant in February, even accounting for 67.2% in the 0.45–5-μm plankton communities in surface water of station CJ4 (Supplementary Fig.S4). Paraphysomonas, the most dominant genus of Chrysophyceae, were more predominant in February, and occupied about 11.6% in the 0.45–5-μm plankton communities. Syndiniales (dominated by Group Ⅰ, Ⅱ, and Ⅲ) occupied over 7.8% of 0.45–5-μm eukaryotic plankton during the study period, and could occupy 46% in the bottom water of station CJ5 in April (Supplementary Fig.S4). The proportion of different groups of Syndiniales also displayed a seasonal variation. Clade 5 of Group Ⅰ was more predominant in January, whereas clade 4 of Group Ⅰ was most predominant in June (Fig. 5e).

3.5 The correlation between eukaryotic communities and environmental factors

The spatio-temporal dynamics of eukaryotic plankton communities were correlated with temperature, salinity, and nutrient conditions, including SiO₄^4--Si, PO₄^3--P, and DIN-N (Fig. 7, PERMANOVA, P < 0.05). Plankton communities were more correlated with PO₄^3--P (P < 0.01) than DIN-N (P < 0.05 in surface samples) and NO₃^--N (P < 0.05 or not significant in 20– 200 μm and 5–20-μm fractions) Spearman's rank correlation tests showed the effect of environmental factors upon particular organisms (Fig. 8). The relative abundance of Bacillariophyta was negatively correlated with water temperature, and positively correlated with the concentration of SiO₄^4--Si and PO₄^3--P. The correlation of eukaryotic organisms and nutrient concentrations were sizedependent. The proportion of Bacillariophyta in netcollected communities was positively correlated with DIN-N and DIN-N꞉PO₄^3--P, but negatively correlated with DIN-N꞉PO₄^3--P in 20–200- and 0.45–5-μm fractions. The relative abundance of Chlorophyta and Dinophyceae did not correlate with the concentration of DIN-N in 20–200-, 5–20-, and 0.45–5-μm fractions, but correlated with DIN-N in net-collected fractions.

At genus level, although Minutocellus, Skeletonema, Thalassiosira, and Actinocyclus were negatively correlated with water temperature, Leptocylindrus and Palmerina were positively correlated with water temperature. The genera of Chlorophyta showed different correlation with water temperature. Micromonas and Mantoniella were positively correlated with water temperature, while Ostreococcus was negatively correlated with water temperature (Fig. 8). In addition, Ostreococcus was negatively correlated with salinity but Micromonas was not.

3.6 Interactions among organisms

The potential interactions among organisms were also detected. The two typical parasites Syndiniales and Ellobiophyceae, played an important role in the dynamic of plankton communities. At community level, Mantel's test showed a significant correlation between these parasites and other eukaryotic plankton during the spatio-temporal variation (Supplementary Fig.S5). At genus level, some particular parasites displayed co-occurrence with their hosts. For example, Syndinium was significantly correlated with Maxillopod (over 98% of which was occupied by copepods in this study, Supplementary Table S3).

3.7 The effect of thermal discharge on the eukaryotic plankton

The difference of temperature among sampling sites was less than 1 ℃, and was smaller in bottom water than in surface water (Table 2). The DO was above 5.8 mg/L during the sample period (Fig. 2c). The differences among bottom of each site were less than 0.2 mg/L (Table 2). Two-way ANOVA revealed that the OTU number, Shannon index at OTU level, and Dinophyceae꞉Bacillaroiphyta changed among study periods, but thermal discharge had no significant effect (Table 3). Two-way PERMANOVA based on Bray-Curtis distance also showed that thermal discharge had no significant effect on community structure, but there was a significant interaction between the effects of sampling time and thermal discharge (Table 3).

4 DISCUSSION 4.1 Biodiversity of different-sized fractions

A higher alpha-diversity in small-size fractions (0.45–5 μm, and 5–20 μm) and a higher beta-diversity in large-size fractions (20–200 μm, and > 20 μm) were revealed in this study. There are two possibilities for higher beta-diversity in large-size fractions. First, at equal sequencing depth, the multiple copies of multicellular organisms might contribute to a higher beta-diversity. Secondly, body size can affect the spatial distribution of passive dispersers, such as phytoplankton, rotifer, and mollusc, whose dispersal ability is positively correlated with body size (Shurin et al., 2009; de Bie et al., 2012).

Filtration might introduce some fragments or cell breakage of large-scale organisms, which could be detected by molecular methods in small-sized fractions (Massana et al., 2015). This bias occurred in this investigation too because a certain amount of metazoans was detected in 0.45–5-μm fraction. Although the higher alpha-diversity of small-sized fractions might be the result of filtration bias, the large amount of unique OTUs in small-sized fractions still showed that the small-sized eukaryotic plankton composition largely differed from those of the largesized. This result supported the idea that piconanoplankton had high alpha-diversity in global sunlit oceans (de Vargas et al., 2015), and indicated that small-size plankton could make a substantial contribution to community diversity in this area.

4.2 Eukaryotic plankton communities of different sizes

Diatoms, together with dinoflagellates, Prasinophytes and Cryptophytes have been reported to be dominant in the coastal waters of SCS (Xiao et al., 2018). Small size diatoms are overlooked but may play an important role in carbon cycle and blooms (Leblanc et al., 2018). Although previous studies have reported the diversity and spatiotemporal variation of diatoms in the Beibu Gulf (Wang et al., 2012; Xu et al., 2015), these results were based on light microscopy and might miss the small size diatoms. The investigation based on Illumina sequencing revealed a remarkable proportion of Minutocellus in 0.45–5- and 5–20-μm fractions, and it was negatively correlated with water temperature and salinity. This result was consistent with the previous research, that pico- and nanoeukaryotic phytoplankton in Beibu Gulf were influenced by nutrients and negatively correlated with salinity (Zhao et al., 2019). The findings of considerable amounts of Minutocellus suggested that the diversity of Bacillariophyta in the Beibu Gulf was underestimated.

Prasinophytes are widely distributed in the global oceans. Previous studies have also shown that there are a large number of Prasinophytes in the coastal area of the SCS (Wu et al., 2014, 2017), but the detailed information of Prasinophytes in the Beibu Gulf is missing. This study firstly reports that Micromonas has a relatively high level among the plankton community in the coast of Beibu Gulf. Although Micromonas was the most dominant Chlorophyta in this area, Ostreococcus was also detected and its vertical distribution and correlations with environmental factors differ from Micromonas. For instance, Ostreococcus is negatively correlated with water temperature and salinity, and more abundant in surface water. Previous research has pointed out that both Micromonas and Ostreococcus have high abundance in the estuary and slope region, but Ostreococcus might be more adapted to estuary ecosystem than Micromonas (Wu et al., 2014). That might help explain why Ostreococcus was negatively correlated with salinity in the study area.

Pigment-based results indicated that pico-and nano-phytoplankton accounted for 50% of total chl a in the Beibu Gulf, and pico-phytoplankton contribute 50% of primary productivity (Cai et al., 2002). Eukaryotic microalgae are important parts of picophytoplankton in the Beibu Gulf (Zhao et al., 2019). The remarkable proportion of Micromonas and Ostreococcus detected in this study helps to explain the contribution of pico-phytoplankton to the primary productivity in this area. Our results also show that Micromonas was positively correlated with temperature, which is consistent with previous observation in the north of SCS (Wu et al., 2014).

Additional to the autotrophic phytoplankton, heterotrophic pico- and nano-eukaryotic plankton have also been detected in the study area. Syndiniales are host-specific parasitic dinoflagellates, and they could infect zooplankton, phytoplankton, fish, and marine crustaceans (Shields, 1994; Guillou et al., 2008; Skovgaard et al., 2009; Coats et al., 2012; Chen et al., 2019; Huang et al., 2019). Ellobiopsidae is also a parasite of copepod (Artüz, 2016). Syndiniales and Ellobiophyceae have widespread distribution in the ocean (Shields, 1994; Guillou et al., 2008), and Amoebophrya-like organisms (Syndiniales) have been found across the coastal waters of China (Chen et al., 2019), but the detailed information about these parasites in the Beibu Gulf is still limited.

In this study, we found a substantial proportion of Syndiniales and Ellobiophyceae in the eukaryotic plankton communities, and the variation of Syndiniales and Ellobiophyceae composition was correlated with the variation of other eukaryotic plankton proportion. The dynamic of parasites is correlated with their hosts. For instance, Syndinium, which could parasitize copepods and lead to the mortality of copepods (Shields, 1994), was significantly correlated with their hosts (Maxillopod) in this study. In addition to parasitic in the hosts, Syndinium and Ellobiopsis also released small-sized spores (Konovalova, 2008), which might explain the fact that Syndinium was correlated with Maxillopod in both net-collected and 0.45–5-μm fractions. Previous studies reported that parasitic dinoflagellates might affect the dynamic of their hosts, including algae and fish (Coats and Bockstahler, 1994; Huang et al., 2019). The co-occurrence of parasites and their hosts indicated that Syndiniales and Ellobiophyceae might play an important role in the spatial and temporal dynamic of eukaryotic plankton in the study area.

Heterotrophic flagellate Paraphysomonas, which typically range from 4 to 8 μm in diameter (Lim et al., 1999), is cosmopolitan in marine environments (Scoble and Cavalier-Smith, 2014). It is omnivorous and can feed on bacteria and algae (Lim et al., 1999), Paraphysomonas is an important component of the autotrophic picoplankton-metazooplankton food chain (Bec et al., 2006). The high proportion of Paraphysomona suggests that Paraphysomonas should not be ignored in the food web and energy flow of ecosystem in this area. In addition, the 0.45– 5-μm eukaryotic community was predominated by Paraphysomonas (heterotrophic) in February, and by Micromonas (autotrophic) in June. In April, November, and January, parasitic dinoflagellates were the most predominant organisms. The variation of predominant organisms displaying different trophic types indicates that the structure of microbial food web and energy flow might display seasonal variations.

4.3 Environmental factors and the spatio-temporal variation of eukaryotic plankton communities

Our results show that the dynamic of plankton communities is correlated with silicate, phosphorus, and nitrogen concentration, and phosphorus played a more important role than DIN-N did. A ratio of Si and N below 1꞉1 limited the growth of diatoms (Turner et al., 1998). In this study, the Si꞉N atomic ratio was below 1꞉1 in April, at three sampling sites in June. These might explain the positive correlation between silicate concentration and the relative abundance of Bacillariophyta. In addition, nitrogen and phosphorus were absorbed by phytoplankton at a ratio of about 16꞉1 (the Redfield ratio) (Redfield, 1958; Zhu et al., 2013). In this study, N꞉P atomic ratio was above 16:1 except for in April. The high N꞉P ratios have also been found in the west coast of Hainan Island by Xu et al. (2015). This might explain that plankton communities were more correlated with PO₄^3--P than with DIN-N, and the relative abundances of Chlorophyta and Bacillariophyta in the eukaryotic communities were positively correlated with phosphorus but DIN-N in < 200-μm fractions, in which Chlorophyta and Bacillariophyta were predominant. Additional to phytoplankton, zooplankton in the northern Beibu Gulf has also been controlled by phosphorus and nitrogen (Wang et al., 2014), which emphasized the influence of phosphorus and nitrogen upon the annual spatio-temporal variation of eukaryotic plankton communities in the study area.

The correlation of eukaryotic organisms and nutrient concentrations were size-dependent. Bacillariophyta is positively correlated with DIN-N in net-collected fraction, which might be explained by the fact that different genera displayed different patterns with environmental factors. For example, Coscinodiscus, which was most predominated in netcollected and 20–200-μm fractions, was positively correlated with DIN-N concentration. In addition, eukaryotic phytoplankton is also influenced by other environmental factors including temperature and salinity. For example, Leptocylindrus(Bacillariophyta) is negatively correlated with silicate, and positively correlated with temperature and salinity. These indicate that the influence of environmental factors upon plankton communities was complex and further studies are needed.

4.4 The effect of thermal discharge on the eukaryotic plankton

The thermal discharge of nuclear power plants could cause serious environmental consequences, including the rise of water temperature and reduced dissolved oxygen (Chen et al., 2010; Prince Prakash Jebakumar et al., 2018). In the present study, temperature rising can hardly be detected near the outfall, suggesting that the effect of thermal discharge to the water temperature was slight, which might be due to the high flow rate (45–123 cm/s) near CNPP (Wang et al., 2012).

For the plankton communities, it was reported that thermal discharge could change the composition or biodiversity (Jiang et al., 2009; Li et al., 2011; Lin et al., 2018). However, in this study, except for the synergistic effects of thermal discharge and seasonal variation, thermal discharge showed little effect on both alpha-diversity (Shannon index and OTU numbers) and community composition of eukaryotic plankton. In addition, an investigation of phytoplankton in a subtropical bay revealed that temperature rising caused by the thermal discharge could induce the shift of dominant species from diatom to dinoflagellates (Li et al., 2011), but the increase of Dinophyceae꞉Bacillariophyta did not appear in this study. The little influence of thermal discharge to plankton community might be a result of the slight increase of water temperature (< 1 ℃) compared with the remarkable temperature rising observed in some other power plants (Supplementary Table S4). The other possibility of the slight effect of thermal discharge might be that CNPP is located near the open sea. Some power plants located inside bays (Jiang et al., 2009; Li et al., 2011; Lin et al., 2018) had negative effect on plankton, but a study in an open sandy tropical coast also found the distribution and abundance of phytoplankton in the coastal water was not affected (Poornima et al., 2005).

5 CONCLUSION

A high diversity of pico- and nano-plankton community in the coast of west Hainan Island was revealed, and the presence of autotrophic phytoplankton, including Micromonas and Minutocellus, and heterotrophic plankton, including Paraphysomonas, Syndiniales, and Ellobiopsis was remarkable. These results indicate that small-sized eukaryotic plankton should not be overlooked as they were an important component of the ecosystem. The community structure displayed significant seasonal variations similar among sampling sites. The spatialtemporal dynamic of plankton was correlated more significantly with phosphorus than with DIN-N in this area. Our results also reveal that the impact of thermal discharge on eukaryotic plankton communities in this open sea area was minor. In the future, the complex environmental factors and interactions among organisms that drive the dynamic changes of eukaryotic plankton communities will be further investigated.

6 DATA AVAILABILITY STATEMENT

The datasets generated and analyzed during the current study are available from the corresponding author on reasonable request.