Six cruises were carried out from October 2018 to March 2019 in Qinzhou Bay (21.55–21.76°N, 108.53–108.72°E). Fifteen stations were set up for sampling during each cruise ( Figure 1 ). Water samples were collected from a depth of approximately 0.5 m at each station using 2.5 L Niskin bottles for the determination of the nutrient concentrations, phytoplankton biomass and species composition.

A multiparameter water quality analyzer kit (JFE AAQ171, Japan) was used to measure water temperature, pH, salinity, and dissolved oxygen (DO). The pH and DO sensors was calibrated before each cruise according to the calibration procedures provided by the manufacturer. Two type of standard solution are used for calibration on pH sensor, including standard 7 and standard 9. When the DO sensor is to be celibrated, fill a dedicated container with half full of water and prepare a 100% saturated water by air pump, then submerge the instrument in it. The steps specified for the calibration were carried out in accordance with the instrument communication software user’s guide.

Water samples were filtered through precombusted Whatman GF/F filters (450°C for 4–5 h) and stored at –20°C. The samples were thawed before they were analyzed, and the concentrations of dissolved inorganic nitrogen (DIN, including NO 3 − , NO 2 − , and NH 4 + ), dissolved inorganic phosphate (DIP), and dissolved inorganic silicate (DSi) were analyzed using a flow injection analyzer (Skalar San++ autoanalyzer, Holland).

One-liter seawater samples were mixed with 2% acid Lugol’s solution and concentrated gradually to 20 mL prior to quantification. A 100 μL concentrated sample was transferred to a 0.1 mL counting chamber. Standard methods (inverted microscopy (Olympus BX51) and Utermöhl were used for phytoplankton species identification and cell enumeration.

500 mL of water sample was poured into a white disk in batches, and lit up by a lamp for counting the density of P. globosa colonies. Colony size was roughly recored by a dividing rule. After the sample concentrating steps, the colonies of P. globosa were broke into pieces when they were absorbed back and forth for several times by a pipette. The cell density of P. globosa was calculated as the sum of solitary cells in the water and cells in the colony debris.

The dominant species of phytoplankton were determined according to the dominance index (Y), which was calculated by the following formula:

where N_i is the cell number of the specific species i, f_i is the frequency of species i appearing in the survey area, and N is the total cell number of phytoplankton at all stations. Species were defined as dominant species when their Y≥0.02.

In addition to the dominance index, the Shannon-Wiener (Shannon, 1948), Margalef (Margalef, 1968), Simpson (Simpson, 1949) and Pielou’s evenness (Pielou, 1966) diversity indices of phytoplankton species were calculated for the study period.

Seawater samples (1000 mL) were filtered onto Whatman GF/F filters and stored at –20°C. The Chl a concentration was measured by the fluorometric method according to Yentsch and Menzel (1963). Ten milliliters of 90% acetone was added to the tubes with Chl a filters; the tubes were stored for 24 h in the dark at 4°C and then centrifuged at 4000 rpm for 10 min. Fluorescence before and after acidification with 3.7% HCl was recorded using a fluorospectrophotometer (Perkin Elmer LS55, Ex = 470 nm, Em = 680 nm).

One-way ANOVA was performed using SPSS 19 software to compare the difference between sets of parameters among cruises using Tukey’s studentized range test (t test), where a p value of 0.05 was used to determine the significance of the results.

Biological and environmental (BioEnv) analysis was performed using R Studio software to quantify the relative contributions of the abiotic variables to the observed patterns of phytoplankton for each cruise, and then the best models with the highest correlations were obtained.

Principal component analysis (PCA) was used to explore the influences of different sampling stations on environmental variables. Redundancy analysis (RDA) was used to evaluate the influence of abiotic variables on the abundances of the dominant species. Only those species with the 10 highest dominance index values (Y) were selected for the ordination analysis. All canonical axes were applied to evaluate the significant variables under analysis with Monte Carlo permutation tests (499 permutations). The above mentioned analyses were performed using Canoco 5 software.

The water temperature decreased from October to January and then increased from January to March ( Figure 2A ; Table 1 ). There were significant differences among the six months (one-way ANOVA, p< 0.01). The area close to the coastline had a lower temperature than the far shore area ( Figure S1 ), which was mainly dependent on latitude. The salinity varied from 15.7 to 31.8 and increased gradually from the inner bay area to the outer area ( Figure 2B ; Figure S2 ). The salinity in October was significantly lower than that in any other month except January (one-way ANOVA, p< 0.05, Table 1 ). The salinity did not show significant variation among the five cruises from November 2018 to March 2019 (one-way ANOVA, p > 0.05). The average pH remained relatively stable throughout the study, varying from 7.97 to 8.14, and the highest average pH was recorded in January ( Table 1 ; Figure 2C ). According to the correlation analysis, pH was significantly positively correlated with temperature and salinity (p< 0.01), while it was negatively correlated with DIN, DIP and DSi (p< 0.01). DO was highest in January and then December, and the DO in these months significantly higher than that in the other four months (one-way ANOVA, p< 0.05, Table 1 ; Figure 2D ). The DO was lowest in October when the highest temperature of the study period occurred. In fact, DO was mainly dependent on temperature, as an opposite trend was found between the two parameters ( Table 1 ; Figure 2D ).

Except for January 2019, the spatial distributions of DIN, DIP, and DSi were consistent during all the cruises. In contrast to salinity, DIN, DIP and DSi decreased gradually from the inner bay to the outer region ( Figures 3 – 5 ). The concentrations of all the inorganic nutrients were lowest in January 2019 ( Table 1 ). The average DIN was the highest in March 2019, significantly higher than that in any other month except October 2018 (p< 0.01), and the highest DIN (49.4 μmol/L) was recorded in the surface water of station 1. DIP showed no significant differences among October, November and December 2018 (t test, p > 0.05), but it decreased sharply in January 2019, i.e., below 0.1 μmol/L at all survey stations ( Figure 4D ), increased gradually from February to March, and reached its highest value in March with an average of 0.57 μmol/L. The DSi in January was significantly lower than that in November, February and March (one-way ANOVA, p< 0.01), while there were no significant differences between January and October or between January and December (one-way ANOVA, p > 0.05).

During the study periods, the DIN/DIP values were higher than 22 in all months except December 2018, in which the average DIN/DIP was 15.4 ± 5.2 ( Table 1 ). The average DSi/DIP ratio ranged from 21.6 to 273.7 during the six cruises, which was the same as the DIN/DIP. The DSi/DIP ratio also had the lowest average value in December 2018, while the ratio was higher than 22 in the other five months ( Table 1 ).

The PCA results ( Figure 6 ) showed that the inorganic nutrient (DIN, DIP and DSi) vectors were positively correlated with PCA axis 1, while the salinity, Chl a, pH, DIN/DIP and DSi/DIP vectors were negatively correlated with PCA axis 1. Moreover, DO was positively correlated with PCA axis 2, and the temperature vector was negatively correlated with PCA axis 2. Samples from the six cruises were distributed along PCA axis 2 with obvious differences, showing that temperature and DO were the main ecological constraints that affected the samples on the temporal scale. The distribution of samples in the same cruise was scattered, and samples near the inner bay (especially stations 01-04) were located closer to the DIN, DIP and DSi vectors, illustrating that these samples were richer in inorganic nutrients.

The highest average Chl a was observed during the cruise in January, and it was significantly higher than that in any other cruise ( Table 1 ) (one-way ANOVA, p<0.01). Distribution tendencies of Chl a were different between cruises ( Figure S3 ). A significant positive correlation was observed between Chl a and the abundance of phytoplankton (p< 0.01).

A total of 193 species (including unidentified species and varieties) of phytoplankton were found during the six cruises, including 141 diatoms and 43 dinoflagellates. The total species number for each cruise varied from 97 to 126, and the highest and the lowest numbers of species were found in December and January, respectively. Over 75% of the taxa were diatoms in the first five months, and 66.7% were diatoms in March. The average total abundance of phytoplankton varied from 0.3 ± 0.2 × 10⁵ cells/L to 43.4 ± 14.3 × 10⁵ cells/L, and phytoplankton abundance was significantly higher in January than in the other five cruises (t test, p< 0.01). The total phytoplankton abundance was higher in the eastern part of the study area than in the western part during all cruises except in January ( Figure 7 ), while it was the highest at station 3, which was closest to the nuclear power plant, in January ( Figure 7D ). Diatoms and P. globosa were previously the most abundant algal groups from October to February ( Figure 8 ). The abundance of diatoms varied from 0.2 ± 0.2 × 10⁵ cells/L to 5.4 ± 3.2 × 10⁵ cells/L, accounting for 12.6 ± 6.8% to 90.2 ± 19.7% of the total phytoplankton abundance during the six cruises. The average cell density of P. globosa ranged from 0 to 3.8 ± 1.3 × 10⁶ cells/L and accounted for 0 to 87.0% of the total abundance during the six cruises. The highest cell densities of both diatoms and P. globosa appeared in January 2019, and they were significantly higher than those in other cruises (t test, p< 0.01).

The succession of the phytoplankton community was very obvious throughout the study period ( Figure 9 ), and there were also monthly variations in temporal and spatial distributions of P. globosa (Figure 10). During the first cruise in October 2018, six dominant species were found in Qinzhou Bay, including P. globosa and five diatom species. Among these species, P. globosa was found at all stations ( Figure 10A ) and had the highest dominance index (0.319), with an average cell density of 8.2 ± l8.0 × 10⁴ cells/L, accounting for 31.6% of the total phytoplankton abundance ( Figure 8 ). The number of dominant species increased to 11 in November; ten of these species were diatoms, among which 9 species had chain-like structures. Although P. globosa was found at only three stations ( Figure 10B ), its high density made it one of the dominant species. Chaetoceros species and P. globosa continued to be the dominant species in December. P. globosa abundance increased dramatically, and this species became the most abundant species in January 2019, accounting for 87.0% of the total phytoplankton, with an average cell density of 3.8 ± 1.3 × 10⁶ cells/L. Only two dominant species were found during this cruise. P. globosa remained one of the dominant species in February, together with Plagioselmis prolonga (Cryptophyta) and Rhizosolenia alata (Bacillariophyta). Dinoflagellate (Dinophyta) species, including Oxyrrhis marina, Prorocentium micans and Gonyaulax spinifera, were found to be dominant species for the first time during the cruise in March 2019, and dinoflagellates replaced P. globosa and had the second-largest abundance following diatoms.

P. globosa colonies were just found at one or two stations in the cruises excluding January, and had low densities (< 5 colonies/L). While colonies were detected at nearly all the stations in January, and the colony density ranged from 0 to 125 colonies/L, with an average of 43.7 ± 22.2 colonies/L. The percentages of colonial cells to total P. globosa cells varied from 0 to 84.5%, with an average of 31.1 ± 22.2%. Both the highest colony density and the highest percentage of colonial cells were recorded in station 12 ( Figure 11 ). All the visable colonies was smaller than 0.4 cm.

Four phytoplankton community diversity indices were calculated to analyze the relationship between the phytoplankton community and the species density and disturbance degree of a colony to environmental factors. During the study period, the trends of all four diversity indices were consistent: the phytoplankton community structure stability was stable during October to December but decreased rapidly in January ( Figure 12 ). The average value of the Shannon-Wiener index ranged from 0.99 to 4.24 throughout the tudy period. The highest average value was observed in November, and the lowest average value was recorded in January. The Margalef diversity index ranged between 1.95 and 3.29. The lowest and highest values were recorded in January and December, respectively. The values of both the Pielou evenness index and the Simpson index were lower than 1. Higher values indicated higher relative diversities and more evenly distributed species. The trend of variation in the Simpson evenness index was almost the same as that in Pielou’s evenness index. P. globosa abundance increased greatly in January, resulting in low values of the diversity indices except for the Margalef index ( Figure 12 ).

BioEnv analysis based on Spearman’s correlation matrix was performed to determine which variable(s) affected phytoplankton abundance and distribution in this study. The best models were selected according to correlations of the variables to phytoplankton communities. The results revealed that the variables that best explained the phytoplankton spatial pattern were variable among the cruises ( Table 2 ). In October, salinity, DSi and pH were the main factors that statistically correlated the variations in the composition of phytoplankton, with correlations of 0.5224. The best model of environmental variables included just temperature and salinity in November, while temperature, salinity, DIN, DIP were contained in January, with high correlations of 0.6114 and 0.6907.

As the length of first DCA was between 3 and 4 (3.08), the relationship between environmental variables and the main phytoplankton species (i.e., those with one of the 10 highest dominance index values (Y) in this study) throughout the study period was analyzed by RDA. The first two axes explained 26.6% and 5.5% of the cumulative variance in the relationship between species and environmental variables, with eigenvalues of 0.266 and 0.055, respectively.

The ordination graph of RDA shows that the succession of phytoplankton was mainly driven by hydrological factors and inorganic nutrients ( Figure 13 ). Instead of having multiple dominant species each month, there were only two dominant species (Y > 0.02) during the entire study period, P. globosa (Y = 0.382) and C. lorenzianus (0.022). The 10 species with the highest dominance index values from all of the samples met the data selection criteria and, consequently, were summed as the species data matrix. One of the predominant species, P. globosa, associated with G. striata and T. rotula, was characterized by high DIN/DIP. However, the other dominant species, C. lorenzianus, had a preference for high pH, and the vectors of C. debilis, C. costatus and L. danicus were close to that of C. lorenzianus in this graph. The other species, C. curvisetus, C. debilis, C. teres, and R. alata showed a particular preference for high pH and high salinity.

Qinzhou Bay in a broad sense includes the inner bay (Maowei Sea) and the outer bay, where most of the stations in this study are located (stations 3 to 15, Figure 1 ). In addition to latitude playing a key role in temperature, the spatial distributions of temperature and salinity in the study area were affected by river input from the Maowei Sea and seawater intrusion from the offshore area ( Figures S1 , S2 ). As a tropical bay, the temperature was obviously different in different months ( Table 1 ; Figure 2A ); the maximum temperature was recorded in October 2018 (25.6 ± 1.0°C), while the minimum temperature was recorded in January 2019 (14.3 ± 0.5°C). Temporal variations in DO were almost the opposite of those of temperature during the six months ( Table 1 ), so the vectors of DO and temperature were in opposite directions in the PCA results ( Figure 6 ). According to the PCA results, we found that temperature and DO were important for the monthly environmental changes during the study period.

Nutrients in Qinzhou Bay were mainly affected by river input, as the Qinjiang River and Maolingjiang River are linked to the inner bay (Lan et al., 2012; Lan, 2012; Hong et al., 2014). The concentrations of DIN, DIP and DSi in all cruises except that in January decreased from the inner bay to the outer region ( Figures 3 – 5 ) and showed a significant negative correlation with salinity during the investigation. Nutrients carried by river inputs in this area were calculated and accounted for more than 78% of the total nutrient input (Xu et al., 2012). Dissolved inorganic nutrient contents (DIN, DIP and DSi) were relatively low in the study area during this study, possibly because the rivers were in a low-flow season during this period, and the water was mixed with “clean water” from the outer sea. In addition, the outer bay is far from the estuary and has less land source sewage discharge and less mariculture area (Xu et al., 2012; Kang et al., 2020a). The samples from all cruises except those in January were dispersed along axis 1 in the PCA result, which was mainly associated with DIN, DIP, DSi and salinity ( Figure 6 ).

In the last several decades, environmental pressure and eutrophication in Qinzhou Bay have been increasing with the rapid development of industrialization and urbanization. DIN was the nutrient that increased most, while DIP remained relatively stable; eutrophication was found only in the inner bay, and oligotrophication was found in the outer bay (Wei and He, 2008; Xu et al., 2012). DIN/DIP ratios are often used to determine whether an aquatic system is experiencing N or P limitation (Redfield, 1958; Redfield et al., 1963). Phytoplankton are considered to become P limited when DIN/DIP is greater than 20-30 (Healey, 1979). A DIP of 0.1 µmol/L is considered to be the minimum threshold for phytoplankton growth; that is, when DIP is lower than this value, phytoplankton are absolutely P limited (Peery and Eppley, 1981; Dortch and Whitledge, 1992). The DSi/DIP ratio was used in combination with the DIN/DIP ratio to determine if there was a P limitation, and when both ratios were greater than 22, the phytoplankton were thought to be in a P-limited state (Justić et al., 1995). Both the DIN/DIP and DSi/DIP in Qinzhou Bay were higher than 22 in this study during five cruises, excluding the cruise in December ( Table 1 ), and the DIP at all the stations was lower than 0.1 µmol/L in January, indicating P limitation in this area. In fact, Qinzhou Bay has long been considered a P-limited area, as the average N/P and Si/P throughout the year were 116 and 152, respectively (Ma et al., 2013; Lin and Zhong, 2014).

A total of 193 species of phytoplankton were identified during the study, and this number was consistent with results obtained in 2010 (Jiang et al., 2012). However, only 148 species were recorded in 2013-2014 (Liu et al., 2020); this variation might be attributed to the difference in sampling stations and survey times. The dominant species changed gradually between cruises. P. globosa was the dominant species in the first five months, and C. curvisetus and C. lorenzianus were additional dominant species that appeared frequently (from October to December, Figure 9 ). Other studies have shown that Chaetoceros was the main dominant species in Qinzhou Bay (Jiang et al, 2012; Liu, 2017). Chaetoceros species were frequently dominant from October to December, the period during which Phaeocystis blooms began in previous years. Abundant Chaetoceros were found during the prebloom period (Peperzak et al., 2000), and many newly formed colonies were found on the setae of Chaetoceros, suggesting that Chaetoceros plays a key role in the inception of Phaeocystis blooms (Rousseau et al., 1994). It was reported that a certain cell density of C. minimus and C. pseudocurvisetus will significantly increase the number of P. globosa colonies and prolong the existence time of the colonies (Li et al., 2010; Liu et al., 2021). In fact, the coexistence of Phaeocystis colonies with diatoms is a well-known and recurring phenomenon in field observations (Rousseau et al., 1994; Peperzak et al., 1998; Wassmann et al., 2010). Diatom species, especially those with thorn, spine or setae structures, are commonly observed (Liu et al., 2021). These cell structures are regarded as beneficial for improving buoyancy or providing a substrate for P. globosa attachment (Rousseau et al., 1994; Sazhin et al., 2007; Liu et al., 2021). Moreover, the allelopathic effect of diatoms plays a key role in the effects of diatoms on P. globosa colony development. For example, the cell-free filtrates of Thalassiosira weissflogii and Ditylum brightwellii stimulate the formation of more and larger P. globosa colonies (Liu et al., 2021; Wang et al., 2021c). Some allelochemicals can support a reallocation of resources involved in forming and maintaining colonies (Mars Brisbin and Mitarai, 2019), while others could provide necessary nutrients and energy to support P. globosa colony formation (Zhang et al., 2020).

The phytoplankton community diversity indices were quite high from October to December, and all four indices decreased in January, which was when P. globosa had the highest density ( Figure 12 ). Similarly, Liu et al. (2020) found that the phytoplankton diversity indices were high and the community was stable in Qinzhou Bay, but this characteristics changed when P. globosa blooms occurred. Significant changes were also found in both the quantity and the variety of phytoplankton-dominant species starting in January, when there were only two dominant species ( Figure 9 ). Those Chaetoceros species that had a chain-like structure and spines on the cells were no longer dominant species after January. There are three possible reasons for the decline in diversity indices in January. One was that the decline in temperature in January led to a decrease in Chaetoceros abundance, as these species are thought to grow better at high temperatures (Wang et al. 2010a). Second, the DIP was relatively low in January, which might limit the growth of Chaetoceros. Third, the rapid reproduction of P. globosa might result in a reduction in Chaetoceros. The phytoplankton community always tends to be uniform during the formation of blooms, and the diversity decreases (He et al., 2008; Xu et al., 2022). Thus, diversity indices were proposed as a predictor of the occurrence of algal blooms (He et al., 2008). Single-cell species appeared successively; for example, P. prolonga (Cryptophyta) was recorded in February and dinoflagellates (O. marina and P. micans) were recorded in March. Colonies of Phaeocystis commonly attach to Plagioselmis species and heterotrophic dinoflagellates before the appearance of Phaeocystis blooms (Peperzak et al., 1998; Shields and Smith, 2008). Moreover, the abundance of heterotrophic dinoflagellates feeding on cells of Phaeocystis increases gradually during blooms (Grattepanche et al., 2011). Thus, it was speculated that the increase in heterotrophic dinoflagellates after the decline in P. globosa in March might be due to feeding on P. globosa cells throughout the investigation period.

According to the results of the BioEnv analysis, the factors regulating the phytoplankton assemblages during the six cruises differed month to month. Temperature, salinity and pH were important factors from October to December when the phytoplankton community was mainly composed of diatom species ( Table 2 ). More variables (temperature, DIP, salinity, DIN) were found to be important by the BioEnv analysis for January, suggesting that the regulating model of environmental factors on phytoplankton assemblages is more complex when a bloom is about to start. When it was analyzed from the perspective of the entire investigation period, the phytoplankton community was mainly influenced by DIN/DIP, temperature, DSi, DIP and DSi/DIP ( Figure 13 ; Table 3 ). In fact, phytoplankton succession is largely determined by the interactions of physical and chemical factors (Sommer and Maciej Gliwicz, 1986). T. rotula, G. striata and P. globosa were distributed in the stations with high DIN/DIP ( Figure 13 ), indicating that these three species may have competitive advantages under P-limited conditions. T. rotula can grow well when the DIN/DIP ratio is 98 (Zhang, 2021), and P. globosa has been found to have competitive advantages when DIP is limited (Wang and Tang, 2006). All Chaetoceros species showed a strong negative relationship with inorganic nutrients, meaning that they were major consumers of these nutrients. Our findings are in agreement with earlier observations showing that hydrological and chemical factors play significant roles in determining the distribution and succession of phytoplankton (Peng et al., 2012; Cui et al., 2018).

Most studies have shown that P. globosa blooms occur during late autumn and early spring, both in China (Guo et al., 2007) and in other countries (Baudoux et al., 2006; Seuront et al., 2006). Previous studies have shown that P. globosa blooms in the coastal waters of Guangxi often occur in autumn or winter, with occurrence having been recorded from November to March (Luo et al., 2016; He et al., 2019a; He et al., 2019b). According to statistics performed by He et al. (2019a), the water temperature range during blooms caused by P. globosa in the coastal of China was 15-27°C, and low and stable temperatures were thought to be favorable for colony formation. At the end of November 2017, a P. globosa bloom began with an average water temperature of 16.9°C in a similar survey area in Qinzhou Bay. The bloom lasted for one and a half months and then dissipated in mid-January 2018 when the average temperature was 15.4°C (He et al., 2019a). The average water temperature was 19.0°C during two P. globosa blooms that occurred in the coastal area of Guangxi (Li et al., 2015). P. globosa was correlated with low temperature according to the RDA results ( Figure 13 ), and it seems that temperature was the key factor regulating the growth of P. globosa in Qinzhou Bay during this study. Previous studies have shown that P. globosa is always one of the predominant species in Qinzhou Bay (Jiang et al., 2012; Xu et al., 2014; Liu et al., 2020). In the present study, P. globosa was found to be a predominant species in Qinzhou Bay from October to February ( Figure 9 ). The bloom occurrence threshold was defined according to the cell size of HABs causative species, and a cell density of 1×10⁷ cells/L was considered to be the bloom occurrence threshold of P. globosa or a corresponding Chl a content of 10 µg/L (Guo et al., 2014). Although the P. globosa density increased sharply and was close to the bloom density in January 2019, it failed to form a bloom, and few colonies were found in the water ( Figure 11 ). During a typical P. globosa bloom, colonies dominate the phytoplankton community with high biomass (Chen et al., 2002; Hai et al., 2010), and the solitary cells were occasionally found or almost nonexistent in these blooms (Rousseau et al., 1994). During a P. globosa bloom occurred in Pearl River Estuary, the average colony density was 675 colonies/L with the highest value of 2000.2 colonies/L, and the diameter of visible colonies ranged from 0.5 to 2.5 cm (Wang et al., 2010b). The relationship between colony formation and density of solitary cells was revealed by Liu et al. (2021). In their study, P. globosa colonies were observed only when the density of solitary cells reached 10⁴ cells cells/mL. Here, although P. globosa abundance increased during January in the study region, cell densities remained below this 10⁴ cells/mL threshold. Even so, colonies were observed in the field, indicating that the factors controlling colony formation are more complex than replicated in laboratory studies.

Qinzhou Bay is mainly affected by northeast monsoon in winter (Gao and Chen, 2014). The date of sampling in January is 22^nd, and this area began to be influenced by strong northeasterly winds two days before the sampling. Three days of strong northeasterly winds resulted in a significant drop in temperature ( Figure 2A ). The average of salinity in January was lower than in December, and the salinity difference between different stations became smaller ( Figures 2B , S2 ). In particular, stations 01-03 near the inner bay had significantly lower salinities than observed at the outer stations in other months, but were all higher than 27 in January. These changes in salinity may imply that the influence of strong northeasterly winds may result in a mixing/upwelling where the outer bay with high salinity waters were mixed with inner bay low salinity nutrient enriched waters. The characteristics of the resultant new water mass (e.g., lower temperature, elevated nutrient concentrations, perhaps additional planktonic components) that drove the increase in P. globosa, and shaped the composition of the phytoplankton assemblage in January. It was also observed that all three surges in colony density occurred around the third day of strong winds (Xu, 2022). These findings suggest there is a correlation between the occurrence of P. globosa blooms and meteorological conditions.

Studies have shown that P plays a key role in the growth and colony formation of Phaeocystis; the growth and life-history transformation of Phaeocystis are inhibited under P deficiency, under which colonies fail to form (Veldhuis et al., 1991; Riegman et al., 1992). A concentration of 0.5-1 μmol/L phosphate is thought to be necessary for the formation of Phaeocystis colonies (Cariou, 1991). The average DIP in this survey was 0.36 ± 0.33 μmol/L, and the DIP was lower than the threshold for colony formation in P. globosa, especially in January. The DIP at all stations in the study area fell below the absolute P threshold for the growth of phytoplankton (0.1 µmol/L, Peery and Eppley, 1981), so there was not enough P for the cells to continue to grow and form a bloom. Many phytoplankton exploit P from organic sources by synthesizing alkaline phosphatase when DIP is depleted, and P. globosa was found to have the ability to utilize DOP under P limitation (Caroline et al., 2015; Qin et al., 2018; Huang et al., 2021). The DOP also declined to 0.15 ± 0.11 μmol/L in January during this study (data not shown). Mo et al. (2022) conducted nutrient enrichment experiments and confirmed that the addition of P to the surface seawater of Qinzhou Bay could stimulate the colony density, cell number and diameter of colonies. Therefore, phosphorus discharge management should be strengthened to avoid the stimulating effect of sudden P pollution on the growth of P. globosa.

Based on the analysis described above, it was speculated that in January, the mixing event that resulted from the strong northeasterly winds created improved environmental conditions for the growth of P. globosa, but that phosphorus was exhausted before higher cell densities could form. Therefore, the blooms failed to develop due to phosphorus limitation. As phosphorus is indispensable element for the formation of P. globosa colonies, it is speculated that the control of phosphorus input is of great significance for the prevention and control of P. globosa blooms and the reduction in disasters caused by its large colonies.