Increased anoxia following species invasion of a eutrophic lake

Species invasions can disrupt aquatic ecosystems by re‐wiring food webs. A trophic cascade triggered by the invasion of the predatory zooplankter spiny water flea (Bythotrephes cederströmii) resulted in increased phytoplankton due to decreased zooplankton grazing. Here, we show that increased phytoplankton biomass led to an increase in lake anoxia. The temporal and spatial extent of anoxia experienced a step change increase coincident with the invasion, and anoxic factor increased by 11 d. Post‐invasion, anoxia established more quickly following spring stratification, driven by an increase in phytoplankton biomass. A shift in spring phytoplankton phenology encompassed both abundance and community composition. Diatoms (Bacillaryophyta) drove the increase in spring phytoplankton biomass, but not all phytoplankton community members increased, shifting the community composition. We infer that increased phytoplankton biomass increased labile organic matter and drove hypolimnetic oxygen consumption. These results demonstrate how a species invasion can shift lake phenology and biogeochemistry.

Anthropogenic impacts and manipulations introduce non-native species to aquatic ecosystems, which can disrupt as well as reshape energy flow and mass transfer in food webs (Boivin et al. 2016;Crespo et al. 2018;Bartley et al. 2019), and can threaten ecosystem stability and services (Vander Zanden et al. 1999;Dudgeon et al. 2006;Lopez et al. 2022).Many studies of species invasions in lakes focus on food web changes, overlooking indirect feedbacks to lake biogeochemistry.The invasion literature on "zoogeochemistry" is mostly focused on nutrient shunting and relocation.A notable example includes the role of dreissenid mussels in shunting carbon, nitrogen, and phosphorus from pelagic to benthic habitats (Ozersky et al. 2015;Li et al. 2021;Vanni 2021).However, there are few examples of food web disruptions that lead to alterations in oxygen dynamics in lakes.The paucity of limnological datasets that involve a species invasion and include both lake biology and biogeochemistry has limited our understanding of how species invasions affect lake biogeochemical processes.
Lake Mendota is a eutrophic lake in Wisconsin, USA, with a long history of limnological observations through the North Temperate Lakes Long Term Ecological Research program (NTL-LTER).Lake Mendota experienced a population irruption of the non-native predatory zooplankton spiny water flea (Bythotrephes cederströmii) in 2009 (Walsh et al. 2016b; note the species name was previously Bythotrephes longimanus).Spiny water flea invasions of temperate North American lakes are widespread, and result in reduced zooplankton diversity (Yan et al. 2011) and disrupted trophic structure (Rennie et al. 2011).In Lake Mendota, Spiny water flea predation reduced the abundance of the lake's dominant zooplankton, Daphnia pulicaria, which is a keystone species in the food web, and a key food item for native fish populations (Johnson and Kitchell 1996;Walsh et al. 2016a;Rani et al. 2022).The reduction in Daphnia abundance caused a water clarity decline due to reduction in Daphnia grazing pressure on phytoplankton (Walsh et al. 2016a).This shortened the duration and intensity of Lake Mendota's spring clearwater phase (Matsuzaki et al. 2020) due to higher diatom biomass (Walsh et al. 2018), and an earlier appearance of Cyanophyta (Cyanobacteria) during clearwater phase (Rohwer et al. 2023a).
Notable shifts in oxygen dynamics in Lake Mendota have also been observed (Ladwig et al. 2022).Ladwig et al. (2021a) applied a mechanistic aquatic ecosystem model that was able to replicate hypolimnetic dissolved oxygen (DO) consumption and bottom-water anoxia from 1995 to 2015, but model performance declined post-2009, with the model overestimating hypolimnetic DO.The model did not include zooplankton grazing on phytoplankton, leading us to wonder about the potential role of food web interactions on DO dynamics.
Past studies have quantified the impacts of trophic cascades on lake ecosystems (Carpenter and Kitchell 1993;Carpenter et al. 2001), including on Lake Mendota (Walsh et al. 2016a).Adding to this knowledge, we explore the aftermath of a trophic cascade by quantifying how the impacts of the spiny water flea irruption in Lake Mendota resulted in an increase in the annual spatial and temporal extent of anoxia using 24 years of long-term data.We hypothesize that the spiny water flea irruption caused an abrupt phenological shift in lake anoxia stemming from increased algal biomass.Mechanistically, increased grazing pressure on planktivorous zooplankton by spiny water flea caused an increase in spring phytoplankton biomass, which resulted in enhanced hypolimentic DO consumption through mineralization of phytoplankton biomass (Fig. 1).

Lake Mendota
Lake Mendota is a 3961 ha, dimictic, eutrophic lake with a maximum depth of 25 m (NTL-LTER 2021).Physical, chemical, and biological data have been collected fortnightly (when ice-free) to monthly (when ice-covered) by the NTL-LTER since 1995 (Magnuson et al. 2006).
All in-lake measurements were collected at the central deep hole during daylight hours (43.0988N, À89.4054 W, 25 m depth) and include: ice duration (Magnuson et al. 2021a), integrated water-column measurements from 0 to 20 m of zooplankton and spiny water flea density measured by zooplankton net tows (Magnuson et al. 2022c), integrated water-column measurements from 0 to 8 m of phytoplankton density and biomass (Magnuson et al. 2022a), and depthdiscrete measurements of DO, water temperature, nitrate/ nitrite (NO 3 À /NO 2 À ), soluble reactive phosphorus (SRP), dissolved reactive silica, and Secchi depth (Magnuson et al. 2021b(Magnuson et al. , 2022b;;Rohwer and McMahon 2022;Magnuson et al. 2023).Discharge data from the Yahara River were obtained from USGS gage 05427718 (U.S. Geological Survey 2022).Precipitation data were obtained from the Dane County Regional Airport (NOAA GHCND:USW00014837) (Menne et al. 2012).
Data analysis to explore hypolimnetic anoxia Annual anoxia was quantified using anoxic factor, AF, which was calculated according to Nürnberg (1995) as: where the daily highest areas, A i , with DO < 1.5 g m À3 (Chapra and Canale 1991) were summed up over time t i (days); A s is surface area (m 2 ).
To identify break points in the time series of anoxic factor, we applied an Ordinary Least Square Cumulative Sum test to quantify the timing of a significant structural change, followed by the "breakpoints" function from strucchange R package (Zeileis et al. 2002).Years were grouped as either pre-irruption (n = 14) or post-irruption(n = 9), with Jan 2010 as the breakpoint, and groupings compared with Wilcoxon rank sum tests.
Biweekly water temperature measurements were temporally interpolated to daily values using linear, constant, and spline interpolation.The transition from mixed to stratified conditions was defined as a density gradient between surface and bottom water layers > 0.1 g m À3 and the water column average temperature > 4 C, and stratification duration as the number of days between stratification onset and offset.Water column stability was quantified using Schmidt stability (Idso 1973).DO measurements were temporally interpolated using spline interpolation.Nutrient data were temporally linearly interpolated to weekly values with the NTLlakeloads R package (Dugan 2023), area-averaged, and labeled as either surface or deep water based on a mean thermocline depth of 13 m.
Clearwater phase intensity (CWP in meter-days) was quantified by integrating Secchi depths between April and June: where z secchi are Secchi depths (m) linearly interpolated to daily values.This method allowed us to quantify year-to-year variability in the intensity of clearwater phase without the need to define a threshold that would arbitrarily correspond to the formation or breakdown of a clearwater phase.Hypolimnetic oxygen consumption fluxes, volumetric, and areal consumption fluxes were calculated from interpolated daily DO data according to the Livingstone and Imboden (1996) model.Only DO concentrations below the bottom metalimnion depth, calculated using the rLakeAnalyzer R package (Read et al. 2011), were used to calculate vertical anoxia height.
Candidate predictor selection to explain interannual variability in anoxia in a multiple linear model was determined with the "Boruta" random forest classifier function from the Boruta R package (Kursa and Rudnicki 2010).Candidate predictor importance was calculated using the relaimpo R package sensu Lindeman et al. (1980).For predictor analysis we included: annual stratification duration, stratification start and break-down date, water column stability, ice cover dates, and duration from the previous winter, summer volumetric, areal, and total oxygen sink, annual days of phytoplankton biomass surpassing 1.0 mg L À1 , annual total Yahara River discharge, annual total precipitation, annual spring clearwater intensity, maximum spring Secchi depth, annual average spiny water flea biomass, annual average diatom biomass, annual average Cyanobacteria biomass, annual average Daphnia biomass of D. mendotae and D. pulicaria, summer stratified average surface and bottom SRP and nitrate concentrations, and summer stratified average silica concentrations.Important candidate predictors were analyzed using a linear regression model.

Phytoplankton and anoxia phenology
Sampling dates were divided annually into four "lake seasons" based on water temperature profiles: 1) ice, 2) spring mixed, 3) stratified, and 4) fall mixed.Phytoplankton biomass within each season and year were averaged to account for uneven sampling and compared between seasons, pre-and post-the 2010 breakpoint.Changes in phytoplankton community composition were further investigated using the vegan R package (Oksanen et al. 2020).Bray-Curtis distances between average annual communities were calculated and analyzed with nonmetric multidimensional scaling.Analysis of similarities (ANOSIM) was applied to the distance matrix to determine if years within an invasion group were statistically more similar to themselves than to all years.Shannon and Simpson diversity were calculated for each year and averaged by invasion group to compare preinvasion and postinvasion diversity.
Oxygen phenology was investigated as the difference in days between stratification development and when the lowest hypolimnion layer dropped to < 1.5 g m À3 DO.

Results
Anoxia increased with spiny water flea irruption Anoxic factor increased from an average of 56 (AE 6.8) days pre-Jan 2010 to 67 (AE 4.4) days post-Jan 2010, concordant with the spiny water flea irruption.Additionally, we observed an increase in total oxygen consumption flux of 0.02 AE 0.02 g m À3 d À1 and a decrease in average spring clearwater intensity by 169 AE 99 m-days per year (Fig. 2a-d).Daphnia mendotae biomass increased by 7.7 AE 9.9 mg L À1 (Fig. 2e, p < 0.05), whereas Daphnia pulicaria biomass did not significantly change (Fig. 2f).Although Walsh et al. (2017) attribute post-spiny water flea water clarity declines to a loss of Daphnea pulicaria, this dynamic is masked in Fig. 2 F by high variability in zooplankton biomass.D. pulicaria proportion declined precipitously in spring (Supporting Information Fig. S1), but we did not include compositional data as linear model input.Phytoplankton-related metrics broadly showed significant change (p < 0.05) between pre-and post-spiny water flea regimes (Fig. 2g,h): average total days with phytoplankton biomass >1.0 mg L À1 increased by 76 AE 20 d per year, and average diatom biomass increased by 0.6 AE 0.6 mg L À1 .Average Cyanobacteria biomass did not significantly change during the summer stratified season (Fig. 2i).Stratification, stability, ice duration, and NO 3 À /NO 2 À did not change significantly (Fig. 2j-l,m).SRP decreased with a significant decline in the surface layer (p < 0.01), and accordingly N : P ratio increased significantly (Fig. 2m,o).Eight predictors were significant predictors for annual anoxia: annual stratification duration and breakdown, annual ice-cover duration, annual duration of phytoplankton biomass over 1.0 mg L À1 , annual average spiny water flea biomass, annual clearwater phase intensity, and annual average surface water layer SRP and bottom water layer SRP concentrations (adjusted R 2 of 0.73, p-value <0.05).Stratification duration and breakdown, phytoplankton biomass and spiny water flea biomass had a positive correlation with anoxia, whereas ice cover duration and both SRP concentrations were inversely correlated.Phytoplankton biomass (30%), stratification breakdown (19%), stratification duration (14%), and clearwater phase intensity (13%) drove interannual variability in anoxic factor, whereas the remaining five predictors accounted for 24%.
An alternate interpretation of Fig. 2a,g is a gradual phytoplankton biomass and anoxic factor increase since 2002.However, we posit that 2002/2003 are outlier years with low precipitation and discharge that resulted in unusual mesotrophic conditions (Lathrop and Carpenter 2014; Supporting Information Fig. S2), rather than the start of a long-term trend.Additionally, in late 2015 zebra mussels invaded the lake, potentially confounding our results.We reanalyzed these same dynamics excluding years after 2015 (Supporting Information Fig. S3).Anoxia, spiny water flea, diatoms and N:P ratios still significantly increased, whereas days of phytoplankton biomass > 1.0 mg L À1 did not.Stratification breakdown (45%), phytoplankton biomass (24%), spiny water flea (15%), diatoms (10%), and bottom SRP (4%) were significant predictors.

Spring phytoplankton changes are coincide with increase in anoxia
After establishing a step-change increase in anoxia after the spiny water flea irruption (Fig. 2a), we observed a concurrent phenological shift in anoxia.Anoxia onset shifted earlier, from July to June (Fig. 3b).However, shifts in anoxia gradient, maximum vertical extent, and breakdown timing were not apparent.Given that the regression model identified phytoplankton biomass as driving anoxia variability, we compared phytoplankton biomass before and after the irruption (Fig. 3a).Post-irruption phytoplankton biomass was elevated prior to the period of anoxia.
Mean phytoplankton biomass during the mixed spring season increased to concentrations typical of the stratified summer season, from 1 AE 1 preinvasion to 3 AE 2 mg L À1 postinvasion (p < 0.005) (Fig. 4a).Similarly, biomass under lakeice increased to concentrations previously typical of spring, from 0.3 AE 0.3 to 2 AE 2 mg L À1 (p < 0.005).In contrast, later in the season no statistically significant change in total biomass was observed during the stratified summer season (p > 0.1) and more modest increases were observed during the fall mixed season (p < 0.05).
Post irruption, the lag between stratification and anoxia onset, when DO in the lowest hypolimnion layer was <1.5 g m À3 , decreased by nearly 2 weeks, from 51 AE 9 d to 39 AE 15 d (p < 0.05) (Fig. 4b), confirming the shift toward earlier anoxia onset.

Diatoms drive phytoplankton increase
Seasonal phytoplankton community composition was broadly consistent before and after the spiny water flea invasion; spring was dominated by Bacillariophyta (diatoms), summer was dominated by Cyanophyta (Cyanobacteria), and fall was dominated by a mix of diatoms and Cyanobacteria (Fig. 5).This follows a typical phenology for a eutrophic lake sensu PEG model (Sommer et al. 1986), which was previously documented in Lake Mendota (Carey et al. 2016).Diatoms were predominantly responsible for the increase in spring biomass, comprising the majority of the phytoplankton community in all years (67 AE 20% and 65 AE 25%, respectively).Diatom biomass in the spring increased two-fold, from 0.9 AE 0.9 to 2 AE 2 mg L À1 (p = 0.08), but diatom proportion to total phytoplankton biomass remained relatively constant.
Although diatoms dominated the spring phytoplankton community, three other phytoplankton divisions also contributed to the increase in spring phytoplankton biomass.
Chlorophyta (green algae) remained at 5-9% of the community but increased four-fold, from 0.04 AE 0.02 to 0.1 AE 0.1 mg L À1 (p < 0.005), Cyanobacteria remained at 5-9% of the community, but increased by six-fold, from 0.03 AE 0.03 to 0.2 AE 0.2 (p < 0.005), and Pyrrhophyta (dinoflagellates) remained at 1-3% of the community but increased three-fold, from 0.02 AE 0.03 to 0.05 AE 0.03 (p < 0.05).Two phytoplankton divisions decreased their relative contribution, thus shifting the community composition.Cryptophyta (cryptophytes) decreased from 17 AE 12 to 9 AE 7% of the community (p = 0.07), and Chrysophyta (golden algae) decreased from 3 AE 2 to 1 AE 1% of the community ( p = 0.05), although the absolute biomass of both taxa remained constant.Excluding  zebra mussel years from the analysis did not change the observed trends, although Cyanobacteria increased by only three-fold (p < 0.05).A comparison of community composition found the communities were more similar during years with the same invasion status than among all years (ANOSIM significance < 0.05), but these changes were modest enough that phytoplankton Shannon and Simpson diversity did not significantly change.

Discussion
Pre-spiny water flea interannual anoxia dynamics in Lake Mendota were predominantly driven by changes in water column stability and stratification duration (Ladwig et al. 2021a), in agreement with global anoxia observations (Jane et al. 2022).We propose a mechanism where the step change increase in anoxia observed after spiny water flea was driven by indirect ecosystem impacts of spiny water flea on Lake Mendota's phytoplankton.Our proposed step-change mechanism is driven by spiny water flea predation on the zooplankton grazer Daphnia pulicaria that enabled spring diatoms to flourish (Walsh et al. 2017(Walsh et al. , 2018)).The link between epilimnetic phytoplankton biomass and elevated hypolimnetic oxygen consumption is well established for eutrophic lakes (Paerl 1988).The increased springtime phytoplankton biomass observed in this study likely increased the settling flux of organic matter and availability of labile substrates for hypolimnetic mineralization.Given that physical factors like stratification did not change following the species invasion, the observed increases in springtime phytoplankton biomass and anoxia, as indicated by the decrease in lag between stratification development and anoxia onset, seems beyond coincidence.
Alongside anoxia changes, epilimnetic SRP concentrations decreased by 65% following the spiny water flea invasion (Walsh et al. 2019).While this pattern runs counter to positive relationships between nutrient availability and phytoplankton biomass (Conley et al. 2009), biophysical processes may provide an explanation.Whiting events occur when phytoplankton blooms raise eplimnetic pH through the uptake of inorganic carbon, thus triggering the precipitation of calcium carbonate and the co-precipitation of SRP (Walsh et al. 2019).Simultaneously, increased phytoplankton biomass may have also reduced surface layer SRP concentrations due to uptake.Nutrients are often considered a bottom-up control on phytoplankton growth; but in this case, we speculate that the top-down processes responsible for increased phytoplankton growth in turn resulted in reduced SRP.This highlights another possible biogeochemical impact of the spiny water flea irruption.
Cascading impacts of species invasions can also extend through time, as disturbed ecosystems may be more vulnerable to future disturbance (Turner et al. 2020;Spear et al. 2021).The susceptibility of Lake Mendota to spiny water flea may stem from a biomanipulation in the 1980s, when piscivorous fish were stocked to improve water clarity.Fewer planktivorous fish opened a trophic niche that spiny water flea filled (Walsh et al. 2017).Like the post-invasion period, smaller-bodied D. galeata mendotae, which reach lower biomass and result in lower water clarity, were the prominent Daphnia species before the 1987 biomanipulation (Lathrop et al. 2002).Consistent DO data prior to 1995 is too scarce to track anoxic factor from the 1970s through the biomanipulation to present.However, historic oxygen profiles taken in 1906/07 and 1961-1971 show that Lake Mendota pre-biomanipulation suffered from similar hypolimnetic oxygen consumption rates compared to present day (Stewart 1976).
In Lake Mendota, zebra mussels invaded in 2015 (Spear et al. 2022), potentially confounding the second half of our post-2009 analysis.In Lake Mille Lacs, USA, a spiny water flea invasion had no net effect on phytoplankton biomass because of a simultaneous zebra mussel invasion (Rantala et al. 2022).In Lake Simcoe, Canada, a simultaneous invasion of dreissenids and spiny water flea improved water quality by lowering hypolimnetic hypoxia and increasing water clarity (Goto et al. 2020).However, in Lake Mendota no water clarity change was observed with the subsequent zebra mussel invasion (Spear et al. 2021), and we did not observe an associated change in anoxia or phytoplankton biomass.Dreissenids are known to shift phytoplankton community composition (Naddafi et al. 2007), and Rohwer et al. (2023a) did observe an earlier seasonal timing of Cyanobacteria onset in the microbial community post-zebra mussels; however, note that our lake season "spring mixed" differs from the "spring" season in Rohwer et al. (2023a) in that "spring mixed" also includes a large portion of clearwater phase.Nonetheless, our finding that Cyanobacteria and green algae increased along with diatoms, by maintaining but not increasing their proportion of total phytoplankton biomass, holds true with or without the zebra mussel years included.Further, a reanalysis of the time series including only years prior to the zebra mussel invasion resulted in similar regression model results.Given our proposed mechanism where higher phytoplankton biomass increases anoxia, it makes intuitive sense that since zebra mussels did not impact water clarity in Lake Mendota, they also did not impact anoxia.We present one plausible mechanistic pathway of decreased phytoplankton grazing leading to increased anoxia, but note that trophic cascades from zooplankton grazers to zooplankton and phytoplankton biomass are complex beyond grazing pressure, food quality and nutrient recycling (McQueen and Post 1988;Vanni and Findlay 1990).
A longer anoxia duration extends the period of reduced fish habitat due to oxythermal stress in Lake Mendota (Magee et al. 2019), as well as the season of cyanotoxin production (Rohwer et al. 2023a).We attribute an anoxic factor increase of 11 d (corresponding to an increase of hypolimnetic anoxic volume of 19 Million m 3 ) to a spiny water flea irruption in Lake Mendota.Climate change impacts, such as shortened ice duration (Sharma et al. 2021), increased water temperature (Woolway et al. 2022), and decreased wind speeds (Magee et al. 2016), as well as freshwater salinization (Ladwig et al. 2021b) are likely to further increase summer water column stability, amplifying habitat constraints.Increased anoxia represents a biogeochemical invasion impact that has not been previously accounted for, and identifies invasive species as interacting drivers of anoxia alongside climate change and road salt.Comprehensive long-term monitoring programs that collect observations of food webs, physical characteristics, and biogeochemistry are essential to studying impacts of these interacting drivers on lake ecosystems.

Fig. 1 .
Fig.1.Consequences of the trophic cascade triggered by spiny water flea in Lake Mendota.(a) Pre-invasion: Diatom blooms after ice melt are grazed on by zooplankton (esp.Daphnia), resulting in a spring clearwater phase that is characterized by deeper Secchi depths (blue shading).After the lake stratifies, hypolimnetic anoxia develops (pink shading) and Cyanobacteria become the dominant phytoplankton.(b) Post-invasion: Spiny water flea graze on Daphnia, in turn reducing grazing pressure on diatoms.The spring diatom bloom extends and intensifies, and the duration and magnitude of the spring clearwater phase decreases.The additional deposition of organic matter from sinking phytoplankton biomass leads to increased hypolimnetic consumption of oxygen.This reduces the lag-time between stratification onset and the formation of hypolimnetic anoxia.

Fig. 2 .
Fig. 2. Long-term dynamics of lake variables.(a) Anoxic factor over time.Breakpoint analysis of anoxic factor identified 2010 as a breakpoint.The dotted vertical line indicates the breakpoint.(b) Modeled hypolimnetic total oxygen depletion flux calculated from observed DO data.(c) Spring clearwater phase intensity (CWP) over time quantified from Secchi depth.(d) Spiny water flea biomass over time.(e) Daphnia mendotae biomass over time.(f) Daphnia pulicaria biomass over time.(g) Days per year with a phytoplankton biomass over 1.0 mg L À1 .(h) Diatom biomass over time.(i) Cyanobacteria biomass over time.(j) Stratification duration over time.The gray ribbon represents the potential uncertainty between sampling points.The red line represents the spline interpolation.(k) Schmidt stability, as a proxy for water column stability, over time.(l) Ice duration over time.(m) Summer stratified average SRP concentrations in the surface water layer (solid line) and bottom water layer (dotted line).Box plot highlights only the surface layer SRP concentrations.(n) Summer stratified average NO 3 À /NO 2 À concentrations in the surface water layer (solid line) and bottom water layer (dotted line).Box plot highlights only the surface layer NO 3 À /NO 2 À concentrations.(o) Annual average N : P ratios in the surface water layer (solid line) and bottom water layer (dotted line).Box plot highlights only the surface layer N : P ratio.

Fig. 3 .
Fig. 3. Phytoplankton biomass and anoxia depth.(a) Annual time series of total phytoplankton biomass before (pre) and after (post) the spiny water flea irruption in late 2009.Gray lines denote the average timing of ice-off and spring stratification.(b) Annual time series of anoxia transition depth (DO < 1.5 g m À3 ).

Fig. 5 .
Fig. 5. Spring phytoplankton biomass composition.(a) Barplots of average annual spring phytoplankton taxa biomass in the spring mixed season.(b)Barplots of average annual spring phytoplankton taxa relative abundances.The other category includes Xanthophyta, Euglenophyta, Haptophyta, and unclassified organisms.The spring mixed season is majority diatoms (Bacillariophyta), which increase along with green algae (Chlorophyta), Cyanobacteria (Cyanophyta), and dinoflagellates (Pyrrhophyta).Despite shifts in spring phytoplankton community composition to include more green algae, Cyanobacteria, and dinoflagellates, diatoms remained most abundant in spring and drove the increase in spring phytoplankton biomass.

Fig. 4 .
Fig. 4. Seasonal phytoplankton biomass and anoxia onset timing before and after spiny water flea.(a) Boxplots of annual averages of phytoplankton biomass in each season.The ice and spring mixed seasons had increased phytoplankton biomass post-spiny water flea ( p < 0.05, Wilcoxon test with Bonferroni correction).(b) Boxplot of lag between stratification development and anoxia onset in days.Anoxia onset occurred sooner after stratification development post-spiny water flea ( p < 0.05).