Impact of Dinophysis acuminata Feeding Mesodinium rubrum on Nutrient Dynamics and Bacterial Composition in a Microcosm

The development of Dinophysis populations, producers of diarrhetic shellfish toxins, has been attributed to both abiotic (e.g., water column stratification) and biotic (prey availability) factors. An important process to consider is mixotrophy of the Dinophysis species, which is an intensive feeding of the Mesodinium species for nutrients and a benefit from kleptochloroplasts. During the feeding process, the nutritional status in the environment changes due to the preference of Mesodinium and/or Dinophysis for different nutrients, prey cell debris generated by sloppy feeding, and their degradation by micro-organisms changes. However, there is little knowledge about the role of the bacterial community during the co-occurrence of Mesodinium and Dinophysis and how they directly or indirectly interact with the mixotrophs. In this study, laboratory experiments were performed to characterize the environmental changes including those of the prey present, the bacterial communities, and the ambient dissolved nutrients during the co-occurrence of Mesodinium rubrum and Dinophysis acuminata. The results showed that, during the incubation of the ciliate prey Mesodinium with its predator Dinophysis, available dissolved nitrogen significantly shifted from nitrate to ammonium especially when the population of M. rubrum decayed. Growth phases of Dinophysis and Mesodinium greatly affected the structure and composition of the bacterial community. These changes could be mainly explained by both the changes of the nutrient status and the activity of Dinophysis cells. Dinophysis feeding activity also accelerated the decline of M. rubrum and contamination of cultures with okadaic acid, dinophysistoxin-1, and pectenotoxin-2, but their influence on the prokaryotic communities was limited to the rare taxa (<0.1%) fraction. This suggests that the interaction between D. acuminata and bacteria is species-specific and takes place intracellularly or in the phycosphere. Moreover, a majority of the dominant bacterial taxa in our cultures may also exhibit a metabolic flexibility and, thus, be unaffected taxonomically by changes within the Mesodinium-Dinophysis culture system.


Introduction
The cosmopolitan dinoflagellate species Dinophysis acuminata is responsible for diarrhetic shellfish poisoning (DSP) events all around the world [1,2]. Okadaic acid (OA) and its derivatives known as dinophysistoxins (DTXs) and/or pectenotoxins (PTXs) are the dominant components in the toxin profile of D. acuminata. As strong inhibitors of serine and threonine protein phosphatases in eukaryotic organisms, OA and DTXs are capable of promoting potent tumors [3], inducing typical diarrhetic symptoms [2,4], and even acting as lethal agents to mammals. Recent transcriptomics analysis also revealed that OA and DTX-1 may induce hypoxia-related pathways or processes, unfolded protein response (UPR), and endoplasmic reticulum (ER) stress [5]. PTXs are generally not responsible for unpleasant gastrointestinal symptoms but are potentially involved in acute toxicity [6].
D. acuminata is a mixotrophic species that primarily requires phototrophic metabolism and plastid retention for long-term maintenance in the laboratory [7][8][9]. The Dinophysis-Mesodinium-cryptophyte is so far the only known food chain for Dinophysis growth. Dinophysis blooms are very much related to the distribution and abundance of Mesodinium [10][11][12]. Therefore, the nutritional status of prey and the surrounding environment may have a critical impact on the growth and toxin production of Dinophysis [13][14][15][16][17]. The feeding process of the latter involves not only the direct uptake of the prey organelles through a feeding peduncle (myzocytosis) and secretion of mucus traps but also the intense lysis of the ciliate cells [18][19][20]. Cell debris and organic substances originating from prey were reported to induce the DSP toxin release from Dinophysis [21]. The suspected harmful compounds (e.g., free polyunsaturated fatty acids) were not the shellfish toxins [22]. Additionally, "sloppy feeding" behavior generates a substantial amount of dissolved and particulate materials in the surrounding environment. This pool of biological organic matter combined with the extracellular toxin fraction may also function as a source of nutrients available to the heterotrophic bacterial community and, in turn, for Dinophysis cells after regeneration [21,23,24] or other biochemical pathways [25]. However, few studies have been conducted to assess the contribution and availability of these nutritional sources.
The role of algal-bacterial interactions during harmful algal bloom (HAB) has received attention in recent years [26][27][28]. The supply of dissolved organic substances through cell exudation or cell lysis is hypothesized to be a major interaction between phytoplankton and the associated bacterial community [24,25,29]. The influence of bacteria on the toxigenic properties of photosynthetic microalgae (mainly Alexandrium spp. producing paralytic shellfish toxins) has been widely examined (Reference [26] and literature therein). The "obligate" relationship between bacteria and mixotrophic Dinophysis species has been explored in terms of cell abundance and carbon equivalents, which show a possible dependence on bacteria-produced vitamin B 12 and, to a lesser extent, the potential of bacterivory for Dinophysis growth [23], which was otherwise confirmed in the case of Mesodinium rubrum [30]. Recently [31], the cluster of Alteromonadales have been identified as the unique prokaryotic microbiome associated with D. acuminata blooms in Northport Harbor, New York. This finding highlighted the importance of biogeochemical conditions in shaping the microbial consortia.
Mixotrophs may become the major players in an aquatic ecosystem due to their substantial contribution to the energy cycles and to nutrient cycles where heterotrophic bacteria control most of the pathways [32,33]. However, more compelling evidence is needed to explain the interactions between specific heterotrophic bacteria and nutrient dynamics mediated by the mixotrophy of Dinophysis species. Therefore, in this study, we focused on the bacterial community associated with the mixotrophic D. acuminata feeding on the mixotrophic M. rubrum in laboratory culture conditions. By tracking the changes of the bacterial assemblages and the nutritional status of the culture medium, we aimed to (i) study the nutrient dynamics mediated by mixotrophy through the different growth phases of Dinophysis feeding Mesodinium and the possible consequences for the ambient microbial community and (ii) identify the prevailing interactions among deterministic factors, the bacterial community, and DSP toxin dynamics during the feeding process. In the context of its mixotrophic nature, we hypothesized that an ingestion-derived nutrient shift combined with the activities of the toxin-producing Dinophysis could lead to niche separation of the microbial community.

Predator-Prey Population Dynamics and Environmental Changes
The simulation started with an initial density of M. rubrum of 6740 ± 1379.3 (mean ± SD, Group A) cells mL −1 and 7902 ± 373.0 (Group B) cells mL −1 (Figure 1). The M. rubrum population developed until the 5th day when Dinophysis cells were inoculated (Group A). The ciliate cell-density gradually declined under an average ingestion rate of 3.25 ± 0.38 prey cells predator −1 day −1 during 16 days while the Mesodinium population in the control group (Group B) doubled in 4 days and then significantly declined to 4833 ± 378.6 cell mL −1 from days 5 to 16. Dinophysis exponential growth lasted 15 days (from T1 to T3) and reached a maximal density of 1302 ± 282.1 cells mL −1 on the 20th day and remained in a plateau phase thereafter ( Figure 1). The bacterial abundance also changed during the feeding process (Group A) and the growth of M. rubrum (Group B).
Environmental characteristics varied differently over the course of the growth curve (  (Table 1 and Figure 1). A separate trend of POC, however, was observed for T3 when this compound continued declining in Group B but gradually accumulated in Group A. In addition, toxins accumulated in the culture medium as the Dinophysis population increased at the expense of Mesodinium in Group A. OA and DTX1 were presented together in the summary and in the following analysis due to the fact that they share the same chemical backbone (Table 1). Total OA + DTX1 and PTX2 contents in the culture medium were shown to increase over the growth curve and reached a concentration of 4109.58 ± 621.79 pg mL −1 and 26.74 ± 0.73 ng mL −1 , respectively, by the end of the experiment (Table 1).

Composition and Structure of the Microbial Community throughout the Growth Curve
Bacteria samples were also harvested six times throughout the entire growth curve. The composition and structure of the microbial community was demonstrated at the class and order levels, according to the DNA results ( Figure 2), and analyzed statistically using unweighted (structure) and weighted (composition) NMDS ( Figure 3) and UniFrac dissimilarity ( Figure 4).
During the molecular analysis, a total of 93 OTUs were observed throughout all the samples after being rarefied to an even depth of 25,396 reads (Table 2). A Good's coverage index of over 0.999 and the high validity of clean tags indicated that the sequencing had covered almost all the species in the samples and the results were convincing. Alpha diversity indexes (Simpson, Shannon Wiener, and Chao1) indicated that there were no significant differences (ANOVA, p > 0.05) of the bacterial community between the two treatments in each crucial time period (Table 2). Proteobacteria (relative abundance = 74.3%) and Bacteroidetes (relative abundance = 21.1%) were the two dominant bacterial phyla in all samples. The majority of Proteobacteria were Alphaproteobacteria (97.8%) and a small fraction of Gammaproteobacteria (2.1%). Bacteroidetes and Sphingobacteria only attributed to the Sphingobacteriales accounted for 14.2% of the total microbial assemblage (Figure 2a). At the order level, Rhodobacterales and Cellvibrionales, which are representative of Alphaproteobacteria and Gammaproteobacteria, respectively, were dominant ( Figure 2b).       The same results were also found in UniFrac dissimilarities analysis, which showed that the composition and structure of the bacterial community significantly changed during Dinophysis growth ( Figure 4a) and M. rubrum decay (Figure 4b) by using the dissimilarity indexes and the generalized, unweighted, and weighted UniFrac distances. As for the difference between the two treatments (Figure 4c), only a minor increase in an unweighted distance analysis was observed and an opposite trend was observed when the abundance of bacterial taxa was considered (generalized and weighted UniFrac distance), which suggests that the presence and proliferation of the Dinophysis population or the consequences of its association with Mesodinium decay may not influence the composition and structure of the dominant bacterial species.
To better identify the interactions between the microbiome and the biotic or abiotic factors characterized in the growth curve, we manually defined the bacteria into three groups, which include the abundant (Ab) group, the moderately abundant (M) group, and the rare (R) taxa group. These groups stand for a relative abundance above 1%, between 0.1% and 1%, and below 0.1%, respectively. In summary, out of the 93 identified OTUs, 5 OTUs were assigned to the abundant (Ab) group, 83 OTUs met the criterion of rare taxa (R), and 5 OTUs belonged to the moderate (M) group. The response of the microbial communities to the environmental changes was further investigated by the Mantel test (Table 3) and the most related factors were selected by the BIOENV procedure ( Table 4). The results showed that all three assemblages were significantly correlated (p < 0.01) to the environmental matrix (Table 1) especially the M and R (r > 0.5), which suggests that the selected parameters had a better interpretation on a relatively lower abundance of bacterial assemblages. Cell density of M. rubrum and D. acuminata were also included in the "environmental matrix" (Table 1) for interpretation through the BIOENV procedure ( Table 4). Out of the 11 parameters, eight were finally selected by BIOENV for the best correlation models to interpret the Bray-Curtis distance matrices of bacterial communities. They were M. rubrum density, Dinophysis cell density, DOC, PO4 3− , NH4 + , POP, OA + DTX1, and PTX2. Similarly, the selected parameters were more powerful in representing the M and R taxa (Table 4). The response of the microbial communities to the environmental changes was further investigated by the Mantel test (Table 3) and the most related factors were selected by the BIOENV procedure ( Table 4). The results showed that all three assemblages were significantly correlated (p < 0.01) to the environmental matrix (Table 1) especially the M and R (r > 0.5), which suggests that the selected parameters had a better interpretation on a relatively lower abundance of bacterial assemblages. Cell density of M. rubrum and D. acuminata were also included in the "environmental matrix" (Table 1) for interpretation through the BIOENV procedure ( Table 4). Out of the 11 parameters, eight were finally selected by BIOENV for the best correlation models to interpret the Bray-Curtis distance matrices of bacterial communities. They were M. rubrum density, Dinophysis cell density, DOC, PO 4 3− , NH 4 + , POP, OA + DTX1, and PTX2. Similarly, the selected parameters were more powerful in representing the M and R taxa (       Furthermore, Bray-Curtis dissimilarities were used to interpret the differences of the three assemblages of microbial communities during the growth phases of Dinophysis (Figure 5a) and Mesodinium ( Figure 6a) and with the surrounding environment in a mixed culture (Figure 5b-i) as well as in a control treatment (Figure 6b-f). In detail, during Dinophysis growth (Group A), the dissimilarity between the M and R taxa varied significantly (p < 0.01) since the culture aged while the abundant taxa showed no differences ( Figure 5a). As for the environmental factors, the M and R taxa also showed a more active response than the Ab taxa. These two assemblages positively correlated     (Table 1), the significant decrease in distance of the 0.1-1% fraction (Figure 7a, r = −0.567, p = 0.0003) may be mostly attributed to the approaching environmental variables. Therefore, the previously detected deleterious effect of Dinophysis activity (Figure 5g and 5h) may be limited to the R taxa fraction in the culture medium. The differences of the selected factors did not appear to be responsible for the differences in the bacterial community (Figures 7c-f) except for the M. rubrum cell density effect with which the M taxa significantly varied (Figure 7b). The impact of M. rubrum cells on the M taxa was partly responsible for the decrease in dissimilarity over culture days (Figure 7a) and the largest differences of ciliate density were actually  (Table 1), the significant decrease in distance of the 0.1-1% fraction (Figure 7a, r = −0.567, p = 0.0003) may be mostly attributed to the approaching environmental variables. Therefore, the previously detected deleterious effect of Dinophysis activity (Figure 5g,h) may be limited to the R taxa fraction in the culture medium. The differences of the selected factors did not appear to be responsible for the differences in the bacterial community (Figure 7c-f) except for the M. rubrum cell density effect with which the M taxa significantly varied (Figure 7b). The impact of M. rubrum cells on the M taxa was partly responsible for the decrease in dissimilarity over culture days (Figure 7a) and the largest differences of ciliate density were actually observed in the early phases of the Mesodinium growth curve (Figure 1). In summary, dissimilarities of the R taxa between the two groups ( Figure 7) were generally higher than those of M and Ab taxas compared to the results within each group (Figures 5 and 6), which leads to our assumption that the influence of Dinophysis and the toxins it produces may affect those R bacterial taxa.
Toxins 2018, 10, x FOR PEER REVIEW 13 of 21 observed in the early phases of the Mesodinium growth curve (Figure 1). In summary, dissimilarities of the R taxa between the two groups ( Figure 7) were generally higher than those of M and Ab taxas compared to the results within each group ( Figures 5 and 6), which leads to our assumption that the influence of Dinophysis and the toxins it produces may affect those R bacterial taxa.

Discussion
The algal-bacterial interaction has been shown to have a critical effect on bloom dynamics, but much of the research on this subject has been limited to photosynthetic algal species [34][35][36]. Due to their mixotrophic nature, Dinophysis species need to ingest Mesodinium, which is the sole determined prey [37,38], for population growth. Field populations of the ciliate usually aggregate in the

Discussion
The algal-bacterial interaction has been shown to have a critical effect on bloom dynamics, but much of the research on this subject has been limited to photosynthetic algal species [34][35][36]. Due to their mixotrophic nature, Dinophysis species need to ingest Mesodinium, which is the sole determined prey [37,38], for population growth. Field populations of the ciliate usually aggregate in the subsurface water layers and perform diurnal vertical migration [39]. At the same time, Dinophysis species normally represent a small proportion of the phytoplankton community with a minor importance in relation to major biogeochemical cycles [25,40,41]. Thus, Dinophysis blooms and the time period from Dinophysis initiation to the M. rubrum decline is difficult to capture [10]. In the current study, we simulated an ideal predator-prey microcosms study to estimate the biogeochemical consequences of their co-occurrence and identify the ecological niche of bacterial communities developed during Dinophysis growth and Mesodinium decline.
We found that the M. rubrum-Dinophysis system intensively altered the biogeochemical status of the culture medium (Table 1). Dinophysis feeding activities accelerated the decline of the M. rubrum population (Figure 1). The differences in microbial communities between Groups A and B were mainly ascribed to the presence of Dinophysis cells (up to 1300 cells per mL −1 ) and the accumulation of toxins or other compounds produced by Dinophysis. In addition, considering that the influence of Dinophysis was restricted to the taxa in relatively low abundance, the dissimilarity may mainly consist of intracellular or phycosphere bacteria given that Dinophysis cells could function as particulate vectors or as hosts of certain bacterial species [42]. This assumption was also shown by the increased particulate organic carbon in this study (Table 1) as well as previous experiments [14]. Furthermore, in a study on bacterial assemblages associated with Dinophysis, it was also found that the addition of Dinophysis culture filtrate caused no significant changes in their relative abundance while the prokaryotic genera directly associated with Dinophysis were found in the >20-µm size fraction [31]. Species-specific grazing or deleterious effects of Dinophysis cells may also gradually shape the structure of the bacterial community. In a recent study [23], the potential contribution of bacterivory and bacterial remineralization to the growth of Dinophysis was calculated. It was concluded that neither of these processes were quantitatively relevant in order to support the increased biomass observed in the study.
Changes in nutrients and carbon content were observed in the M. rubrum-D. acuminata culture (Group A) ( Figure 1 and Table 1), but the pattern was quite similar to of M. rubrum population dynamic itself (Group B, the control group), which indicates that ingestion by Dinophysis was unable to integrate most of the nutrients and carbon compounds derived from the M. rubrum population. Enhancement of secondary metabolite (DSP toxins) production was verified in Dinophysis cells [21] and in heterotrophic microbes [32] when exposed to M. rubrum living cells as well as cell lysate, but no direct evidence has emerged to decipher the associations among Dinophysis, microbes, and detritus in terms of nutrient cycling [37,43]. Tong et al. [14] estimated that D. acuminata was able to assimilate 65% and 25% of the particulate nitrogen and phosphate through predation, respectively. In the current study, the ingestion rate (ca. 3.25 prey cells predator −1 day −1 ) approached the highest level of the calculated growth rate where the predator: prey ratio was considered saturated [44]. The ingestion rate may not able to reach a higher value due to the low growth rate of the prey. However, the amount of ingested carbon calculated from a recently published paper [15] was about 2466.7 pg C cell −1 d −1 , which is already beyond Dinophysis needs to maintain growth [8]. Field studies also found that the Dinophysis population preferred savaging on their prey in a short period of time [10,45].
After the M. rubrum population collapsed, the nutritional status of the culture medium was remarkably affected, which is expressed by the elevated levels of prey ammonium and fluctuant DOC and forms a hotspot of biogeochemical activities for the Dinophysis [40,46] and also the heterotrophic bacteria therein [47][48][49][50]. In this case, ammonium became the key component in the nitrogen cycle Mesodinium decline. NH 4 + concentration remained constant until the addition of Dinophysis (T1, Table 1). A sharp increase was noticed due to a possible consequence of bacterial ammonification when M. rubrum decayed (Figure 1). Then NH 4 + was significantly low at the late plateau phase of Dinophysis (T4 and T5, Table 1, Mann-Whitney Rank Sum Test, p = 0.041) when compared to the Mesodinium control. This finding is consistent with previous results showing that ammonium may lead to an increase of field populations of Dinophysis [51] and could be assimilated, which enhances the growth of Dinophysis at certain levels [23]. As for the other forms of nitrogen, nitrate with a moderate level (~50 µM) was dominant in our initial culture system but was used up in the first five days (Table 1, phase T1). When a high concentration of nitrate, e.g., up to 200 µM, is available in the culture medium, the availability of ammonium at lower concentrations (<2 µM) may be undetected or compensated [14]. The decay of the M. rubrum population seems unrelated to the nutrient limitation given the availability of both nitrogen (ammonium and nitrate) and phosphate in the culture medium [32,37]. Characterized as a plastidic-specific non-constitutive mixotroph [52], Dinophysis cells mainly retain chloroplasts from their prey and perform photosynthesis by those kleptoplastids as a carbon source. At the same time, marine heterotrophic bacteria play a major role in incorporating, respiring, and degrading dissolved organic carbon. However, changes of DOC concentration in this study hardly demonstrated the dissimilarities of the microbial community ( Figures 5-7), which suggests that the heterotrophic bacteria assemblages in our culture system may exhibit a metabolic versatility at least within the range of our DOC measurement. Moreover, bacterial communities at large phylogenetic group levels may exhibit general outcomes when exposed to high DOC concentrations [24]. Therefore, variation of DOC during phytoplankton dynamics may not taxonomically drive the shift of the major microbial community. Bacterial assemblages of M taxa (>0.1% and <1% of relative abundance) seem to be the most sensitive portion to changing the environmental nutrient conditions (Figures 5-7). High-throughput sequencing revealed that this portion is composed of bacterial taxa assigned to the Rhodospirillaceae, Cytophagaceae, Flavobacteriaceae, and CHAB-XI-27 at family-level resolution. According to the 16s rDNA sequencing, Proteobacteria (Alphaproteobacteria-Rhodobacterales and Gammaproteobacteria-Cellvibrionales dominated) and Bacteroidetes (Sphingobacteria-Sphingobacteriales dominated) constituted more than 90% of the relative abundance of the microbial community in the culture medium cumulatively. These results are not surprising since only a limited number of heterotrophic bacterial lineages dominate those eukaryotic phytoplankton-associated communities [48]. Furthermore, these bacteria lineages cover those groups responsible for both monomer (such as amino acids) and polymer (such as chitin and protein) degradation in the ocean [53]. Less than 100 bacterial OTUs were assigned in our culture system, which is far less than previous field studies [45,54]. Bacterial community results from laboratory cultures show a lower diversity of bacterial assemblages because long-term maintenance may have eliminated those species that had already been overwhelmed and laboratory studies could avoid invasion of accidental species that are common in field studies. This is also the reason for which some studies exclude those rare species in their analysis [54]. The only study addressing the interaction between a Dinophysis bloom and the microbial community revealed that, even during the peak of a Dinophysis acuminata bloom (cell density~1300 cell mL −1 ), the Dinophysis cells only accounted for 29% of the phytoplankton community [31]. Thus, the comparison of the bacterial community between our two groups could merely be attributed to the influence of the D. acuminata population. Our finding that only a low abundance of bacterial species was altered during Dinophysis intensive feeding activities indicated that, even though intense mixotrophy could remarkably drive biogeochemical dynamics, changes of the phytoplankton population may not be reflected by changes in those abundant bacterial phylotypes or in metabolic generalists. Moreover, connections may exist between specific species of bacteria and Dinophysis cells. Locating these connections by using metabolite and meta-transcriptome analysis may give us a further understanding for how these organisms interact with each other. Establishing an axenic culture of Dinophysis and comparing physiologies with non-axenic cultures over long-term periods may offer more robust evidence of the dependence of Dinophysis on bacteria [25]. However, attempts to purify Dinophysis and M. ruburm cells by using antibiotic treatment were not successful either in our laboratory or after efforts devoted by other groups [31]. Novel approaches to generate axenic algae cultures have been tested on the cyst-forming species Gymnodinium catenatum starting from resting cysts [55] and freshwater species by using fluorescence-activated cell sorting [56]. Yet, the methods proposed may not be universal and transferable considering that Dinophysis and M. rubrum cells are oddly shaped and fragile. More efforts are needed in the future to come up with appropriate approaches to initiate axenic Dinophysis cultures or elegant methods to directly target specific interactions between Dinophysis and associated bacteria.

Cultures
A unicellular algal culture of D. acuminata (DAYS01) was established from cells previously isolated from Xiaoping Island (121.53 • E 38.83 • N), the Yellow Sea, China in July 2014 [13]. The ciliate M. rubrum (AND-A0711) and the cryptophyte, Teleaulax amphioxeia (AND-A0710) were isolated from coastal waters off Huelva, Southern Spain in 2007 [57]. All cultures were routinely inoculated based on the cryptophyte-M. rubrum-Dinophysis food chain [7,13] in f/6-Si medium, which was prepared with 1/3 nitrate, 1/3 phosphate, 1/3 metals, and 1/3 of the vitamins concentrations in the f/2-Si medium. Cultures were maintained at 15 • C under a light intensity of 3000 lux and a 14 h light:10 h dark photo cycle. Dinophysis and/or the ciliate subsamples were taken every 2 or 3 days and fixed with 3% (v/v) formalin solution for microscopic enumeration in a Sedgewick-Rafter counting chamber at 100× magnification. For bacterial analysis, the formalin-preserved samples (1 mL) were stained with 2 µL of 4 ,6-diami-dino-2-phenylindole (DAPI) solution (1 mg mL −1 ) and filtered onto a black polycarbonate filter (pore size: 0.22 µm, diameter: 25 mm, Millipore, Burlington, MA, USA). Then the filters were gently removed onto a glass slide and observed at 600× by using fluorescence microscopy (DMi8, Leica Microsystems, Buffalo Grove, IL, USA) under UV excitation.

Nutrient Sample Collection and Preparation
The growth curve was manually divided into six different growth phases. Six sampling spots were set up to collect nutrient and toxin samples, which are, hereafter, referred to as T0-at the very beginning of the incubation, T1-early phase following inoculation of Dinophysis, T2-the middle of the exponential growth of Dinophysis, T3-the end of the exponential growth of Dinophysis, T4-the depletion of M. rubrum, and T5-the end of Dinophysis growth.
For nutrients, 30-mL culture medium were filtered through pre-combusted GF/F filters (25 mm, Whatman, Maidstone, UK) for particulate organic carbon (POC) and particulate organic phosphate (POP) collection, respectively. The filters for POC were dried in a 60 • C oven for 24 h, stored at −20 • C, and analyzed on an elemental Analyzer (EA3000, EuroVector S.p.A, Milan, Italy). The particulate phosphate was converted to orthophosphate (PO 4 3− ) by first hydrolyzing it by the addition of 5 mL of 5% potassium persulfate and 10 mL of Milli-Q water and then autoclaving it (121 • C) for 20 min.
The filtrate was used for quantifying the dissolved inorganic nutrients (DIN and DIP) and dissolved organic carbon (DOC). Nitrate, ammonium, and phosphate were analyzed by SKALAR SAN ++ Autoanalyser (SKALAR, Breda, The Netherlands). DOC concentration was determined by using a TOC Analyzer (Multi C/N 3100, Analytik Jena, Jena, Germany). All analyses were conducted following the protocols of the manufacturers.

Toxin Analysis
Dinophysis cells and culture medium were separated for toxin analysis. Between 10 and 30 mL of culture medium was harvested and analyzed for toxins since the inoculation day (T1) of Dinophysis cells. The medium was kept in 50-mL centrifuge tubes and stored at −20 • C before extraction. Solid-phase extraction (SPE) was employed [13,58] for the extraction of cells or medium samples. The SPE column (Oasis HLB 60 mg, Waters, Milford, MA, USA) was preconditioned with 6 mL of methanol and 6 mL of Milli-Q water. Once the cells or medium samples were loaded, the cartridge was washed with 3 mL of Milli-Q water and then blow-dried and eluted with 1 mL of methanol to collect the toxins into an HPLC vial. Eluates from the samples were then heated at 40 • C in a heating block (HP-S016SY), dried under a stream of N 2 , and re-suspended in 1 mL of 100% methanol for toxin analysis.
Toxin analysis was performed on an UltiMate 3000 LC (Thermo Scientific TM Dionex TM , Waltham, MA, USA) and an AB 4000 mass spectrometer system (SCIEX, Framingham, MA, USA) with electrospray ionization. PTX2 was analyzed in positive mode, while OA and DTX1 were analyzed in negative mode. Chromatographic separation was performed by using a Waters XBridgeTM C18 column (3.0 × 150 mm, 3.5-µm particle size) (Milford, MA, USA) at 40 • C for a negative mode. The mobile phase consisted of phase A, 0.05 v/v % ammonia hydroxide in water, and phase B, 0.05 v/v % ammonia hydroxide in 90% acetonitrile with a flow rate of 0.4 mL min −1 and 10 µL injection. A linear gradient elution from 10% to 90% B was run for 9 min, held for 3 min at 90% B, decreased to 10% B for 2 min, and held at 10% B for 4 min to equilibrate at the initial conditions before the next run. In a positive mode, a Waters XBridgeTM C18 column (2.1 × 50 mm, 2.5-µm particle size) at 25 • C was used for chromatographic separation. A linear gradient from 10% to 80% acetonitrile containing a constant concentration of buffer (2 mM ammonium formate and 50 mM formic acid) was run between 0 min and 9 min and held at 80% acetonitrile for 2 min at a flow rate of 0.3 mL min −1 . Standards for OA, DTX1, and PTX2 were purchased from the National Research Council, Canada.

DNA Extraction and Illumine Sequencing
At each time spot, a 120 to 150 mL culture medium was retrieved and filtered by using a 0.22-µm cellulose filter. Then the filters were folded and stored at −80 • C. Total genome DNA was extracted by using GenJET Genomic DNA Purification Kits (Thermo Scientific, Waltham, WA, USA) following the protocols of the manufacturer. Bacterial amplicons were produced by targeting the 16S V3-V4 hypervariable region with universal primers 343F (5 -TACGGRAGGCAGCAG-3 ) and 798R (5 -AGGGTATCTAATCCT-3 ). The amplicon quality was visualized by using gel electrophoresis, which was purified with AMPure XP beads (Agencourt) and amplified for another round of PCR. After being purified with the AMPure XP beads again, the final amplicon was quantified by using Qubit dsDNA assay kits. Equal amounts of purified amplicons were pooled for subsequent sequencing. The amplicon libraries were then sequenced on an Illumina MiSeq platform (Shanghai OE Biotechnology Co., Ltd., Shanghai, China).

Bioinformatics Analysis
Paired-end reads were preprocessed by using Trimmomatic software [59] to detect and cut off the ambiguous bases (N), barcodes, primers, and low-quality sequences. After trimming, paired-end reads were assembled by using FLASH [60]. Valid tags were subjected to clustering to generate operational taxonomic units (OTUs) using VSEARCH software (Version 2.4.2) at a 97% similarity setting [61]. The representative sequence of each OTU was selected by using the QIIME package. Representative reads were annotated and blasted against the Silva database (Version 123) using the ribosome database project (RDP) classifier with a confidence threshold of 70%. Sequences were submitted to the NCBI Sequence Read Archive with the accession number SRR6048156.
Package phyloseq [62] in R (http://www.Rproject.org, v. 3.3.3) was used to perform alpha and beta diversity calculations and to visualize the results of dimensional reduction approaches (nonmetric multidimensional scaling, NMDS). UniFrac distance matrices (weighted, unweighted, and generalized) were calculated with the R package GUniFrac [63] based on the OTU table and the phylogenetic tree. Note that the OTU table was square-root transformed and the explanatory matrix was z-score-transformed in R before the statistic procedures. To characterize the potential functions of different bacteria groups, we defined "abundant" (Ab) and "rare" (R) OTUs with the criteria of the average relative abundance across all the samples above 1% and below 0.1%, respectively [64]. The rest of the bacterial lineages (>0.1% and <1%) were assigned to a "moderately abundant" (M) group. The BIOENV procedure was implemented to identify the subset of a set of explanatory variables (nutrients profile, toxin content, and mixotrophs density). The Euclidean distance matrix of this correlates maximally with the Bray-Curtis compositional dissimilarity matrix of the OTU table. The Mantel test was also run with 999 permutations to check whether the subset of explanatory variables was able to capture the variation of bacterial communities from the three fractions. The BIOENV procedure and Mantel test on the three matrices (abundant, moderately-abundant, and rare fraction) were achieved by using relevant functions in the R package vegan [65]. Bray-Curtis dissimilarities of abundant (Ab), moderately abundant (M), and rare taxa (R) among samples were plotted against differences of explanatory variables (subtraction between samples) selected via BIOENV.