The cyanobacterium Gloeotrichia echinulata increases the stability and network complexity of phytoplankton communities

Changes in the abundance of a taxon can have large effects on communities, particularly if that taxon is a strong interactor. These changes may arise as a consequence of environmental change, recruitment from dormant stages, or quirks of population dynamics, and have effects that ripple through a community interaction network. We hypothesized that cyanobacteria, which are increasing in many freshwater lakes globally, may be strong interactors because they can exert large and persistent effects on the biomass and composition of other phytoplankton. To test this hypothesis, we evaluated how the phytoplankton community responded to different densities of Gloeotrichia echinulata, a large colonial cyanobacterium increasingly observed in low-nutrient lakes in northeastern North America, in an in situ mesocosm experiment. We observed that many phytoplankton taxa, especially diatoms and green algae, responded primarily to increased nutrient availability (a result of Gloeotrichia's nitrogen fixation and translocation of phosphorus from the sediments), while a few taxa (two euglenophytes, one dinoflagellate, and one cyanobacterium) responded to both the direct and indirect effects of Gloeotrichia. Surprisingly, Gloeotrichia reduced the compositional variability of the phytoplankton community relative to the non-Gloeotrichia control treatment; there was no effect on the aggregate temporal variability of total non-Gloeotrichia biovolume. Moreover, experimentally increased densities of Gloeotrichia coincided with increasing complexity of the phytoplankton community in network analyses of taxon co-occurrences, as indicated by significantly greater network density and transitivity and shorter path lengths. Taken together, these findings suggest that Gloeotrichia may be a strongly interacting species in low-nutrient lakes, with the potential to increase the resilience of phytoplankton communities to future disturbance by increasing compositional stability and network complexity.


INTRODUCTION
Ecological communities can be described and visualized as networks of interacting populations (Dunne et al. 2002). Fluctuations in the relative abundance of certain taxa within communities can trigger rapid changes in population interactions and thus alter network structure, relative abundance, and aggregate properties such as total biomass (Micheli et al. 1999). For example, if the addition of a new strongly interacting species (sensu Paine 1980), or an increase in the total or relative abundance of a previously rare one, has a large impact on resource availability, it may change not only which other species can persist in the community, but also their abundance and interactions with one another (Paine 1980, Power et al. 1985, Soule et al. 2005. Increases in strong interactors may occur as a result of environmental change, recruitment from dormant stages, or changes in population dynamics, and have effects that ripple through the community network. Because disentangling the effects of strong interactors on their community requires experimental manipulation of species assemblages, most studies on the effects of these taxa have focused on macroscopic species. Less is known about the consequences of increases in potential strong interactors for microscopic communities, especially those that have the potential to become macroscopic problems from the perspective of the public, such as large surface aggregations of cyanobacteria in freshwater ecosystems (blooms). Cyanobacterial blooms have the potential to substantially affect the community structure of, and interactions among, other phytoplankton (Cottingham et al. 2015), and could be of increasing importance in the future, since their duration, magnitude, and geographic extent appear to be increasing globally (Paerl andHuisman 2008, Brookes andCarey 2011).
There are multiple mechanisms by which cyanobacteria may act as strongly interacting taxa and alter phytoplankton community structure and variability. Most previous work has focused on negative interactions due to the ecophysiological adaptations that allow cyanobacteria to outcompete other phytoplankton (reviewed by Carey et al. 2012b). For example, some cyanobacteria excrete toxins and allelopathic chemicals that inhibit the growth, division, and metabolism of target phytoplankton (Leflaive and Ten-Hage 2007), and some can create surface scums that decrease light availability for other phytoplankton (Reynolds 2006). However, cyanobacteria may also have positive effects on other phytoplankton: Some taxa are able to fix nitrogen (N) and some are able to acquire phosphorus (P) from sources that other taxa cannot access (reviewed by Cottingham et al. 2015). These "new" nutrients have the potential to benefit other phytoplankton because even live cyanobacterial cells tend to leak N and P in forms available for uptake (Cottingham et al. 2015). Moreover, cyanobacteria can also excrete stimulatory allelochemicals that increase growth rates and division of some phytoplankton (e.g., Suikkanen et al. 2005).
We tested the hypothesis that cyanobacteria exert strong and persistent effects on phytoplankton community structure and variability, even after blooms have ended, by analyzing data on phytoplankton communities before, during, and after experimental blooms in mesocosms. We chose the large, colonial cyanobacterium Gloeotrichia echinulata (J. E. Smith) P. Richter (hereafter, Gloeotrichia) as our focal species for two reasons. First, the large size of its colonies (1 to >3 mm in diameter) allows its densities to be directly manipulated. Second, it is of practical relevance, as reports of Gloeotrichia in oligotrophic and mesotrophic lakes across the northeastern USA and Canada are increasing (Carey et al. 2008, 2012a, Winter et al. 2011, and both ecologists and the general public are concerned about its potential effects on aquatic food webs and water quality. We used a suite of community analyses to identify the effects of Gloeotrichia on other phytoplankton taxa, focusing particularly on the effects mediated by altered N and P concentrations, and the concordance of results across metrics. These included differentiating the non-nutrient vs. nutrient effects of Gloeotrichia on other phytoplankton using path analysis ( Fig. 1) and calculating aggregate and compositional temporal variability (sensu Micheli et al. 1999) for the rest of the phytoplankton community. We also used network analysis to examine how taxon co-occurrences varied with Gloeotrichia density.
We predicted that the experimental Gloeotrichia blooms would have both nutrient and non-nutrient effects on other phytoplankton taxa and that these effects would increase variability, thereby decreasing stability, in phytoplankton community composition relative to no-bloom controls. Here, we define stability as the inverse of variability, following Grimm and Wissel (1997). These expectations were derived from our previous work suggesting that Gloeotrichia can exert both strong positive and negative effects on the biomass of other phytoplankton, depending on lake conditions (Carey et al. 2014a, b). On the positive side, Gloeotrichia is an N-fixer (Stewart et al. 1967) and takes up and stores large quantities of P in excess of its immediate needs (Istv anovics et al. 1993); these nutrients can be released to the water column in available forms, alleviating nutrient limitation (Pitois et al. 1997, Nõges et al. 2004). In the laboratory, Gloeotrichia stimulates growth of other phytoplankton by excreting stimulatory allelochemicals (Carey and Rengefors 2010). However, Gloeotrichia can also create surface scums that reduce light (Carey et al. 2014b). Consequently, we expected that both the aggregate and compositional temporal variability of the phytoplankton community during and after a Gloeotrichia bloom would be greater than in a no-bloom control and that Gloeotrichia would exert larger effects on other phytoplankton at higher densities, with effects not mediated by N and P being more important at higher densities. Finally, because Gloeotrichia can alter the biomass and composition of phytoplankton at the division level (Carey et al. 2014a), we expected that the phytoplankton community experiencing the highest density of Gloeotrichia blooms would exhibit greater network complexity than the no-bloom control.

In situ mesocosm experiment
As previously reported (Carey 2012, Carey et al. 2014a), we conducted an in situ mesocosm experiment in oligotrophic Lake Sunapee, New Hampshire, USA (43°24 0 N, 72°20 0 W), in July 2008 to examine the effects of different densities of Gloeotrichia on phytoplankton communities. The timing of this experiment was motivated by previous observations of Gloeotrichia scums in Lake Sunapee in early July (Carey et al. 2014c). Addition of Gloeotrichia to the mesocosms simulated the natural process of bloom formation and the translocation of N and P within the Gloeotrichia colonies from the sediments into the water column (Cottingham et al. 2015). Here, we summarize the major aspects of the experiment that pertain to this study; for complete details, see Appendix S1.
Just prior to the experiment, we collected Gloeotrichia colonies from Lake Sunapee, processed them into aliquots in the laboratory after a careful cleaning and inspection procedure, and stored them in incubators. Only the largest, buoyant colonies with all trichomes intact were used in the experiment. The day before the experiment began, we deployed 16 clear, 50 L mesocosms in a littoral cove of Lake Sunapee. The mesocosms were suspended from frames with the openings raised above the lake surface to prevent waves from overtopping the rims and kept covered with mesh to prevent insect, bird, or zooplankton immigration. The mesocosms were filled with unfiltered water from the upper 0.3 m of the lake, which included only small zooplankton (rotifers and nauplii) as most large zooplankton avoid visual predators during the day (Lampert 1989).
We randomly assigned the mesocosms to four Gloeotrichia density treatments: 0 (control), 25, 50, and 400 Gloeotrichia colonies/L, with four replicates each. These treatments reflect observed Gloeotrichia densities in low-nutrient lakes in the northeastern United States (0-450 colonies/L; Carey et al. 2012a); the highest littoral Gloeotrichia density observed in Lake Sunapee to date is 80 colonies/L (K. L. C., unpublished data).
The first morning of the experiment, we sampled each mesocosm to establish baseline conditions. Immediately after this initial sampling, we added Gloeotrichia colonies to the mesocosms to create the treatments, simulating the recruitment of different densities of colonies leaving the sediments and entering the water column to form blooms. We sampled each mesocosm again 24 h after Gloeotrichia addition and then every four days until day 13; at the end of the experiment, all Gloeotrichia had senesced and sunk to the bottom of the mesocosms. This pattern of senescence in the mesocosms over 13 days is similar to patterns of Gloeotrichia density observed in natural systems (Carey et al. 2008). On each sampling day, we used an integrated tube sampler to collect water for phytoplankton samples (250 mL, preserved with Lugol's iodine solution in opaque bottles), total N and P (100 mL, kept frozen until analysis), and soluble nutrients (100 mL, filtered through 0.7-lm Whatman GF/F filters (GE Healthcare Life Sciences, Marlborough, Massachusetts, USA) and frozen until analysis for nitrate, NO 3 À ; ammonium, NH 4 + ; and soluble reactive P, SRP). We also sampled Gloeotrichia by filtering 3-5 L of water through 80-lm mesh, preserving the collected sample, and returning the filtrate to the mesocosms.
We analyzed the nutrient and phytoplankton samples according to the standard protocols. Throughout the experiment, NO 3 À , NH 4 + , and SRP concentrations were consistently below method detection limits and thus cannot be used to evaluate the effects of Gloeotrichia on soluble nutrient pools (Appendix S1). At least 50 mL of each phytoplankton sample was concentrated and settled for 3 days before C.C.C. enumerated all phytoplankton cells, generally to genus, at 4009 magnification according to Uterm€ ohl (1958) on an inverted Nikon MSD microscope (Nikon Inc., Melville, New York, USA). Four phytoplankton genera (Aphanocapsa, Chroococcus, Merismopedia, and Trachelomonas) found in many samples were differentiated to sub-genus levels. We calculated non-Gloeotrichia phytoplankton biovolume (lm 3 /mL) by approximating the cells to geometric shapes. In addition, we counted Gloeotrichia colonies on an Olympus SZH10 dissecting microscope (Olympus America, Center Valley, Pennsylvania, USA).
We determined the effects of the Gloeotrichia treatments on the biovolume of each phytoplankton taxon using one-way repeated-measures ANOVA using SAS PROC MIXED (SAS v.9.2; SAS Institute, Cary, North Carolina, USA; Wolfinger and Chang 1999). We selected the best error structure from three different covariance matrices (compound symmetry, Huynh-Feldt, and unstructured) using the Bayesian Information Criterion.

Path analysis
To evaluate the relative strengths of the direct (non-nutrient) and indirect (nutrient-mediated) effects of Gloeotrichia on the biovolume of the taxa that were present in >5% of the observations, we fit saturated path analysis models to the four post-treatment sampling dates using PROC CALIS in SAS on the correlation matrices among variables. Specifically, we evaluated the direct effect of observed Gloeotrichia density and the indirect effect of total N and P (TN and TP) on each taxon (Fig. 1). To estimate the indirect effect, we used the first principal component from a principal components analysis of all observed TN and TP values, due to the correlation between TN and TP concentrations (r = 0.39) and our focus on the combined effects of increased nutrients rather than the separate effects of N and P. These path analyses were bidirectional: They were not examining causal effects of Gloeotrichia density and nutrients on phytoplankton taxa, but rather the strength of the interaction between each phytoplankton taxon and Gloeotrichia/nutrients because the phytoplankton taxa themselves may have altered nutrient concentrations. We fit these models only for the 29 taxa which had approximately linear relationships with the predictors and evaluated the statistical significance of the estimated effect coefficients using a ❖ www.esajournals.org Bonferroni-corrected a of 0.0009 (0.05 divided by 58). Because the models are saturated, we do not report traditional fit metrics here because they all suggest a very good fit to the data.

Aggregate and compositional variability
We tested the prediction that Gloeotrichia would increase variability in the phytoplankton community by first comparing aggregate temporal variability (hereafter, aggregate variability) among treatments by calculating the coefficient of variation (CV) of total non-Gloeotrichia biovolume in each mesocosm over time, then comparing CVs across treatments with one-way ANOVA. Second, to evaluate compositional temporal variability (hereafter, compositional variability), we used multivariate dispersion (Anderson 2006), which was originally proposed as a test for b-diversity in space or a measure of multivariate homogeneity of variance, but is also an effective test of temporal variability when applied to time series (Brown and Lawson 2010). We used the function betadisper in the vegan package in the R statistical environment (R Development Core Team 2015) to calculate the dispersion metric for each experimental unit, and one-way ANOVA with Tukey's post hoc multiple comparison tests to examine the differences in dispersion between the Gloeotrichia density treatments. To evaluate whether the observed treatment differences were driven more by changes in species presence or by changes in relative abundance, we performed the analysis using both the Jaccard and Bray-Curtis distance metrics.

Network analysis
We used network analysis to examine the statistical and structural characteristics of the community of phytoplankton genera (nodes) and the connections between them (edges) (Newman 2003). Network analysis is increasingly being used to explore interactions among taxa in microbial and phytoplankton communities (e.g., Barberan et al. 2012, Kara et al. 2012, Patrick et al. 2014), but experimental applications remain rare, especially for aquatic communities.
We constructed separate co-occurrence networks for the phytoplankton genera observed in each of the four Gloeotrichia density treatments on the four post-treatment sampling dates. Co-occurrence networks are the collective interconnections of nodes based on their paired presence within a community (Proulx et al. 2005). Given the small number of observations, we pooled all 16 observations collected within a treatment after Gloeotrichia addition (4 replicates per treatment 9 4 sample days) to create one observed network and 1000 simulated networks for each of the four Gloeotrichia density treatments. The simulated networks came from Monte Carlo simulations drawing 12 observations from the 16 total observations, and were used to compare structural network characteristics among the treatments. Within each of the 1000 simulations for a treatment, we used the number of co-occurrences of each pair of taxa during the 12 observations to calculate probabilities of each taxa pair's co-occurrence. We created 1000 simulated undirected Erd} os-R enyi networks based on the number of observed nodes, with the probability of undirected edge connections among them derived from the probabilities of cooccurring taxa pairs in each Monte Carlo simulation (following Barberan et al. 2012).
For the observed and simulated networks, we calculated structural properties that assessed network complexity. These characteristics included the number of edges and nodes (i.e., taxa richness), network density (the proportion of potential connections in a network that are actually connected), mean path length (the mean shortest path length between all nodes), and transitivity (also called clustering coefficient), a metric of how connected nodes are in the network (Wasserman andFaust 1994, Newman 2003). Complete transitivity, indicated by a value of 1, implies that if taxon x cooccurs with y, and y co-occurs with z, then x and z also co-occur. We compared these structural metrics among Gloeotrichia treatments for the simulated networks using permutational ANOVA in R. All network analyses were carried out in R using sna (Butts 2008) and igraph (Csardi and Nepusz 2006) packages. Networks were visualized with the software program Cytoscape v.2.8.3 (Smoot et al. 2011). After examining the structural properties of the networks, we also compared node degree, the number of edges connected to each taxon, among the taxa that were present in the phytoplankton community in all four Gloeotrichia density treatments. Taxa with a high node degree are referred to as "hubs" and play an important role in ecological networks because of their high connectedness to other taxa (e.g., Watson et al. 2011).

In situ mesocosm experiment
As previously reported (Carey et al. 2014a), there were no significant differences in Gloeotrichia density, non-Gloeotrichia phytoplankton biovolume, TN, or TP among treatments prior to Gloeotrichia addition (Appendix S1). However, additions of Gloeotrichia colonies increased Gloeotrichia density, TN, and TP (Appendix S1: Table S1, Fig. S1). Additions of 50 or 400 colonies of Gloeotrichia/L resulted in surface scums, which lasted for~12 d; after 13 d, the Gloeotrichia colonies senesced and sunk to the bottom of the mesocosms. The total biovolume of the phytoplankton community, excluding Gloeotrichia, increased in response to Gloeotrichia addition and exhibited different temporal dynamics among treatments (Appendix S1: Table S1). Mean non-Gloeotrichia phytoplankton biovolume increased steadily beginning one day after Gloeotrichia addition in the 50 and 400 Gloeotrichia colonies/L treatments, but was consistently low in the 0 and 25 colonies/L treatments (Appendix S1: Fig. S1).
At the genus and sub-genus level, there were substantial shifts in community composition in response to the Gloeotrichia treatments: 15 of the 59 identified taxa showed significant Gloeotrichia or Gloeotrichia 9 time interactions in the RM-ANOVA (Appendix S2: Table S1). For example, the diatoms Asterionella, Stauroneis, Stephanodiscus, Synedra, and Tabellaria showed both statistically and biologically significant increases following the 400 colonies/L Gloeotrichia treatment, with more modest increases in the 50 colonies/L Gloeotrichia treatment (Appendix S2: Fig. S1). Mougeotia, Oedogonium, and Zygnema increased following the addition of 400 Gloeotrichia colonies/L, suggesting a shift toward filamentous green algae dominance relative to the non-Gloeotrichia controls (Appendix S2: Fig. S1); however, only the responses in Mougeotia reached statistical significance (Appendix S2: Table S1). Similarly, the filamentous cyanobacteria Anabaena (taxonomically revised to Dolichospermum) and Oscillatoria increased immediately following Gloeotrichia addition, but declined to baseline levels thereafter (Appendix S2: Fig. S1). Although there were some changes in the biovolume of other taxa, there were no statistically significant treatment effects on any of the individual cryptophyte, chrysophyte, or dinoflagellate taxa (Appendix S2: Table S1).

Path analysis
Individual phytoplankton taxa responded to both Gloeotrichia's direct non-nutrient effects and its indirect nutrient-mediated effects in our path analysis models (Table 1, Fig. 1). Altogether, more taxa responded significantly to Gloeotrichia's indirect, nutrient effects than to direct, nonnutrient effects. Gloeotrichia's direct, non-nutrient effects were positive, while its indirect, nutrientmediated effects were both positive and negative; the magnitude and direction of the nutrient and non-nutrient effects varied among phytoplankton taxa. The direct, non-nutrient effects were strong and positive for four taxa: one cyanobacterium, one dinoflagellate, and two euglenophytes. Of the 19 taxa with significant responses to Gloeotrichia's indirect, nutrientmediated effects, 10 taxa (primarily diatoms and green algae) had positive responses and nine taxa (primarily cryptophytes, cyanobacteria, euglenophytes, and dinoflagellates) had negative responses.

Aggregate and compositional variability
Despite an increase in the total biovolume of the phytoplankton community and major shifts in individual taxa, aggregate phytoplankton variability did not change with Gloeotrichia treatment, contrary to expectations ( Fig. 2a; one-way ANOVA on CV for each mesocosm, F 3,12 = 0.04, P = 0.99). Moreover, compositional variability decreased with increasing Gloeotrichia density (Fig. 2b).
While the results were robust to the choice of distance metric (Jaccard vs. Bray-Curtis), there were differences in the strength of pairwise comparisons of temporal variability among treatments. Based on a Jaccard presence/absence distance metric, dispersion of the non-Gloeotrichia control was significantly different from all +Gloeotrichia treatments: 25 colonies/L treatment (95% CI, 0.0087-0.094, P = 0.017); 50 colonies/L treatment (CI, 0.017-0.10, P = 0.006), and 400 colonies/L Gloeotrichia treatment (CI, 0.065-0.15, P = 0.00003). There were also significant differences between the 25 colonies/L treatment and the 400 colonies/L treatment (CI, 0.014-0.099, P = 0.009) and between the 50 colonies/L and the 400 colonies/L treatment (CI, 0.005-0.091, P = 0.026). There was no significant difference between the 25 and the 50 colonies/L treatments.
Differences in temporal variability across treatments calculated using the Bray-Curtis metric showed the same pattern in terms of strength of differences between treatments, but the differences were generally smaller and only three comparisons-control vs. 25 colonies/L, control vs. 50 colonies/L, and control vs. 400 colonies/Lwere significantly different in temporal variability. The difference in dispersion results between the two metrics indicates that changes in community variability were driven primarily by changes in the presence and absence of species, rather than changes in relative abundance.

Network analysis
The phytoplankton communities in the highest Gloeotrichia density treatment exhibited significantly higher network complexity than communities in the non-Gloeotrichia control (Fig. 3), as indicated by more nodes (i.e., greater richness) and edges, higher network density and transitivity, and lower mean shortest path lengths in the simulated Erd} os-R enyi networks ( Fig. 4; permutational ANOVA, all F 3, 3996 ≥ 911.85, P < 0.001). We observed significantly more taxa co-occurrences in the 400 colonies/L treatment vs. control (Fig. 4), with the identity of the highest degree nodes changing with the density of the Gloeotrichia treatment (Fig. 3). In general, the structural Notes: The estimated regression coefficients with standard errors (SE) and P-values are given for the direct, non-nutrient and indirect, nutrient effects for each phytoplankton taxon, sorted by phytoplankton division. Effects that are significant after Bonferroni correction are highlighted in boldface. characteristics of the observed networks closely followed the patterns of the simulated networks (Fig. 4). Tukey's post hoc tests showed significant differences in all network structural characteristics among treatments.
Interestingly, node degree in the phytoplankton networks may be potentially related to a taxon's sensitivity to Gloeotrichia's indirect, nutrientmediated effects and direct, non-nutrient effects. Two of the phytoplankton taxa with the highest node degree ("hubs") in the 400 colonies/L treatment (diatoms Asterionella and Synedra) exhibited strong positive responses to Gloeotrichia's indirect, nutrient-mediated effects. In contrast, two of the taxa with the lowest node degree in the 400 colonies/L treatment (the dinoflagellate Gymnodinium and euglenophyte Trachelomonas) exhibited some of the largest positive responses to Gloeotrichia's direct, non-nutrient effects.

DISCUSSION
Studies with macroscopic organisms have demonstrated that increases in a strongly interacting taxon can have widespread effects on the structure and variability of ecological communities (Paine 1980, Power et al. 1985. Here, we show that there are also microscopic strong interactors: This study demonstrates that the cyanobacterium Gloeotrichia echinulata can be a strong interactor with the ability to increase biomass, stabilize community composition, and increase network complexity of the phytoplankton community. This experiment simulated the natural blooms of Gloeotrichia that occur in lownutrient lakes, which translocate large quantities of N and P from the sediments into the nutrientlimited pelagic zone; these nutrients can then become available through multiple mechanisms to increase other phytoplankton. Interestingly, some of Gloeotrichia's effects were counter to expectations, in regard to both our own predictions and conventional ecological wisdom on the effects of cyanobacterial blooms on phytoplankton communities. First, our data contradict the expectation that cyanobacteria generally have negative effects on eukaryotic algae, and add to the small but growing body of literature that demonstrates stimulatory effects by cyanobacteria, especially in natural in situ phytoplankton communities (e.g., Pitois et al. 1997, Nõges et al. 2004, Suikkanen et al. 2005. Second, contrary to our initial predictions, Gloeotrichia decreased the compositional variability of the phytoplankton community, despite increasing the biovolume, richness, evenness, and diversity (as Temporal variability of the phytoplankton community after Gloeotrichia addition. Aggregate variability (a) was calculated as the coefficient of variation (CV) of total phytoplankton biovolume (BV). Compositional variability (b) was calculated using multivariate dispersion (Anderson 2006), here visualized as a principal coordinates (PCoA) ordination of phytoplankton community composition. Each point in b represents the mean value of a treatment on a given date. In both panels, error bars are standard errors across the four replicates. Fig. 3. Observed networks of phytoplankton co-occurrences in the four Gloeotrichia density treatments (0, 25, 50, and 400 colonies/L). Each colored circle represents a node, or taxon in the phytoplankton community, with the line (or edge) connecting taxa representing an observed co-occurrence. The color of the node denotes the division of the phytoplankton taxon, and the size of the node is scaled to the node degree, the number of edges connected to each taxon (hubs are the largest nodes in each network). previously reported by Carey et al. 2014a). Consequently, while we generally think of cyanobacteria as increasing phytoplankton variability, Gloeotrichia's stimulation of other phytoplankton resulted in a more predictable community composition in this experiment.
The positive effects of Gloeotrichia on phytoplankton community stability correspond with increased phytoplankton network complexity at higher Gloeotrichia densities. This increase in network complexity is driven by increases in taxa richness (as indicated by the number of nodes) as well as the number of co-occurrences among taxa: At high Gloeotrichia densities, the community had more taxa that consistently co-occurred. This follows the Bray-Curtis and Jaccard dispersion results, which suggest that Gloeotrichia's primary mechanism of increasing network complexity in the 400 colonies/L treatment was by stimulating many rare and moderately abundant taxa, not by changing the abundance of the dominant species.
The close correspondence in structural metrics between the empirical and simulated networks emphasizes that the phytoplankton communities in the mesocosms within a treatment were much . Phytoplankton community network complexity increased with Gloeotrichia density, as indicated by a greater number of nodes and edges, decreased mean shortest path length, higher network density, and higher transitivity. The boxplots show the structural network characteristics of the simulated 1000 networks; the red triangle denotes the observed network attributes. more similar to each other than among treatments over the course of the experiment, even though the simulated networks were generated from resampling observations both during and after the Gloeotrichia blooms. In addition, the similarity between the empirical and simulated networks demonstrates that the simulated networks, derived from 1000 Monte Carlo simulations, successfully captured the underlying network properties of the observed phytoplankton communities within each Gloeotrichia treatment.
We observed that the mean shortest path length between taxa was much lower, and transitivity and network density were much higher, in the 400 colonies/L treatment than in the non-Gloeotrichia control. In particular, transitivity was 0.62 in the 400 colonies/L observed network (Fig. 4), which is extremely high for an ecological community (e.g., Kara et al. 2012, Patrick et al. 2014) and indicates a very dense network. Importantly, because transitivity is robust to the number of nodes, our data suggest that Gloeotrichia increases phytoplankton network complexity not only by increasing richness but also by increasing co-occurrences, as illustrated by the very high node degree of some taxa (Fig. 3). In the 400 colonies/L treatment, the network hubs exhibited large positive responses to Gloeotrichia's nutrient stimulation, without decreases in evenness (Carey et al. 2014a).
Across treatments, Gloeotrichia at high densities may increase the "small-world" nature of the phytoplankton community relative to non-Gloeotrichia controls. Small-world networks, as defined by Watts and Strogatz (1998), are characterized by a low shortest mean path length and high transitivity, as well as by having hubs. Small-world networks tend to be robust to perturbations and the random loss of nodes, but are vulnerable to targeted removal of hub taxa (Sole andMontoya 2001, Olesen et al. 2007), which has implications for the stability of a community network. Consequently, it may be possible that cyanobacterial blooms increase the resilience of phytoplankton communities to future disturbance after a bloom ends. While network analysis is in its early stages of application for phytoplankton communities, it holds much promise for analyzing relationships of co-occurring taxa. In particular, integrating measurements of community stability (via dispersion) and network complexity may provide insight into how communities respond to changes in the abundance of strongly interacting taxa.
Gloeotrichia's indirect nutrient stimulation effects acted in concert with its direct, non-nutrient effects to constrain the phytoplankton community, decreasing compositional variability but not aggregate variability. Micheli et al. (1999) described some communities with low compositional variability as exhibiting synchrony, which is often driven by parallel responses of different taxa to environmental drivers such as nutrients, as we observed in the high Gloeotrichia treatments relative to the control (Appendix S2: Fig. S1).
The decrease in compositional variability is likely also due to additional direct mechanisms that enhance Gloeotrichia's positive nutrient effects, such as the excretion of bioactive secondary metabolites that other phytoplankton can use for their own metabolism, or antibacterial or antifungal compounds that promote phytoplankton growth (Suikkanen et al. 2005, Leflaive andTen-Hage 2007). While we did not measure these compounds, other studies in both laboratory and field settings have quantified the positive effects that cyanobacteria can exert on other phytoplankton via biochemical interactions (Suikkanen et al. 2005, Leflaive andTen-Hage 2007). We note that there are likely other mechanisms by which Gloeotrichia may modify the environment (e.g., by altering light availability) for which we do not have data but may have had important effects (Carey et al. 2014b). Regardless of which mechanisms were occurring in the mesocosms, our results support the hypothesis that both facilitation and the increase in weak pairwise interactions can increase community stability and resilience to disturbance, as observed in other communities (McCann et al. 1998, Downing et al. 2014. Based on the path analysis models, the magnitude of Gloeotrichia's indirect, nutrient and direct, non-nutrient effects varied among phytoplankton taxa. Diatoms, green algae, and chrysophytes were generally more sensitive to nutrientmediated effects, whereas the cyanobacteria, dinoflagellate, and euglenophyte taxa responded to both the non-nutrient-and nutrient-mediated effects of Gloeotrichia density (Table 1). This result may reflect the ecological and life history characteristics of these phytoplankton groups. For example, many diatom and green algae taxa are able to respond quickly to increased nutrient availability (reviewed by Reynolds 2006), which may explain their rapid increase in biovolume in the mesocosms with high Gloeotrichia densities. We think it unlikely that the increases were due to phytoplankton attached to Gloeotrichia colonies because the colonies were cleaned very thoroughly, and the growth rates (r) of these taxa within the mesocosms are consistent with published phytoplankton growth rates (Reynolds 2006, Carey et al. 2014a). Carey and Rengefors (2010) observed that Gloeotrichia stimulated a diverse array of cryptophytes, cyanobacteria, diatoms, and dinoflagellates, even under high-nutrient conditions. There, too, the majority of taxa exhibited positive responses to Gloeotrichia presence.
The sensitivity of individual phytoplankton taxa to Gloeotrichia's nutrient and non-nutrient effects has implications for the larger phytoplankton community network. The diatoms Asterionella and Synedra, network hubs in the experiment, exhibited very strong positive responses to Gloeotrichia's nutrient-mediated effects. Studies examining networks composed of taxa across multiple trophic levels find that hubs play a critical role in maintaining network connectedness and providing resilience to disturbance (e.g., Sole andMontoya 2001, Olesen et al. 2007); we expect hubs play similar roles in networks composed of taxa on the same trophic level. Thus, Gloeotrichia's indirect nutrient stimulation effects may be more important than its direct, non-nutrient effects for altering overall network complexity and community structure.
In conclusion, as a strongly interacting species, Gloeotrichia and its nutrient subsidies determined the structure (i.e., co-occurrence networks) and compositional variability of the phytoplankton community through both direct and indirect effects. In contrast to conventional ecological wisdom, blooms of this cyanobacterium appear to act as a stabilizing force in the phytoplankton community of a low-nutrient lake, as demonstrated by increased temporal stability in community composition and increased network complexity. Importantly, these stabilizing effects persisted even after Gloeotrichia decreased. Are these effects a general pattern among strongly interacting species, or are they a consequence of Gloeotrichia's ability to simultaneously affect many members of the phytoplankton community through both direct and indirect mechanisms? Our results suggest that strongly interacting species such as Gloeotrichia can have unanticipated and profound effects, substantially altering the temporal dynamics and stability of complex communities.