Evolutionary adaptation to steady or changing environments affects competitive outcomes in marine phytoplankton

The interplay of phytoplankton competition and adaptation affects how phytoplankton, and ultimately marine ecosystems, respond to global warming. However, current ecosystem models that are run under global warming scenarios do not include both processes simultaneously. To fill this gap, we developed an innovative ecosystem model for the Baltic Sea that simulates competition between three phytoplankton functional groups and allows for adaptation to changing temperatures. As adaptation can be affected by the resuspension of dormant resting cells from the sediment, we explicitly implemented this mechanism. We found that resuspension tends to slow down adaptation, and that competition and adaptation influence each other. The outcome of the competition‐adaptation interplay depends on environmental conditions. In a steady environment, competition drives adaptation to individual temperature niches to reduce competition pressure. In a changing environment, adaptation allows inferior competitors to mitigate the dominance of preadapted superior competitors. Our results demonstrate that by neglecting adaptation, models can systematically overestimate warming‐related changes in taxa dominance. Ecosystem models should include both competition and adaptation to accurately simulate phytoplankton responses to global warming. Our model is ideally suited to integrate emerging evolutionary data based on long‐term data series (e.g., from sediment archives) to further improve projections of future ecosystem change.

The functioning of marine ecosystems is affected by the dominant phytoplankton taxa and the functional traits that these taxa express (Litchman et al. 2015).Due to global warming, both taxa dominances and functional traits are changing (Klais et al. 2011;Irwin et al. 2015).Even though competition between phytoplankton taxa and adaptive trait changes affect phytoplankton responses, and ultimately ecosystem responses to global warming (Litchman et al. 2015), current ecosystem models do not consider both processes simultaneously (Laufkötter et al. 2015;Munkes et al. 2021).To study how the interplay of competition and adaptation influences community, population, and trait changes of phytoplankton in response to global warming, we developed an ecosystem model for the Baltic Sea that includes three competing phytoplankton functional groups and allows for adaptation to changing temperatures.
Global warming has already caused ecological responses of phytoplankton, which include shifts in bloom timing and changes in species dominance (Winder and Sommer 2012).As phytoplankton form the base of the marine food web (Fenchel 1988), these responses can alter food web structures and eventually lead to ecosystem-level changes (Edwards and Richardson 2004).In the Baltic Sea, an extension of the phytoplankton growing season (Wasmund et al. 2019), a shift from diatom to dinoflagellate dominance during spring bloom (Klais et al. 2011), as well as an increase in cyanobacterial summer biomass (Suikkanen et al. 2007) have been observed over the past 50 yr.Especially the increasing cyanobacteria represent a threat for the Baltic Sea ecosystem due to their toxicity for higher trophic levels (Chorus and Welker 2021) and their contribution to ocean deoxygenation (Long et al. 2021;Munkes et al. 2021).Cyanobacteria enhance primary production through nitrogen fixation (Hense 2007) and hence promote hypoxia by increasing the microbial oxygen demand in bottom waters as dead biomass is decomposed.In turn, spreading hypoxia enhances the release of iron-bound phosphate from sediments (Conley et al. 2002), promoting cyanobacterial summer blooms and creating a positive feedback mechanism that maintains hypoxia (Vahtera et al. 2007).
Changes in phytoplankton phenology and species dominance do not necessarily result from ecological processes alone but can also be affected by evolutionary adaptation.Owing to their large population sizes and short generation times, phytoplankton possess a high potential to adapt to new environments through mutation and selection.Laboratory experiments demonstrated that phytoplankton can adapt to environmental changes within 200-600 generations, corresponding to a few years in nature (Jin and Agustí 2018).Consequently, phytoplankton already show adaptive changes in temperature-dependent functional traits in response to global warming (Irwin et al. 2015;Hinners et al. 2017).
Both competition and adaptation shape phytoplankton community composition and influence each other.Experiments with bacteria (Lawrence et al. 2012) and phytoplankton (Collins 2011) demonstrated that interspecific competition has a major impact on adaptation to a novel environment.Likewise, a competition-evolution experiment with coccolithophores and diatoms revealed that strong selection pressures affect competitivity and species dominance (Listmann et al. 2020).Thus, to realistically estimate phytoplankton responses to global warming, we need to understand how competition and adaptation influence each other in different environments.
So far, no long-term experiments exist that allow us to study the interplay between phytoplankton competition and adaptation under global warming on realistic time scales.Numerical models can bridge this knowledge gap by simulating species-and ecosystem-level responses to global warming at realistic rates of environmental change.
In the Baltic Sea, ecosystem models are used to perform both hindcasts and future climate projections.Prominent examples include ERGOM (Neumann et al. 2002), SCOBI (Eilola et al. 2009), BALTSEM (Savchuk 2002), BFM (Fransner et al. 2018), and ERSEM (Vichi et al. 2004).While all these biogeochemical models consider physical processes, only ERGOM, SCOBI, and BFM are coupled to full ocean circulation models.The Baltic Sea ecosystem is represented by three to four phytoplankton functional groups, which compete for three to four different nutrients, and are grazed by one to four zooplankton groups.While hindcasts were performed with all models presented here, a multimodel ensemble of ERGOM, SCOBI, and BALTSEM was used for climate projections (Meier et al. 2014).Even if the results of climate projections are relevant for stakeholders, neither the models presented here nor similar models like CEMBS (Dzierzbicka-Głowacka et al. 2013) or ECOSMO (Daewel and Schrum 2013) consider evolutionary adaptation of phytoplankton.Instead of adaptation, however, BFM, ERSEM, and an updated version of ERGOM (Neumann et al. 2022) allow for variable nutrient stoichiometry.Whether variable stoichiometry or the interplay of competition and adaptation is more important for realistic climate projections is beyond the scope of this article.Nonetheless, the absence of adaptation in the above models remains a major uncertainty factor that calls into question the validity of their climate projections.
Contrary to the models presented above, there are ecosystem models that consider adaptation.Phytoplankton cell size is frequently used as adaptive master trait, including sizerelated physiological cell properties and trade-offs (Daines et al. 2014;Sauterey and Ward 2022).Other ecosystem models simulate adaptative changes in the optimum growth temperature (Beckmann et al. 2019), the encystment rate (Hinners et al. 2019), two different functional traits (Le Gland et al. 2021), or overall fitness (Walworth et al. 2020).However, none of the previously mentioned ecosystem models considers competition between different functional groups and adaptation to global warming simultaneously.Current eco-evolutionary models that do so, in turn, lack ecosystem dynamics (De Mazancourt et al. 2008;Northfield and Ives 2013;Barab as and D'Andrea 2016).
Although adaptation can be affected by the resuspension of dormant resting cells from the sediment (Rengefors et al. 2017), current ecosystem models largely ignore this process.So far, resuspension has either been represented in a simplified manner, e.g., by disabling the background mortality of resting cells (Daines et al. 2014), or in a completely conceptual framework (Sundqvist et al. 2018).
Here, we present for the first time an ecosystem model that simultaneously considers competition and adaptation of three major phytoplankton functional groups (dinoflagellates, diatoms, and diazotrophic cyanobacteria) in a global warming scenario, including an explicit representation of resuspension.We apply our model to the Baltic Sea, an ecosystem that is heavily influenced by anthropogenic activities, including above-average levels of warming (Reusch et al. 2018).For each phytoplankton functional group, we simulate the life cycle of one representative species or of a complex that contains multiple species.In addition, our model allows for two flexible temperature-dependent traits: cell size, which responds plastically to changing temperatures, and optimum growth temperature, which can adapt through mutation and selection.We assess future climate-driven community, population, and trait changes, and disentangle to which degree these changes result from competition or adaptation.Taking advantage of our detailed representation of phytoplankton life cycles, we explicitly simulate the resuspension of dormant resting cells from the sediment and study potential effects of resuspension on adaptive changes in the optimum temperature.

Model concept
Our model uses a 0-dimensional agent-based approach (Supporting Information Text S1) to simulate the life cycles of three major phytoplankton functional groups (dinoflagellates, diatoms, and diazotrophic cyanobacteria).Compartments for bulk zooplankton, nutrients, and detritus complete the ecosystem (Fig. 1).We simulate only one nutrient, nitrogen, which represents the most limiting nutrient in coastal ecosystems (Howarth and Marino 2006).Being the first to disentangle the separate effects of competition and adaptation on warmingrelated phytoplankton responses, we decided to restrict competition to a single nutrient, which produces clear results that are straightforward to interpret.Including more nutrients (e.g., phosphorus and silica) would unnecessarily increase the complexity of our model without helping us to answer our research question.
In addition to nitrogen, phytoplankton growth depends on photosynthetically available radiation and temperature.Dissolved inorganic nitrogen is taken up by all actively growing phytoplankton cells except for vegetative cyanobacteria cells that grow in filaments with nitrogen-fixing heterocysts.The nitrogen content of all dead phytoplankton and zooplankton cells is added to the detritus pool, which is remineralized back into bioavailable nitrogen at a constant rate.Due to constant sinking of detritus and stochastic burial of phytoplankton resting cells, nitrogen is lost from the system.Nitrogen loss is counteracted by the resuspension of previously buried resting cells and cyanobacterial nitrogen fixation.The probability of buried resting cells to be resuspended from the sediment decreases exponentially with time and is only possible between September and April due to vertical mixing.Thus, the nitrogen cycle in our model system is open, meaning that the amount of nitrogen in the system can change over time.All sources and sinks of nitrogen, as well as the nitrogen fluxes between the different ecosystem components, are shown in Fig. 1.A detailed model description is available in Supporting Information Text S2.
Our model simulates the dynamics of three important phytoplankton functional groups in the Baltic Sea, each represented by one common taxon or by a complex of common taxa.Following Hinners et al. (2019), we model the life cycle of a cold-water dinoflagellate species of the genus Apocalathium.We distinguish between two different life cycle stages, actively growing vegetative (VEG) cells and resting cysts.Resting cysts germinate after a dormancy period of several months (Kremp 2000); we simplify this process by prescribing germination to a certain period in early spring (between day 44 and day 60).The transfer from VEG cells back to cysts depends on temperature, with the encystment rate increasing strongly around a temperature threshold of 6 C. For diatoms, we chose a cold-water species of the genus Thalassiosira.Again, we distinguish between a resting stage (RES; spores) and a growing stage (VEG cells).Like Warns (2013), we assume that the transfer between stages depends on the actual growth conditions, i.e., the realized growth rate.For cyanobacteria, we consider a complex that represents the dominant nitrogen-fixing genera in the Baltic Sea, Nodularia, Aphanizomenon, and Anabaena (Karlsson et al. 2005).The cyanobacteria life cycle comprises three different stages given by VEG cells without heterocysts (nitrogen-limited), VEG cells with heterocysts (not nitrogenlimited), and resting cells (akinetes).We parameterize the transfer between the cyanobacteria life cycle stages by combining two modeling approaches (Hense andBeckmann 2006, 2010): Similar to diatoms, the transfer between stages depends on the actual growth conditions, which are either measured by the realized growth rate or by the severity of different growthlimiting factors.For VEG cells, we separate growth conditions into nitrogen limitation and limitation by temperature and light.While nitrogen limitation leads to heterocyst formation, unfavorable temperature and light conditions initiate the transfer from VEG cells to resting cells.Following other models (Hense and Beckmann 2006;Lee et al. 2018), we assume that cyanobacteria are non-grazeable due to toxicity, while dinoflagellates and diatoms are equally grazed by zooplankton.Further details on the simulated taxa and their life cycles can be found in Supporting Information Text S2.
In addition to phytoplankton life cycle dynamics, we consider two flexible temperature-dependent phytoplankton traits.The first flexible trait is the optimum temperature, which determines temperature limitation of growth.The optimum temperature is subject to random mutations, which we implement using an agent-based modeling approach (Supporting Information Text S1).In our model, each agent represents the same amount of biomass.As the three functional groups differ in for nitrogen (N), detritus (D), and zooplankton (Z), as well as agent-based life cycles of dinoflagellates (din), diatoms (dia), and cyanobacteria (cya).
Each life cycle comprises a resting (RES) and a growing stage (VEG cells).
For cyanobacteria, we consider a second, nitrogen-fixing growing stage (VEG cells with heterocysts, HET).Also shown are the nitrogen fluxes between the different ecosystem components, and the sources and sinks of nitrogen (cyanobacterial nitrogen fixation, resuspension and burial of resting cells, and sinking of detritus).The figure was created with Bio Render.com.
their cell sizes (Supporting Information Table S1), the number of cells per agent, i.e., resolution, differs among groups.Based on Beckmann et al. (2019), we assume that resolution affects the probability that a mutation occurs when an agent divides.
We calculate the mutation probability for one agent by multiplying the number of cells within the agent by the mutation probability for one cell, for which we use an experimentally derived value of 2.5 Â 10 À3 (Lenski and Travisano 1994).As our agents combine more than 2.5 Â 10 3 cells, mutations occur at every division, with the number of mutations being proportional to the number of cells within the agent.Thus, an agent combining more cells experiences more mutations per division.However, due to the larger number of cells within the agent, these mutations require the same time to be fixed as in smaller agents with less cells and less mutations.For this reason, we use the same mutational standard deviation (or mutational step size) for all agents independent of their resolution (0.1 C; Beckmann et al. 2019).These assumptions are similar to a modeling study by Merico et al. (2014), which uses a fixed mutational step size independent of population size, with larger populations having more mutations in total.
The second flexible trait, the cell size, responds plastically to temperature.The relation between cell size and temperature is inversely proportional, meaning that a temperature increase of 1 C causes a cell size decrease of $ 2.5% (Atkinson et al. 2003).
In our model, plastic responses occur during cell division and affect the newly produced daughter cell.Following Beckmann et al. (2019), we assume that a cell divides after having grown to a critical size.As Beckmann et al. (2019) do not consider cell size plasticity, they use the same critical size of 2b 0 for all cells, with b 0 being the initial cell size right after division.To account for cell size plasticity, we calculate the critical size of the newly produced daughter cell depending on environmental temperature.We assume that the daughter cell divides after having grown to 2b 0,T , with b 0,T being the initial size that the daughter cell would have if cell size responded instantaneously to temperature.We calculate b 0,T for the temperature during daughter cell production after Atkinson et al. (2003).In this way, an increase in temperature leads to a decrease in the daughter's critical size, reducing the initial size of the following generation.In addition, we consider that changes in cell size affect metabolic cell properties including maximum nitrogen uptake rate (Ward et al. 2017), half saturation constant for nitrogen (Litchman et al. 2007), basal cellular nitrogen requirement (Ward et al. 2017), maximum nitrogen storage capacity (Marañ on et al. 2013), and theoretical maximum metabolic rate (Ward et al. 2017).In our model, cell size thus determines the nitrogen-limited growth rate, which we calculate from internal nitrogen quotas using a variable-internal-stores model by Grover (1991) (Supporting Information Text S2).

Model setup and environmental forcing
For our simulations, we use a 0-dimensional model setup, which corresponds to a well-mixed tank.We assume that the sides of the tank are closed, while both top and bottom are open.This means that we neglect immigration and emigration of phytoplankton and zooplankton but still allow for changes in total mass due to interactions with the atmosphere (nitrogen fixation) and the benthic zone (sinking of detritus, burial of resting cells, and resuspension of resting cells).
At the beginning of a simulation, we initialize each phytoplankton functional group with the same biomass and overlapping generations.To create overlapping generations, we sample the biomass of each cell randomly from a uniform distribution between minimum cell biomass right after division (b 0,T , Supporting Information Table S1) and maximum cell biomass right before division (2b 0,T ).In contrast to cell biomass, optimum temperature is initially the same for all individuals of a functional group.In this way, we ensure that the observed adaptation results from mutations alone and not from initial intraspecific diversity.
Following other models (Hense and Beckmann 2006;Hinners et al. 2019), we simulate 12 equal months for each year (30 d per month) with a time step of 1 h.For a steady seasonal forcing, the model requires a spin-up period of $ 65 yr until phenology, taxa abundance, and mean traits have largely stabilized (Supporting Information Figs. S1 and S2).The steady seasonal forcing used for spin-up represents present-day conditions in the Gulf of Finland.For the seasonal variation of irradiance, we use the same function as Hinners et al. (2019), which was originally adapted from Stramska and Zuzewicz (2013).The function provides higher irradiance levels during summer and lower irradiance levels during winter.For the seasonal temperature forcing, we use a sinusoidal fit to sea surface temperature data from the Gulf of Finland.Prior to fitting, we averaged the data over the whole gulf area and over the last 10 yr (2011-2021) to compensate for regional differences and special temperature events such as heat waves.Sea surface temperature data were downloaded from the Copernicus database (https://resources.marine.copernicus.eu/products).The seasonal temperature and irradiance forcings are visualized in Supporting Information Fig. S3a; the corresponding equations are available in Text S2.

Model scenarios
To disentangle the separate effects of competition and adaptation on climate-driven phytoplankton responses, we simulate six different model scenarios (Table 1).
The first two model scenarios represent control scenarios (C: control; CA: control and adaptation), which we force with the present-day seasonal forcing for the Gulf of Finland over 100 yr.Both C and CA serve as spin-up for two global warming scenarios W (warming) and WA (warming and adaptation).To simulate global warming, we add a continuous temperature increase of 0.3 C per decade to the seasonal temperature forcing, which corresponds to the IPCC scenario SSP3-7.0 (Allan et al. 2021).This configuration is run for further 100 yr to simulate global warming during the next century (see Supporting Information Fig. S3b).While adaptation in the optimum temperature is disabled in C and W, it is enabled in CA and WA.This means that in C and W, we simulate ecological phytoplankton responses only, with competition (i.e., selection) as the single controlling factor.In CA and WA, on the contrary, phytoplankton responses to global warming are controlled by the combined effects of competition and adaptation (i.e., mutation and selection), meaning that they are eco-evolutionary.
To add further complexity to the representation of evolutionary processes in our model, we consider another essential mechanism that can both enhance and slow down adaptation (Rengefors et al. 2017): the resuspension of resting cells from previous blooms that were buried in the sediment.To study potential effects of resuspension on our model phytoplankton community, we simulate two additional scenarios termed CAR (control and adaptation and resuspension) and WAR (warming and adaptation and resuspension), where both adaptation and resuspension are enabled.For each of the six model scenarios (Table 1), we perform seven different simulations and average the output to ensure robust results.

Model validation for present-day conditions
To validate our model for present-day conditions, we use the two control scenarios C and CA; CAR is evaluated in a separate section (Effects of RES resuspension).Both C and CA produce the same seasonal succession of functional groups (Fig. 2): The blooming season starts with the spring bloom of dinoflagellates and diatoms, which reach maximum abundances in March and May, respectively.The spring bloom is followed by a summer bloom of cyanobacteria, which are most abundant in July.In autumn, diatoms form a second but weaker bloom with maximum abundances at the beginning of November.Our simulated seasonal succession of functional groups results from their individual responses to environmental growth conditions, including nitrogen concentration, temperature, and irradiance (Supporting Information Text S2 and Figs.S4-S6).
However, even if C and CA produce the same seasonal succession of functional groups, the timing of the individual blooms may differ by several days.Bloom timing, in this context, refers to the time when the bloom reaches its peak.In C, the spring bloom occurs later than in CA, with dinoflagellates and diatoms peaking $ 11 d and $ 6 d later, respectively (Table 2).The diatom autumn bloom, on the contrary, occurs $ 4 d earlier in C than in CA.In addition, C shows lower bloom amplitudes, especially for cyanobacteria, which Table 1.Overview of the six model scenarios that we analyze in this study.We perform seven different simulations for each scenario, with each simulation being run over 100 yr."Control" corresponds to a present-day seasonal temperature forcing for the Gulf of Finland, to which "warming" adds a steady temperature increase of 0.3 C per decade (IPCC scenario SSP3-7.0,Allan et al. 2021).produce a $ 42% weaker summer bloom.According to a t-test, the mentioned differences between C and CA are statistically significant (we use a significance level of 0.05 for t-tests, Supporting Information Table S2).Despite these differences, the two control scenarios C and CA agree reasonably well with recent monitoring data from the Baltic Sea (Hjerne et al. 2019), except for slight deviations in spring bloom timing (see Deviations in spring bloom timing in the Discussion section).Moreover, we find that in CA, the optimum temperatures of both dinoflagellates and diatoms track seasonal changes in environmental temperature following an initial acclimation period after germination (Supporting Information Fig. S7).

Adaptation and resuspension
Dinoflagellates show the highest adaptation rate with a change in mean optimum temperature by $ 0.18 C within 1 month.This simulated adaptation rate agrees well with observed adaptation rates for Chaetoceros tenuissimus (Jin and Agustí 2018).

Ecological phytoplankton responses to global warming
In the two warming scenarios W and WA, our phytoplankton community shows different ecological responses to global warming, including shifts in bloom timing and changes in bloom amplitude (Fig. 2).Shifts in bloom timing are only notable for dinoflagellates and diatoms but not for cyanobacteria.For dinoflagellates and diatoms, blooms are shifted toward the winter period, with the shift being stronger in W than in WA.Dinoflagellates bloom $ 8 d earlier in W compared with C, while the shift between WA and CA only amounts to $ 3 d (Table 3).Diatoms show an even stronger shift in bloom timing than dinoflagellates: while the spring bloom is shifted by $ 25 d and $ 17 d, the autumn bloom is shifted by $ 17 d and $ 8 d in W and WA, respectively.
Apart from shifts in bloom timing, W and WA produce increasing bloom amplitudes for all functional groups, except for dinoflagellates.Again, the observed changes are stronger in W than in WA.Cyanobacteria show the strongest amplitude increase with $ 223% in W and $ 57% in WA, followed by diatoms.For diatoms, the autumn bloom increases more than the spring bloom.Dinoflagellates, however, increase only slightly in W, while they decrease in WA; this is in line with the overall stronger amplitude increase in W.
According to a t-test, bloom timing and bloom amplitudes differ significantly between W and WA for all functional groups; the only exception is the timing of the cyanobacterial summer bloom (Supporting Information Table S3).

Adaptive phytoplankton responses to global warming
Already in the control scenario CA, we observe notable adaptation of all functional groups.Dinoflagellates adapt to notably lower temperatures in CA (Fig. 3), with most of the adaptation occurring within the first $ 60 yr of the simulation period (Supporting Information Fig. S2).Throughout the entire CA simulation period of 100 yr, dinoflagellates decrease their mean optimum temperature (T opt ) by 5.8 C, from an initial T opt of 10.8 C (Hinners et al. 2017) to a new T opt of 5.0 C (Table 2).Diatoms and cyanobacteria also show a tendency towards lower temperatures; this tendency is considerably weaker than for dinoflagellates, though.Both diatoms and cyanobacteria reduce their initial optimum temperatures of 12 C (Spilling 2011) and 28.5 C (Collins and Boylen 1982;Lehtimäki et al. 1997) by 0.7 C throughout CA.At the end of CA, diatoms and cyanobacteria thus grow at new optimum temperatures of 11.3 C and 27.8 C, respectively.While most of the adaptation of diatoms occurs during the first $ 65 yr of the simulation period, which is similar to dinoflagellates, cyanobacteria show steady adaptation over the 100 yr of simulation (Supporting Information Fig. S2).Overall, adaptive trait changes in CA correlate temporally with ecological changes in bloom timing and taxa abundance (Supporting Information Figs.S1 and S2).
Table 2. Average bloom timing, peak abundance, and optimum temperature (T opt ) of all taxa for the three control scenarios C (control), CA (control and adaptation), and CAR (control and adaptation and resuspension), including standard deviations.For each scenario, averages were calculated from the last simulation year of seven different simulations.

C
CA CAR In the warming scenario WA, we observe notable adaptation of all functional groups to the increasing temperatures.Dinoflagellates increase their mean T opt by 1.4 C throughout WA (Table 3).The resulting T opt of 6.4 C is still notably lower than the initial T opt of dinoflagellates at the beginning of CA (10.8 C, Table 2).Diatoms and cyanobacteria, on the contrary, increase their mean T opt above their initial optimum temperatures.Throughout WA, the mean T opt of diatoms and cyanobacteria increases by 1.5 C and 1.2 C, respectively, giving final optimum temperatures of 12.7 C and 29.1 C. Thus, the final optimum temperatures of diatoms and cyanobacteria at the end of WA are by $ 0.7 C and $ 0.6 C higher than their initial optimum temperatures at the beginning of CA, respectively.

Effects of RES resuspension
When we enable the resuspension of resting cells from the sediment, we observe only slight phenological changes in  the control scenario CAR compared with CA (Fig. 2 and Table 2).While taxa abundances are comparable, autumn diatoms bloom $ 2 d earlier in CAR; this shift in bloom timing is statistically significant according to a t-test (Supporting Information Table S2).Similarly, we found autumn diatoms to bloom $ 4 d earlier in the control scenario without adaptation C when compared with CA (see section Model validation for present-day conditions).
The slight phenological changes in CAR are accompanied by reduced thermal adaptation.While in both CA and CAR, optimum temperatures evolve toward lower values, the final optimum temperatures at the end of CAR are between 0.14 C and 0.39 C higher than at the end of CA (Table 2).Changes in adaptation are, however, only statistically significant for cyanobacteria, for which adaptation is reduced the most (Supporting Information Table S2).
Under warming conditions, we observe stronger phenological effects of resuspension as part of a general weakening of adaptation-related effects.The results of the WAR scenario resemble more closely the results of the W than the WA scenario (Fig. 2).While the two diatom blooms tend to be shifted $ 1-2 d more toward the winter period, the cyanobacterial summer bloom is by $ 18% stronger in WAR than in WA.The difference in cyanobacteria biomass is statistically significant (Supporting Information Table S3).
The stronger effects of resuspension under warming conditions are accompanied by a greater slowdown of adaptation (Fig. 3 and Table 3).In WAR, the three functional groups increase their mean optimum temperatures between 0.20 C and 0.54 C less than in WA.The resulting differences in optimum temperature are, however, only statistically significant for dinoflagellates (Supporting Information Table S3).In addition, we observe a slight narrowing of the T opt trait distribution in WAR for both dinoflagellates and cyanobacteria, for which adaptation is slowed more than for diatoms.

Discussion
In this study, we investigated the interplay between phytoplankton competition and adaptation under global warming using an eco-evolutionary ecosystem model.We found that adaptation reduces warming-related ecological phytoplankton responses in terms of changes in bloom timing and amplitude.Resuspension of resting cells can slow down adaptation and hence enhance changes in bloom timing and amplitude.To the best of our knowledge, this is the first ecosystem model that simulates competition and adaptation of multiple phytoplankton functional groups under global warming, allowing us to disentangle how the interplay of both processes affects warming-related phytoplankton responses on realistic time scales.Due to the simplicity of the model design, we will not discuss exact magnitudes in detail but focus on qualitative differences.

Deviations in spring bloom timing
Contrary to recent monitoring data from the Baltic Sea (Hjerne et al. 2019), dinoflagellates bloom earlier in spring than diatoms in our simulations.One possible explanation for this could be that our model does not consider re-stratification during spring, which plays an important role during spring bloom formation in the Baltic Sea (Klais et al. 2011).While motile dinoflagellates require stratified conditions to form a bloom (Margalef et al. 1979), nonmotile diatoms are favored under turbulent mixing (María Trigueros and Orive 2001).In addition, the monitoring data do not distinguish between species but only functional groups, while we simulate a single species from each functional group (except for cyanobacteria).The traits of our chosen species do not necessarily coincide with the bulk traits of the functional group that the species belongs to.Consequently, the species we chose may not bloom at the same time as the bulk of species in the corresponding functional group.Indeed, lower optimum temperatures were measured for the cold-water dinoflagellate Apocalathium malmogiense than for the cold-water diatom Thalassiosira baltica (Spilling 2011;Hinners et al. 2017).As we simulate a dinoflagellate species of the genus Apocalathium and a diatom species of the genus Thalassiosira, we use these optimum temperatures to initialize our model.The earlier bloom of dinoflagellates in our simulations is hence a result of the physiological properties of Apocalathium, including the initial optimum temperature and in particular the strong increase in encystment at low temperatures ($ 6 C).Thus, we can conclude that for our chosen phytoplankton species, our model produces a realistic seasonal cycle.

Eco-evolutionary dynamics in a steady environment
The results of control scenario C differ notably from the two control scenarios CA and CAR, which themselves are comparable.Contrary to C, both CA and CAR include adaptation, meaning that the observed differences must result from adaptation-related effects.In comparison to C, CA and CAR show slightly shifted blooms of higher amplitude.Adaptation in CA and CAR allows each functional group to improve its growth conditions by utilizing its niche more efficiently.
Dinoflagellates improve their growth conditions by adapting to lower temperatures.Additional CA simulations using only dinoflagellates reveal that this adaptation pattern is not driven by interspecific competition with diatoms but by intraspecific competition for nitrogen, which selects for early bloomers (Supporting Information Fig. S8).Even in the absence of interspecific competition, dinoflagellates do not extend their niche towards higher temperatures due to restrictions in their life cycle.Resurrection experiments demonstrated that the temperature threshold for encystment of A. malmogiense remained constant at around 6 C over the last century of global warming (Hinners et al. 2017), prohibiting a shift toward later blooming.As our simulations cover a similar time span, we implemented the encystment threshold as a fixed trait.
Contrary to dinoflagellates, adaptation of diatoms is controlled by interspecific competition with cyanobacteria.As inferior competitors under warm, nitrogen-limited conditions, diatoms lower their optimum temperature to avoid competition for nitrogen.This hypothesis is confirmed by additional CA simulations using only diatoms.These simulations show that without interspecific competition, diatoms adapt to higher temperatures in a steady environment to merge their spring and autumn blooms (Supporting Information Fig. S9).Thus, we can conclude that interspecific competition between diatoms and cyanobacteria drives niche separation.This finding goes in line with a modeling study by Barab as and D'Andrea (2016), who found that competing species place themselves more evenly across the trait axis if they can adapt.
Contrary to dinoflagellates and diatoms, adaptation of cyanobacteria is neither driven by intra-nor by interspecific competition for nitrogen, which is not surprising given their ability to fix atmospheric nitrogen.Instead, the difference between environmental and optimum temperature is the driving factor, with the initial T opt of cyanobacteria being almost 9 C higher than maximum environmental temperature in summer.Despite the large temperature difference, cyanobacteria reduce their mean optimum temperature by only 0.7 C throughout CA.This is most likely caused by the wide, plateau-shaped thermal reaction norm of cyanobacteria (Supporting Information Fig. S4), which does not lead to a strong selection pressure on the optimum temperature even if environmental temperature differs notably from the optimum.Indeed, the growth rate of initial genotypes is only 17% lower than the growth rate of genotypes that are adapted optimally to environmental temperature.Our hypothesis is further supported by the T opt trait distribution of cyanobacteria, which is the broadest among the functional groups in our model (Supporting Information Fig. S10).Simulations by Beckmann et al. (2019) showed that trait diversity is reduced under strong selection pressure.Thus, in our simulations, selection pressure indeed seems to be weaker for cyanobacteria than for dinoflagellates and diatoms.The weak selection pressure on cyanobacteria prevents them from completing adaptation to control conditions within the 100 yr of simulation, meaning that we still observe transient trait dynamics at the end of CA.Dinoflagellates and diatoms, on the contrary, complete adaptation in CA within < 100 yr.As neither dinoflagellates nor diatoms are implemented with larger mutation rates or larger mutational step sizes than cyanobacteria, their faster adaptation must result from the stronger selection pressure they experience.We propose that the stronger selection pressure on dinoflagellates and diatoms results from their narrower thermal reaction norms (Supporting Information Fig. S4).For comparison, cyanobacteria can reach 90% of their maximum growth rate over a temperature range of 9.4 C. For dinoflagellates and diatoms, this range is considerably narrower with 4.7 C and 4.2 C, respectively.To conclude, selection pressure and adaptation speed are not only affected by the difference between optimum and environmental temperature but also by the width and shape of the thermal reaction norm.
From the above discussion, we can conclude that adaptation allows all functional groups to improve their individual growth conditions to a certain degree, with the drivers of adaptation differing among groups.However, adaptation of one group may also affect the other groups through changes in competition.For example, diatom adaptation reduces competition with cyanobacteria, which is beneficial for both sides.Reduced competition with diatoms leaves more nitrogen for the nitrogen-limited vegetative cyanobacteria life cycle stage, which can hence initiate a larger summer bloom.Indeed, additional CA simulations without cyanobacterial adaptation reveal that reduced competition with diatoms causes 71% of the amplitude increase in CA compared with C, while only the remaining 29% result from adaptation of cyanobacteria to colder temperatures.As a result, reduced competition with diatoms indirectly enhances the bloom of vegetative cyanobacteria with heterocysts and thus increases the amount of fixed atmospheric nitrogen.The higher nitrogen input into the system in summer allows for stronger blooms of dinoflagellates and especially diatoms, which can directly take up the newly available nitrogen in autumn.In conclusion, adaptation in a steady environment induces positive feedback in our simulations: While diatoms actively improve their growth conditions through niche separation, cyanobacteria benefit passively from reduced competition with diatoms.The higher nitrogen availability due to enhanced cyanobacterial nitrogen fixation feeds back positively on both diatoms and dinoflagellates.

Eco-evolutionary dynamics in a changing environment
The two warming scenarios with and without adaptation WA and W produce ecological phytoplankton responses of notably different intensity, which nevertheless follow the same trends.Both scenarios show a shift of dinoflagellates and diatoms toward the winter period, as well as a notable increase in cyanobacteria biomass in summer.As these findings are consistent with trends observed over the past 50 yr of global warming in the Baltic Sea (Suikkanen et al. 2007;Wasmund et al. 2019), we believe them to be realistic.
By comparing the relative magnitude of warming-related changes between W and WA, we can disentangle the effects of competition and adaptation in a warming environment.In a steady environment, we found adaptation to enhance bloom amplitudes; in a warming environment, however, adaptation has the opposite effect.For increasing temperatures, amplitudes are lower in WA, where adaptation is enabled, especially for cyanobacteria.To understand if the weaker cyanobacteria bloom in WA is mainly caused by cyanobacterial thermal adaptation or adaptation-related differences in competition, we performed additional simulations in which we disabled the thermal adaptation of cyanobacteria.The simulations reveal that cyanobacterial thermal adaptation leads to an increase in bloom amplitude of 2% and hence only has a minor effect.Thus, the observed differences in cyanobacterial bloom amplitude must result from adaptation-related differences in competition with diatoms.Due to lacking adaptation in W, diatoms (and dinoflagellates) can only compensate the increase in temperature by shifting their blooms toward the winter period, when temperatures are closer to their fixed optimum temperature.Thus, diatoms are less abundant and less competitive in late spring and leave more nitrogen for the nitrogen-limited vegetative cyanobacteria, which can initiate a larger summer bloom.In WA, on the contrary, diatoms increase their optimum temperature above their fixed T opt in W.However, this shift is lower than the shift in environmental temperature.In the absence of interspecific competition, diatoms track changes in environmental temperature more closely as confirmed by additional WA simulations using only diatoms (Supporting Information Fig. S9).This suggests that cyanobacteria restrict diatom adaptation to increasing temperatures, which matches findings from a modeling study by De Mazancourt et al. (2008). De Mazancourt et al. (2008) showed that the presence of preadapted species (in our case cyanobacteria) restricts adaptation of other species.Even if cyanobacteria restrict diatom adaptation in WA, diatoms manage to increase their optimum temperature above their fixed T opt in W, which reduces the shift in bloom timing, and increases the abundance and competitivity of diatoms in late spring.The enhanced competition with diatoms reduces nitrogen availability for nitrogen-limited vegetative cyanobacteria and leads to a weaker cyanobacterial summer bloom.
In conclusion, our warming simulations reveal that adaptation reduces ecological phytoplankton responses to global warming in terms of changes in bloom timing and amplitude.This finding matches findings from previous theoretical modeling studies.For example, Barab as and D'Andrea (2016) found that communities that can evolve are much more robust against environmental change.Similarly, Northfield and Ives (2013) found that if global warming leads to conflicting interests between species (i.e., stronger competition), coevolution reduces warming-related effects and hence changes in species abundance.

Effects of RES resuspension
In our simulations, the resuspension of resting cells (or propagules) tends to slow down adaptation to global warming, which is expressed by slightly stronger shifts in bloom timing and a stronger cyanobacterial summer bloom.We can explain the observed slowing effect of resuspension by the re-introduction of past-adapted resting cells to the actively growing population.A similar slowing effect was observed for lake copepods (Hairston Jr and De Stasio Jr 1988).However, some studies suggest that propagule banks can aid phytoplankton to cope with environmental change as they contain a high diversity of genotypes (Ribeiro et al. 2013;Kremp et al. 2016).The effect of resuspension on adaptation seems to depend on whether the resuspended genotypes are a random or nonrandom sample of the total gene pool (Hairston Jr and De Stasio Jr 1988;Rengefors et al. 2017).For layered sediments, this suggests that the resuspended sample of the total gene pool, and hence the effect of resuspension on adaptation, are affected by the adaptation history of the population.
In our study, propagule banks are formed in an initially steady environment that transforms into a steadily warming environment.Under steady conditions, all taxa adapt to lower temperatures, with adaptation being weakened by resuspension.Thus, when temperatures begin to increase, initial optimum temperatures are higher in the warming scenario with resuspension WAR than in the warming scenario without resuspension WA.Due to the previous adaptation to colder temperatures, the sediment additionally contains genotypes that are adapted to even higher temperatures than the actively growing population at the beginning of WAR.This raises the question of why resuspension slows down adaptation to higher temperatures instead of enhancing it.During the 100 yr of steady environmental conditions in CAR, most adaptation to colder temperatures occurs during the first $ 70 yr (Supporting Information Fig. S11), meaning that warm-adapted genotypes are buried deep in the sediment.As resuspension probability decreases exponentially with time in our model, the number of resuspended warm-adapted genotypes is negligible compared with the number of cold-adapted genotypes that are resuspended.Thus, the slowing effect of resuspension dominates in our warming simulations.Indeed, the higher optimum temperatures in the beginning of WAR are outpaced after $ 60-80 yr by the faster adaptation in WA, leading to the lag in adaptation that we observe after 100 yr of warming (Supporting Information Fig. S11).The lag in adaptation in WAR leads to a stronger shift of the two diatom blooms toward winter, which reduces the competition between spring diatoms and cyanobacteria and hence allows for a stronger cyanobacterial summer bloom.
In addition, resuspension narrows the T opt trait distributions of both dinoflagellates and cyanobacteria, which are much more affected by slowed adaptation than diatoms.We suggest that the lesser effect of resuspension on diatoms results from the two diatom blooms per year, which may allow for faster adaptation.Indeed, a comparison of the additional CA simulations using only dinoflagellates with those using only diatoms reveals that diatoms have the potential to adapt much faster than dinoflagellates when not restricted by interspecific competition (Supporting Information Fig. S12).During the first 4 yr of simulation, when the diatom spring and autumn blooms have not merged yet, diatoms adapt more than two times faster than dinoflagellates.
To conclude, our results show that resuspension of resting cells tends to slow down adaptation.Our results, however, strongly depend on adaptation history.We suggest for future work to further investigate the effect of adaptation history on resuspension.

Model biases and suggestions for future work
Below, we describe the major assumptions and simplifications that we made in our model and discuss which model biases they may imply.In addition, we give suggestions for future work.
It is a major challenge to parameterize the transfer between different life cycle stages, for which we need information on triggers and rates of transition.As both parameters are poorly constrained, we must estimate them though model calibration.However, the parameters estimated for present-day conditions may not be applicable to a global warming scenario.As revealed by resurrection experiments, the cold-water dinoflagellate A. malmogiense showed a decrease in the encystment rate over the past century of global warming, while the temperature threshold for encystment remained unchanged (Hinners et al. 2017).Consequently, different life cycle traits seem to respond differently to environmental change, meaning that further experimental research is needed to improve the model parameterization of phytoplankton life cycles.
Another potential model bias affects our mutational algorithm.We simulate trait changes due to random mutations (and sexual recombination) by randomly sampling the modified trait value from a Gaussian distribution centered at the original trait value.The width of this distribution, i.e., the standard deviation of mutations, cannot be measured directly in the laboratory and must hence be estimated through model calibration.For our simulations, we used the same value as Beckmann et al. (2019) (see section Model concept).In our warming scenario WA, dinoflagellates, diatoms, and cyanobacteria show an average increase in T opt of $ 0.14 C per decade for a steady temperature increase of 0.3 C per decade.As evolution experiments can neither replicate the complexity of marine ecosystems nor the time scales of global warming, we compare our simulated adaptation rates to observations.Evaluating 15 yr of observational data for 67 species, Irwin et al. (2015) found a mean increase in the realized niche of 0.45 C for a mean temperature increase of 0.73 C, giving an average adaptation rate of 0.3 C per decade.Thus, our simulations show adaptation rates of the same order of magnitude.The slightly slower adaptation in our simulations may result from the slower temperature increase (0.3 C per decade vs. 0.45 C per decade), leading to a lower selection pressure.For higher selection pressures, e.g., over the course of the seasonal cycle, we observe higher adaptation rates in the range of $ 0.18 C within 1 month for dinoflagellates (Supporting Information Fig. S7), which suggests that our chosen mutational standard deviations are not artificially limiting adaptation to temperature changes.
Being the focus of our study, phytoplankton adaptation is simulated with a high degree of complexity, while adaptation of zooplankton is not taken into account.As we consider zooplankton merely for the sake of grazing mortality, we use a simplistic representation, where we ignore potential zooplankton adaptation, life cycle dynamics, and effects of irradiance and temperature on grazing.Thus, the zooplankton in our model is entirely controlled by phytoplankton, which leads to comparable responses to global warming in terms of increasing abundances and shifts in bloom timing.Observations showed that the abundance of herbivorous zooplankton is strongly controlled by their phytoplankton prey (Richardson and Schoeman 2004), while phenological shifts synchronized with phytoplankton are mainly associated with fast-growing zooplankton (Adrian et al. 2006;Dam 2013).This suggests that our representation of zooplankton may be reasonable for fast-growing taxa like Daphnia, while our model is inappropriate to represent slow-growing taxa with more complex life cycles like copepods (Adrian et al. 2006).We expect that an explicit simulation of the adaptation of fast-growing zooplankton would not change our results notably.Adaptation of slow-growing zooplankton, however, may have a larger effect, but we still do not expect qualitative changes.Future work can expand on our model and include a more complex representation of zooplankton, including light-and temperaturedependent grazing, life cycle dynamics, and adaptation.
Another bias is that our model cannot capture the functional diversity of the real-world Baltic Sea phytoplankton community.Our model community consists of three major functional groups (dinoflagellates, diatoms, and cyanobacteria), each represented by one taxon or by an averaged complex.Thus, we do not account for other important groups such as ciliates or other flagellates (Thamm et al. 2004) and for the functional diversity within each group.Most likely, across-species functional diversity within a functional group greatly enhances the group's potential to adapt to environmental changes (Litchman et al. 2015).With that in mind, we might be missing other key competitors like summer-blooming dinoflagellate or diatom species.Based on our own findings, the absence of these competitors most likely influences both phenology and adaptation of our focal phytoplankton taxa.For the above reasons, we do not interpret our results quantitatively but focus on analyzing qualitative differences between scenarios and identifying general eco-evolutionary mechanisms.We suggest that future work builds on our findings and investigates which role diversity within functional groups plays in the Baltic Sea phytoplankton community regarding both competition and adaptive potential.
Finally, as mentioned in section "Deviations in spring bloom timing," our model is currently not coupled to a physical ocean model.Thus, we ignore physical processes like vertical mixing and stratification, which may, e.g., affect our simulated species succession during phytoplankton spring bloom.In addition, our 0-dimensional model setup leads to an unrealistic adaptation pattern of cyanobacteria.Under warming conditions, preadapted cyanobacteria adapt to even higher temperatures and hence evolve away from environmental temperature.We can explain this apparent maladaptation by our 0-dimensional representation of competition.Nitrogen concentration is low at the end of the diatom spring bloom when vegetative cyanobacteria start to grow.As warm conditions are favorable for cyanobacteria, they grow faster and become nitrogen-limited more quickly than under control conditions.With nitrogen recovering over time due to nitrogen fixation and remineralization, the environment selects for later-blooming vegetative cyanobacteria with higher optimum temperatures.We expect that by resolving the spatial distribution of nitrogen in a 1D or 3D physical environment, our model would simulate a more realistic adaptation pattern of cyanobacteria.Previous 1D models showed that diatoms may shift their spring bloom to deeper water layers when nitrogen gets depleted at the surface (Warns 2013;Lee et al. 2018).As cyanobacteria can grow more efficiently than diatoms on low nitrogen due to their smaller size (Litchman et al. 2007), they could take over growth at the surface on the remaining nitrogen without competing directly with diatoms.As a result, cyanobacteria would no longer experience an unrealistic selection pressure toward higher temperatures.Still, we think that the coupling to a physical model would not change our main results and conclusions.
Based on the above discussion, we can conclude that a physical environment would make our model more realistic.In turn, however, our model would also be a valuable addition to existing coupled biological-physical models.Existing models for the Baltic Sea project a notable increase in cyanobacteria biomass in the future if nutrient loads are not reduced (Meier et al. 2012;Hense et al. 2013).As these models ignore phytoplankton adaptation, they may systematically overestimate warming-related changes in cyanobacteria biomass.Integrating our model into a coupled biological-physical model could hence notably improve climate projections.As long-term evolutionary data (e.g., from sediment archives) become available, these will moreover represent a valuable source for model calibration.As improving climate projections is relevant for political decision making, we suggest for future work to test our or a similar approach in a coupled biologicalphysical environment.

Conclusions
Our modeling study demonstrates that the combined effects of phytoplankton competition and adaptation may differ from their separate effects.Both processes influence each other and shape phytoplankton communities through their interplay.The outcome of this interplay depends on environmental conditions.
In a steady environment, each functional group occupies a fixed niche in which it can improve its growth conditions through adaptation.The main driver of adaptation is given by the most limiting factor for growth.For nitrogen fixers, temperature appears to be the main driver, while it is nitrogen for nitrogen-limited groups.Nitrogen-limited groups maximize their growth rate by adapting to the temperature that coincides with the highest nitrogen concentration within their niche.If groups are inferior competitors, they adapt to avoid interspecific competition for nitrogen.If groups are superior competitors or if no competitor is present, adaptation is driven by intraspecific competition for nitrogen.In a changing environment, however, niches are no longer fixed and thus, groups can no longer avoid competition.Adaptation allows inferior competitors to increase their competitivity and hence to mitigate the dominance of superior competitors.
In conclusion, the role of adaptation cannot be neglected in changing environments, otherwise warming-related changes in taxa dominance (and bloom timing) can be systematically overestimated.To realistically simulate phytoplankton responses to environmental changes, future ecosystem models must consider both competition and adaptation.Moreover, our results demonstrate that adaptation can be affected by dormant resting cells that are resuspended from the sediment.Our study highlights that eco-evolutionary ecosystem models represent a powerful tool to estimate phytoplankton responses and even ecosystem responses to global warming.Future work can build on our model and expand on our representation of competition, e.g., by including multiple limiting nutrients and/or a physical environment.In addition, future models can integrate emerging evolutionary data based on long-term data series (e.g., from sediment archives).The integration of such data will allow us to further improve projections of future ecosystem changes in response to anthropogenic environmental change.

Fig. 1 .
Fig. 1.Ecosystem components of our model including compartments

Fig. 2 .
Fig. 2. Biomass of phytoplankton growing stages during the last simulation year for all six model scenarios (C: control; W: warming; CA: control and adaptation; WA: warming and adaptation; CAR: control and adaptation and resuspension; WAR: warming and adaptation and resuspension).For each model scenario, seven different simulations were averaged.The dashed lines in each warming scenario indicate the timing of the bloom peaks in the corresponding control scenario.Please note that the biomass of both cyanobacteria growing stages (VEG cells without heterocysts and VEG cells with heterocysts) is summarized.

Fig. 3 .
Fig. 3. Evolution of the optimum temperature through mutation and selection.Shown are the initial value, the preliminary trait distribution after 100 yr of present-day seasonal forcing, and the final trait distribution after 100 yr of warming.Upper panels: adaptation only, i.e., CA (control and adaptation) and WA (warming and adaptation); lower panels: adaptation plus resuspension, i.e., CAR (control and adaptation and resuspension) and WAR (warming and adaptation and resuspension).For each model scenario, trait distributions were merged from monthly histograms of the last simulation year of seven different simulations.Please note that the temperature range on the x-axis differs between functional groups.

Table 3 .
Average warming-related changes in bloom timing, peak abundance, and optimum temperature (T opt ) for all taxa and model scenarios, including propagated errors.Scenario abbreviations: C = control; W = warming; CA = control and adaptation; WA = warming and adaptation; CAR = control and adaptation and resuspension; WAR = warming and adaptation and resuspension.Negative changes in bloom timing indicate a shift toward earlier in the year.Please note that changes in taxa abundance are not given as absolute values like in Table2but as relative changes.According to a t-test, all changes are statistically significant at the 0.05 level.For details, see Supporting Information TableS4.