Factors controlling the competition between Phaeocystis and diatoms in the Southern Ocean

The high-latitude Southern Ocean phytoplankton community is shaped by the competition between Phaeocystis and silicifying diatoms, with the relative abundance of these two groups controlling primary and export production, the production of dimethylsulfide, the ratio of silicic acid and nitrate available in the water column, and the structure of the food web. Here, we investigate this competition using a regional physical-biogeochemical-ecological model (ROMS-BEC) configured at eddypermitting resolution for the Southern Ocean south of 35◦ S. We extended ROMS-BEC by an explicit parameterization of 5 Phaeocystis colonies, so that the model, together with the previous addition of an explicit coccolithophore type, now includes all biogeochemically relevant Southern Ocean phytoplankton types. We find that Phaeocystis contribute 46% and 40% to annual NPP and POC export south of 60◦ S, respectively, making them an important contributor to high-latitude carbon cycling. In our simulation, the relative importance of Phaeocystis and diatoms is mainly controlled by the temporal variability in temperature and iron availability. The higher light sensitivity of Phaeocystis at low irradiances promotes the succession from Phaeocystis to 10 diatoms in more coastal areas, such as the Ross Sea. Still, differences in the biomass loss rates, such as aggregation or grazing by zooplankton, need to be considered to explain the simulated seasonal biomass evolution.


Introduction
Unused nutrients from the Southern Ocean (SO) fuel global primary and export production (e.g. Sarmiento et al., 2004;Palter et al., 2010), and the amount and stoichiometry of these laterally exported nutrients is determined by multiple types of 15 phytoplankton with different nutrient requirements. As on-going climate change alters the relative contribution of different phytoplankton groups to total net primary production (NPP; IPCC, 2014; Constable et al., 2014;Deppeler and Davidson, 2017), this will hence have ramifications for global biogeochemical cycles and food web structure (Smetacek et al., 2004).
Today, the SO phytoplankton community is largely dominated by silicifying diatoms (e.g. Swan et al., 2016), but the contributions of calcifying coccolithophores and dimethylsulfide (DMS) producing Phaeocystis are substantial in the subantarctic 20 (Balch et al., 2016;Nissen et al., 2018) and in the high latitudes, respectively (Smith and Gordon, 1997;Arrigo et al., 1999;Poulton et al., 2007), with a likely increase in the relative importance of the latter two types in a warming world (Bopp et al., 2005;Winter et al., 2013;Rivero-Calle et al., 2015). On a global scale, Phaeocystis has been suggested to be a major player in the marine cycling of DMS (e.g. Keller et al., 1989;Liss et al., 1994) and to contribute 6-65% to total phytoplankton carbon biomass (Buitenhuis et al., 2013b). Yet, the contribution of Phaeocystis to the export of particulate 25 2 Methods

ROMS-BEC with explicit Phaeocystis colonies
We use a quarter-degree SO setup of the Regional Ocean Modeling System ROMS (latitudinal range from 24 • S-78 • S, 64 topography-following vertical levels, time step to solve the primitive equations is 1600 s; Shchepetkin and McWilliams, 2005;Haumann, 2016), coupled to the biogeochemical model BEC (Moore et al., 2013), which was recently extended to include an 100 explicit represenation of coccolithophores and thoroughly validated in the SO setup (Nissen et al., 2018). BEC resolves the biogeochemical cycling of all macronutrients (C, N, P, Si), as well as the cycling of iron (Fe), the major micronutrient in the SO.
The model includes four phytoplankton functional types (PFT) -diatoms, coccolithophores, small phytoplankton/SP, and N 2fixing diazotrophs -and one zooplankton functional type (Moore et al., 2013;Nissen et al., 2018). Here, we extend the version of Nissen et al. (2018) to include an explicit parameterization of colonial Phaeocystis antarctica, which is the only species 105 of Phaeocystis occurring in the SO . For the remainder of this manuscript, we will refer to the new PFT as "Phaeocystis". Generally, model parameters for Phaeocystis are chosen to represent the colonial form of Phaeocystis whenever information is available in the literature (see e.g. review by Schoemann et al., 2005). By only simulating the colonial form of Phaeocystis, we assume enough solitary cells of Phaeocystis to be available for colony formation at any time as part of the SP PFT. As for the other phytoplankton PFTs, growth by Phaeocystis is limited by surrounding temperature, nutrient, 110 and light conditions as outlined in the following (for a complete description of the model equations describing phytoplankton growth, see Nissen et al., 2018).
As the new PFT in ROMS-BEC represents a single species of Phaeocystis, we use an optimum function to describe its temperature-limited growth rate µ PA (T ) (d −1 , Schoemann et al., 2005): (1) 115 In the above equation, the maximum growth rate (µ PA max ) is 1.56 d −1 at an optimum temperature (T opt ) of 3.6 • C and the temperature interval (τ ) is 17.51 • C and 1.17 • C at temperatures below and above 3.6 • C, respectively. With these parameters, the simulated growth rate of Phaeocystis in ROMS-BEC is zero at temperatures above ∼8 • C (in agreement with laboratory experiments with Phaeocystis Antarctica, see Buma et al., 1991) and higher than that of diatoms for temperatures between ∼0-4 • C (Fig. A1a). We acknowledge that the range of temperatures for which the growth of Phaeocystis exceeds that of diatoms is 120 possibly too small, as the temperature-limited growth rate by diatoms in ROMS-BEC is too high at low temperatures compared to available laboratory data (see Fig. A1a). Yet, we note that temperature-limited growth by diatoms in the model is tuned to fit the data at the global range of temperatures, in particular for the competition with coccolithophores at subantarctic latitudes (Nissen et al., 2018).
Half-saturation constants for macronutrient limitation are scarce for P. antarctica , and macronutri-125 ent limitation of Phaeocystis is therefore chosen to be identical to that of diatoms in ROMS-BEC (Table 1). As the availability of the micronutrient Fe generally limits phytoplankton growth in the high-latitude SO (Martin et al., 1990a, b) and accordingly in ROMS-BEC (Fig. S1), this choice is not expected to significantly impact the simulated competition between diatoms Table 1. BEC parameters controlling phytoplankton growth and loss for the five phytoplankton PFTs diatoms (D), Phaeocystis (PA), coccolithophores (C), small phytoplankton (SP), and diazotrophs (N). Z=zooplankton, P=phytoplankton, PI=photosynthesis-irradiance. If not given in section 2.1, the model equations describing phytoplankton growth and loss rates are given in Nissen et al. (2018).

Parameter
Unit  Nissen et al. (2018), the k Fe of diazotrophs in ROMS-BEC is higher than for all other PFTs, consistent with literature reporting high Fe requirements of Trichodesmium (Berman-Frank et al., 2001). Furthermore, the maximum grazing rate on diazotrophs is lowest in the model (Capone, 1997).
Still, diazotrophs continue to be a minor player in the SO phytoplankton community, contributing <1% to domain-integrated NPP in ROMS-BEC. ‡ The temperature-limited growth rate of Phaeocystis is calculated based on an optimum function according to Eq. 1 (see also Fig. A1a).
and Phaeocystis in this area. In contrast, differences in the half-saturation constants with respect to dissolved Fe concentrations (k Fe ) of Phaeocystis and diatoms possibly critically impact the competitive success of Phaeocystis relative to diatoms 130 throughout the year (see e.g. Sedwick et al., 2000Sedwick et al., , 2007. Here, due to their larger size, we assume a higher k Fe for Phaeocystis Table 1). We note however, that the k Fe of Phaeocystis has been reported to vary over one order magnitude depending on the ambient light level (0.045-0.45 µmol m −3 , see Fig. A1b and Garcia et al., leads to a better agreement of the k Fe used in ROMS-BEC with values suggested for large SO diatoms by Timmermans et al. (2004), but we acknowledge that the chosen value here is still at the lower end of their suggested range (0.19-1.14 µmol m −3 ).
We note that we currently do not include any luxury uptake of Fe by Phaeocystis into their gelatinous matrix (Schoemann et al., 2001). Serving as a storage of additional Fe accessible to the Phaeocystis colony when Fe in the seawater gets low, this luxury uptake is thought to relieve it from Fe limitation when Fe concentrations become growth limiting (see discussion in Schoemann et al., 2005). We therefore probably overestimate the Fe limitation of Phaeocystis growth in ROMS-BEC.

145
P. antarctica blooms are typically found where and when waters are comparatively turbulent and the mixed layer is comparatively deep (in comparison to blooms dominated by diatoms, see e.g. Arrigo et al., 1999;Alvain et al., 2008), suggesting that Phaeocystis is better in coping with low light levels than diatoms (e.g. Arrigo et al., 1999). In agreement with laboratory experiments (Tang et al., 2009;Mills et al., 2010;Feng et al., 2010), we therefore choose a higher α PI , i.e. a higher sensitivity of growth to increases of photosynthetically active radiation (PAR) at low PAR levels, for Phaeocystis than for diatoms in 150 ROMS-BEC (see Table 1). Our value (0.63 mmol C m 2 (mg Chl W s) −1 ) corresponds to the average value compiled from available laboratory experiments .
In addition to environmental conditions directly impacting phytoplankton growth rates, loss processes such as grazing, non-grazing mortality, and aggregation impact the simulated biomass levels at any point and time. Grazing on Phaeocystis varies across zooplankton size classes, as a consequence of Phaeocystis life forms spanning several orders of magnitude in 155 size (few µm to cm, Schoemann et al., 2005). Furthermore, Phaeocystis colonies are surrounded by a skin (Hamm et al., 1999), potentially serving as protection from zooplankton grazing. While small copepods have been shown to graze less on Phaeocystis once they form colonies, other larger zooplankton appear to continue grazing on Phaeocystis colonies at unchanged rates (Granéli et al., 1993;Schoemann et al., 2005;Nejstgaard et al., 2007). Based on a size-mismatch assumption of the single grazer in ROMS-BEC and Phaeocystis colonies, we assume a lower maximum grazing rate on Phaeocystis than on diatoms 160 (3.6 d −1 and 3.8 d −1 , respectively, see γ g,max in Table 1). Upon grazing, we assume the fraction of the grazed phytoplankton biomass that is transformed to sinking POC via zooplankton fecal pellet production to be higher for larger and ballasted cells. Consequently, the fraction of grazing routed to POC increases from grazing on SP or diazotrophs to coccolithophores, Phaeocystis, and diatoms (r g in Table 1). Consistent with Nissen et al. (2018), we keep a Holling Type II ingestion functional response here and compute grazing on each prey separately. We refer to Nissen et al. (2018) for a discussion of the relative 165 merits and pitfalls for using Holling Type II versus III.
Non-grazing mortality (such as viral lysis) has been shown to increase under environmental stress for Phaeocystis colonies, causing colony disruption and ultimately cell death (van Boekel et al., 1992;Schoemann et al., 2005). To account for processes causing colony disintegration and for grazing by higher trophic levels not explicitly included in ROMS-BEC, Phaeocystis in ROMS-BEC experience a higher mortality rate than diatoms (0.18 d −1 and 0.12 d −1 , respectively, see γ m,0 in Table 1).

170
Thereby, the chosen non-grazing mortality rate of Phaeocystis assumed in the model is still lower than the estimated rate of viral lysis for Phaeocystis in the North Sea by van Boekel et al. (1992, 0.25 d −1 ), but we note that data on non-grazing mortality of P. antarctica are currently lacking . Furthermore, based on the assumption that for a given biomass concentration, larger cells are more likely than smaller cells to form aggregates and subsequently sink as POC, we use a higher quadratic mortality rate for Phaeocystis (0.005 d −1 ) than for diatoms (0.001 d −1 ) in the model (see γ a,0 in Table 1).
In summary, the spatio-temporal variability of the relative importance of Phaeocystis and diatoms in ROMS-BEC is controlled by the interplay of the environmental conditions and loss processes, which differentially impact the growth and loss rates of these two PFTs and consequently their competitive fitness in the model. In the following, we will describe the model setup and the simulations that were performed to assess the competition between Phaeocystis and diatoms throughout the year in the high-latitude SO. The simulations include a set of sensitivity experiments, with the aim to assess the impact of choices 180 of single parameters or parameterizations on the simulated Phaeocystis biogeography.

Model setup and simulations
With few exceptions, we use the same ROMS-BEC model setup as described in detail in Nissen et al. (2018): At the open northern boundary, we use monthly climatological fields for all tracers (Carton and Giese, 2008;Locarnini et al., 2013;Zweng et al., 2013;Garcia et al., 2014b, a;Lauvset et al., 2016;Yang et al., 2017), and the same data sources are used to initialize 185 the model simulations. At the ocean surface, the model is forced with a 2003-normal year forcing for momentum, heat, and freshwater fluxes (Dee et al., 2011). Satellite-derived climatological total chlorophyll concentrations are used to initialize phytoplankton biomass and to constrain it at the open northern boundary in the model (NASA-OBPG, 2014b), and the fields are extrapolated to depth following Morel and Berthon (1989). Due to the addition of Phaeocystis, the partitioning of total chlorophyll onto the different phytoplankton PFTs is adjusted compared to Nissen et al. (2018): 90% is attributed to small 190 phytoplankton, 4% to diatoms and coccolithophores, respectively, and 1% to diazotrophs and Phaeocystis, respectively. This partitioning is motivated by the phytoplankton community structure at the open northern boundary at 24 • S, where small phytoplankton typically dominate and P. Antarctica are only a minor contributor to phytoplankton biomass (see e.g. Schoemann et al., 2005;Swan et al., 2016). Phaeocystis is initialized with a carbon-to-chlorophyll ratio of 60 mg C (mg chl) −1 (same as small phytoplankton and coccolithophores), whereas diatoms are initialized with a ratio of 36 mg C (mg chl) −1 (Sathyendranath 195 et al., 2009).
We first run a 30 year long physics-only spin-up, followed by a 10 year long spin-up in the coupled ROMS-BEC setup. Our Baseline simulation for this study is then run for an additional 10 years, of which we analyze a daily climatology over the last 5 full seasonal cycles. i.e. from 1 July of year 5 until 30 June of year 10. Apart from having added Phaeocystis and adjusted the parameters of the other PFTs as described in section 2.1, the setup of the Baseline simulation in this study is thereby identical 200 to the Baseline simulation in Nissen et al. (2018). We will evaluate the model's performance with respect to the simulated phytoplankton biogeography in section 3.1 and in the supplementary material.
Furthermore, we perform seven sensitivity experiments, in order to assess the sensitivity of the simulated Phaeocystis biogeography and the competition of Phaeocystis and diatoms to chosen parameters and parameterizations (Table 2). To do so, we set the parameters and parameterizations of Phaeocystis to those used for diatoms in ROMS-BEC (runs 1-6 in Table 2). Generally, 205 the differences in parameters between Phaeocystis and diatoms affect either the simulated growth rates (runs TEMPERATURE, ALPHA PI , and IRON) or loss rates (runs GRAZING, AGGREGATION, and MORTALITY). By successively eradicating the differences between Phaeocystis and diatoms, these simulations allow us to directly quantify the impact of differences in parameters on the simulated relative importance of Phaeocystis for total phytoplankton biomass. To assess the impact of ironlight interactions on the competitive success of Phaeocystis at high SO latitudes, we ultimately run a simulation in which the 210 half-saturation constant of iron (k Fe ) of Phaeocystis is a function of the light intensity, following a polynomial fit of available laboratory data (VARYING_kFE, see Fig. A1b and Garcia et al., 2009). All sensitivity experiments use the same physical spin-up as the Baseline simulation and start from the end of year 10 of the Baseline simulation. Each simulation is then run for an additional 10 years, of which the average over the last 5 full seasonal cycles is analyzed in this study. We compare the simulated spatio-temporal variability in phytoplankton biomass and community structure to available observations of phytoplankton carbon biomass concentrations from the MAREDAT initiative (O'Brien et al., 2013;Leblanc et al., 2012;Vogt et al., 2012), satellite-derived total chlorophyll concentrations (Fanton d'Andon et al., 2009;Maritorena et al., 2010), DMS measurements (Curran and Jones, 2000;Lana et al., 2011), the ecological niches suggested for SO phytoplankton 220 taxa (Brun et al., 2015), and the CHEMTAX climatology based on high performance liquid chromatography (HPLC) pigment data (Swan et al., 2016). The latter provides seasonal estimates of the mixed layer average community composition, which we compare to the seasonally and top 50 m averaged model output of each phytoplankton's contribution to total chlorophyll biomass. The CHEMTAX analysis splits the phytoplankton community into diatoms, nitrogen fixers (such as Trichodesmium), pico-phytoplankton (such as Synechococcus and Prochlorococcus), dinoflagellates, cryptophytes, chlorophytes (all three combined into the single group "Others" here), and haptophytes (such as coccolithophores and Phaeocystis). As noted in Swan et al. (2016), the differentiation between coccolithophores and Phaeocystis in the CHEMTAX analysis is difficult and prone to error.

Analysis framework
Possibly, this is due to the large variability in pigment composition of Phaeocystis as a response to varying environmental conditions, especially regarding light and iron levels (Smith et al., 2010;Wright et al., 2010). Coccolithophores have been reported to only grow very slowly at low temperatures (below ∼8 • C, Buitenhuis et al., 2008), and in the SO, their abundance in the high 230 latitudes south of the polar front is very low (Balch et al., 2016). Therefore, whenever the climatological temperature in the World Ocean Atlas 2013  is below 2 • C at the time and location of the respective HPLC observation, we re-assign data points identified as "Hapto-6" (hence e.g. Emiliania huxleyi) in the CHEMTAX analysis to "Hapto-8" (hence e.g. Phaeocystis Antarctica). Throughout the manuscript, this new category ("Hapto-8 re-assigned") is indicated separately in the respective figures, and leads to a better correspondence of the functional types included in the CHEMTAX-based climatology mixed layer photosynthetically active radiation (MLPAR; W m −2 ). Subsequently, we compare the simulated ecological niche to that observed for example taxa from the SO of each model PFT (such as Phaeocystis Antarctica, Fragilariopsis kerguelensis, Thalassiosira sp., or Emiliania huxleyi, see Brun et al., 2015). While this analysis informs on possible links between the competitive fitness of a PFT and the environmental conditions it lives in, the assessment is hindered due to difficulties in comparing a model PFT to individual phytoplankton species, a sampling bias towards the summer months and the low latitudes, 245 and the neglect of loss processes such as zooplankton grazing to explain biomass distributions. As a consequence, the ecological niche analysis does not allow for the assessment of any temporal variability in PFT biomass concentrations.
In order to assess the simulated seasonality and the seasonal succession of Phaeocystis and diatoms, we identify the bloom peak as the day of peak chlorophyll concentrations throughout the year. Besides the timing of the bloom peak, phytoplankton phenology is typically characterized by metrics such as the day of bloom initiation or the day of bloom end (see e.g. Soppa et al.,250 2016). In this regard, the timing of the bloom start is known to be sensitive to the chosen identification methodology (Thomalla et al., 2015). At high latitudes, the identification of the bloom start based on remotely sensed chlorophyll concentrations is additionally impaired by the large number of missing data in all seasons (even in the summer months, a large part of the SO is sampled by the satellite in less than 5 of the 21 available years, see Fig. S2), complicating any comparison of the high-latitude satellite-derived bloom start with output from models such as ROMS-BEC. To minimize the uncertainty due to the low data 255 coverage in the region of interest for this study, and as the seasonal succession of Phaeocystis and diatoms in the high-latitude SO is mostly inferred from the timing of observed maximum abundances in the literature (e.g. Peloquin and Smith, 2007;Smith et al., 2011), we focus our discussion of the simulated bloom phenology on the timing of the bloom peak. To evaluate the model's performance, we compare the timing of the total chlorophyll bloom peak in the Baseline simulation of ROMS-BEC to the bloom timing derived from climatological daily chlorophyll data from Globcolor (climatology from 1998-2018 based on the daily 25 km chlorophyll product, see Fanton d' Andon et al., 2009;Maritorena et al., 2010).
The release of dimethylsulfoniopropionate (DMSP) via zooplankton grazing, cell lysis, and exudation and the subsequent transformation to DMS by bacterial activity are the major sources of DMS in seawater (e.g. Stefels et al., 2007). With Phaeocystis being the major DMSP producer in the SO (Keller et al., 1989;Liss et al., 1994), the timing of observed peak seawater DMS concentrations (Curran et al., 1998;Curran and Jones, 2000) will allow us to assess the simulated seasonality of Phaeo-265 cystis in the model. Though not explicitly including the biogeochemical cycling of sulphur, we can nevertheless use model output from ROMS-BEC to obtain an estimate of DMS production by Phaeocystis through a simple back-of-the-envelope calculation. To this aim, we use the model-based Phaeocystis biomass loss rates via zooplankton grazing and non-grazing mortality to get the DMSP release from this PFT (integrated annually over the top 10 m; we neglect exudation here), a molar DMSP:C ratio for Phaeocystis of 0.011 (Stefels et al., 2007), and a DMSP-to-DMS conversion efficiency between 0.2-0.7 (the 270 DMS yield depends on the local sulphur demand of bacteria, Stefels et al., 2007;Wang et al., 2015). By comparing the resulting model-based estimates to previously published global estimates of marine DMS emissions (Lana et al., 2011), we obtain an estimate of the importance of SO Phaeocystis for global sulphur cycling.

Assessing phytoplankton competition throughout the year
In ROMS-BEC, phytoplankton biomass P i (mmol C m −3 , i ∈ {P A, D, C, SP, N }) is the balance between growth (µ i ) and loss 275 terms (grazing by zooplankton γ i g , non-grazing mortality γ i m , and aggregation γ i a , see appendix in Nissen et al. (2018) for a full description of the model equations). Here, in order to disentangle the factors controlling the relative importance of Phaeocystis and diatoms for total phytoplankton biomass throughout the year, we use the metrics first introduced by Hashioka et al. (2013) and then applied to assess the competition of diatoms and coccolithophores in ROMS-BEC in Nissen et al. (2018). Same as in Nissen et al. (2018), the relative growth ratio µ ij rel of phytoplankton i and j (e.g. diatoms and Phaeocystis) is defined as the 280 ratios of their specific growth rates (µ i , d −1 ), which in turn depends on environmental dependencies regarding the temperature T , nutrients N , and irradiance I, following: In the above equation, the specific growth rate µ i of each phytoplankton i is calculated as a multiplicative function of a 285 temperature-limited growth rate (f D (T) · µ D max for diatoms and µ PA T for Phaeocystis, see Eq. 1), a nutrient limitation term (g i (N), limitation of each nutrient is calculated using a Michaelis-Menten function, and the most-limiting one is then used here), and a light limitation term (h i (I), Geider et al., 1998). For a detailed description of the underlying functions f (T ), g(N ), and h(I), the reader is referred to the appendix in Nissen et al. (2018). At high-latitudes south of 60 • S, the ratio of the nutrient limitation of growth β N corresponds to that of the iron limitation β Fe in our model (Fig. S1). Consequently, environmental 290 10 https://doi.org/10.5194/bg-2019-488 Preprint. Discussion started: 31 January 2020 c Author(s) 2020. CC BY 4.0 License. conditions regarding temperature, iron, and light decide whether the relative growth ratio is positive or negative at a given location and point in time, i.e., which of the two phytoplankton types has a higher specific growth rate and hence a competitive advantage over the other regarding growth.
Similarly, the relative grazing ratio γ ij g,rel of phytoplankton i and j (e.g. diatoms and Phaeocystis) is defined as the ratio of their specific grazing rates (γ i g , d −1 ) following: In ROMS-BEC, grazing on each phytoplankton i is calculated using a Holling Type II ingestion function (Nissen et al., 2018).
As described in section 2.1, Phaeocystis and diatoms in ROMS-BEC do not only differ in parameters describing the zooplankton grazing pressure they experience, but in parameters describing their non-grazing mortality and aggregation losses as well.
Therefore, in accordance with the relative grazing ratio defined above, we define the relative mortality ratio (γ ij m,rel ) and the rel-300 ative aggregation ratio (γ ij a,rel ) of phytoplankton i and j (e.g. diatoms and Phaeocystis) as the ratio of their specific non-grazing mortality rates (γ i m , d −1 ) and aggregation rates (γ i a , d −1 ), respectively, following: Since the total specific loss rate (γ ij total , d −1 ) of phytoplankton i is the addition of its specific grazing, non-grazing mortality, and aggregation loss rates, the relative total loss ratio γ ij total,rel of phytoplankton i and j (e.g. diatoms and Phaeocystis) is defined as If γ DPA total,rel is positive, the specific total loss rate of Phaeocystis is larger than that of diatoms (and accordingly for the indi-310 vidual loss ratios in Eq. 3-5), and loss processes promote the accumulation of diatom biomass relative to that of Phaeocystis.
While the maximum grazing rate on Phaeocystis is lower than that of diatoms, their non-grazing mortality and aggregation losses are higher (see section 2.1 and Table 1). Ultimately, at any given location and point in time, the interaction between the phytoplankton biomass concentrations (impacting the respective loss rates) and environmental conditions (impacting the respective growth rate) will determine the relative contribution of each phytoplankton type i to total phytoplankton biomass.

Phytoplankton biogeography and community composition in the SO
In ROMS-BEC, total summer chlorophyll is highest close to the Antarctic continent (>10 mg chl m −3 ) and decreases north-320 wards to values <1 mg chl m −3 close to the open northern boundary (Fig. 1a). While this south-north gradient is in broad agreement with remotely sensed chlorophyll concentrations (Fig. 1b), our model generally overestimates high-latitude chlorophyll levels, which has already been noted for the 4-PFT setup of ROMS-BEC (Nissen et al., 2018 Buitenhuis et al. (2013aBuitenhuis et al. ( , 2002Buitenhuis et al. ( -2016 d Monthly output from a biogeochemical inverse model (Schlitzer, 2004) and a data-assimilated model (DeVries and Weber, 2017). and surface chlorophyll are overestimated by a factor 1.8-4.4 and 1.8, respectively (Table 3), and the bias is likely due to a combination of underestimated high-latitude chlorophyll concentrations in satellite-derived products   centrations quickly decline to levels <0.1 mmol C m −3 north of 50 • S (Fig. 1c), a direct result of the restriction of Phaeocystis growth to temperatures < ∼8 • C in the model (Fig. A1a). This is in broad agreement with in situ observations, which suggest highest concentrations (>20 mmol C m −3 ) south of ∼75 • S, and concentrations <5 mmol C m −3 north of ∼65 • S (Fig. 1c & Fig. S3a & b). As a response to the addition of Phaeocystis to ROMS-BEC, the simulated high-latitude diatom biomass concentrations decrease compared of the 4-PFT setup of the model (Nissen et al., 2018). In the 5-PFT setup, the model simu-  (Nissen et al., 2018).
Taken together, the model simulates a phytoplankton community with substantial contributions of coccolithophores and Phaeocystis in the subantarctic and high-latitude SO, respectively (Fig. 2a). CHEMTAX data generally support this latitudinal trend (CHEMTAX is based on high performance liquid chromatography pigment data, see

360
In the 4-PFT setup of ROMS-BEC, the simulated summer phytoplankton community south of 60 • S was often almost solely composed of diatoms ( Fig. S4 and Nissen et al., 2018), suggesting that the implementation of Phaeocystis led to a substantial improvement in the representation of the observed high-latitude community structure (Fig. 2). Concurrently, as the distribution of silicic acid and nitrate is directly impacted by the relative importance of silicifying and non-silicifying phytoplankton, such as Phaeocystis, in the community, the addition of Phaeocystis to the model led to an improvement in the simulated high- and Garcia et al., 2014b), this can certainly be expected to affect the competitive fitness of individual phytoplankton types in 370 the subantarctic and possibly at lower latitudes, which we did not assess further in this study. Overall, our model suggests that Phaeocystis is an important member of the high-latitude phytoplankton community. In the remainder of the manuscript, we will therefore explore the temporal variability in the relative importance of diatoms and Phaeocystis and its implications for SO carbon cycling in more detail.

375
Maximum total chlorophyll concentrations are simulated for the first half of December across latitudes in ROMS-BEC (solid blue line in Fig. 3a). At high SO latitudes south of 60 • S, this is 1-2 months earlier than suggested by satellite estimates (black line in Fig. 3a). Yet, compared to the 4-PFT setup (dashed blue line in Fig. 3a), the simulated timing of peak chlorophyll levels improved in this study, with peak chlorophyll delayed by on average a week in the model upon the implementation of Phaeocystis. The simulated physical biases (i.e., generally too high temperatures and too shallow mixed layer depths, both 380 favoring an earlier onset of the phytoplankton bloom, see Nissen et al., 2018) only partially explain the bias in the simulated timing of maximum chlorophyll levels (see red and green dashed lines in Fig. S7a), suggesting that biological factors must explain the difference between ROMS-BEC and the satellite product. As diatoms dominate the phytoplankton community at peak total chlorophyll concentrations everywhere in the model domain (compare their bloom timing in Fig. 3c to Fig. 3a and to the simulated community composition in Fig. 2b-d), the mismatch in timing is likely related to the representation of this PFT 385 in the model, and is possibly at least partly caused by their comparatively high growth rates at low temperatures (see Fig. A1a).
In contrast to diatoms, maximum chlorophyll concentrations of Phaeocystis are simulated for late November or early December across most latitudes in the model (only around 70 • S a peak in late January is simulated, Fig. 3b). Overall, the timing of simulated peak Phaeocystis chlorophyll levels corresponds well to the suggested timing of observed maximum seawater DMSP concentrations (peak in November/December in Curran et al., 1998;Curran and Jones, 2000) and the delayed maximum at-390 mospheric DMS concentrations (January/February, e.g. Nguyen et al., 1990;Ayers et al., 1991). This further corroborates the The difference in the timing of the bloom peak between the two PFTs is largely <10 days when averaged zonally, but locally 395 exceeds 30 days when looking at individual grid cells in the model (Fig. S7b), in broad agreement with observations, which suggest up to 2 months between the peak chlorophyll concentrations of Phaeocystis and diatoms in the Ross Sea (see e.g. Peloquin and Smith, 2007;Smith et al., 2011). Subsequently, we will assess how environmental conditions and biomass loss processes interact to control the SO phytoplankton biogeography and in particular the competition between Phaeocystis and diatoms at high SO latitudes. ∼20 W m −2 and ∼30 W m −2 for Phaeocystis, diatoms, and coccolithophores, respectively, Fig. 4d-f). Yet, the 5-20 W m −2 higher model-based estimates are consistent with the mixed layer depth bias in ROMS-BEC (see Nissen et al., 2018). While the simulated niche of diatoms regarding temperature and nitrate is in broad agreement with that observed for example SO taxa 425 (Fig. 4b), the simulated niches are wider for Phaeocystis (Fig. 4a) and generally centered around lower temperatures for coccolithophores (Fig. 4c). Acknowledging the difficulties comparing a model PFT to individual phytoplankton taxa, a sampling bias towards temperate and tropical species/strains, and the overall low data coverage in the high-latitude SO, we conclude that the simulated ecological niches of Phaeocystis, diatoms, and coccolithophores are largely in agreement with available observations. Since the analysis of the summer average ecological niches does not inform about sub-seasonal environmental 430 variability and neglects the role of top-down factors in controlling the simulated distributions of phytoplankton types, we will assess these in more detail for the simulated competition between Phaeocystis and diatoms in the following.

High-latitude competition of Phaeocystis and diatoms: An assessment of bottom-up and top-down factors
The temporal evolution of the relative growth ratio, i.e., the ratio of the specific growth rates of diatoms and Phaeocystis (Eq. 2), informs about the competitive advantage of one PFT over the other throughout the year. In ROMS-BEC, the relative 435 growth ratio is negative throughout spring and fall between 60-90 • S (µ PA is on average 5%/6% larger than µ D in spring/fall) and only positive between early December and early February (µ D is on average 5% larger than µ PA in summer, Fig. 5a &   c). Hence, bottom-up factors promote the accumulation of Phaeocystis (diatom) relative to diatom (Phaeocystis) biomass in Calculated from non-log-transformed ratios, i.e., e.g. black bar corresponds to 10 µ DPA rel (see Eq. 2). e) & f) Relative total loss ratio (black) of diatoms vs. Phaeocystis, with contributions of the relative grazing ratio (blue), relative non-grazing loss ratio (red), and relative aggregation ratio (yellow, see Eq. 3-6). g) & h) Seasonally averaged percent difference between diatoms and Phaeocystis in the total specific loss rate (black), specific aggregation rate (yellow), specific grazing rate (blue), and specific mortality rate (red), calculated from non-log-transformed ratios. i) & j) Phaeocystis (blue) and diatom (red) surface carbon biomass concentrations [mmol C m −3 ]. For all metrics, the left panels are averaged over 60-90 • S and those on the right for the Ross Sea.
spring/fall (summer) in this area in our model. In comparison, the relative growth ratio is positive for a shorter time period in the Ross Sea (µ D >µ PA between mid-December and mid-January, Fig. 5b). In fact, averaged seasonally, the specific growth rate 440 of Phaeocystis is larger than that of diatoms for all seasons (µ DPA rel <0), being up to 38% larger in spring (Fig. 5d), suggesting that bottom-up factors are particularly favorable for the accumulation of Phaeocystis biomass in the Ross Sea.
To understand the controls of the spatio-temporal differences in the specific growth rates, the relative growth ratio can be broken down into the environmental factors contributing to the simulated specific growth rate of each phytoplankton type at any point in time (Eq. 2, colors in Fig. 5a-d). As a direct result of the higher k Fe of Phaeocystis in the model (Table 1) limitation of their growth is stronger than that of diatoms, and iron availability, i.e., β Fe in Eq. 2, is an advantage for diatoms both between 60-90 • S and in the Ross Sea at all times (blue area is positive, up to 14% stronger iron limitation of Phaeocystis in both areas in summer, see Fig. 5a-d). Yet, the two subareas differ in the simulated temperature and light limitation of growth of Phaeocystis and diatoms. Overall, temperature is limiting diatom growth more than Phaeocystis growth in both subareas throughout the year (β T <0; up to 16% and 19% stronger growth limitation between 60-90 • S and in the Ross Sea, respectively, see red areas in Fig. 5a-d). While the growth advantage of Phaeocystis regarding temperature is rather constant over time in the Ross Sea (Fig. 5b & d), their advantage is substantially smaller in summer between 60-90 • S (5%, Fig. 5a & c). Additionally, the difference in light limitation between diatoms and Phaeocystis is rather small in this area (3-4% throughout the year, yellow areas in Fig. 5a & c), implying that their differences in α PI (Table 1) are balanced by differences in photoacclimation in ROMS-BEC (Geider et al., 1998;Nissen et al., 2018). Consequently, the combination of the comparatively large advantage in 455 iron limitation and the rather small disadvantage in temperature limitation results in the higher specific growth rates of diatoms in summer between 60-90 • S. In contrast, in the Ross Sea, differences in light limitation between diatoms and Phaeocystis are large, with a maximum advantage for Phaeocystis in spring (their growth is 32% less light limited; Fig. 5b & d). Therefore, the difference in light limitation predominantly controls the seasonality of the relative growth ratio (Fig. 5b) and promotes the dominance of Phaeocystis over diatoms early in the growing season in this area in our model (Fig. 5j), which is not simulated 460 when averaging over 60-90 • S (Fig. 5i). Altogether, in ROMS-BEC, differences in growth between diatoms and Phaeocystis are mostly controlled by seasonal differences in iron/temperature (60-90 • S) and iron/light conditions (Ross Sea), respectively.
Still, given the simulated growth advantage of Phaeocystis throughout much of the growing season in both subareas, bottom-up factors alone cannot explain why Phaeocystis only dominates over diatoms temporarily (Fig. 5i & j), implying that top-down factors need to be considered to explain their biomass evolution in our model.

465
In both subareas, the simulated relative total loss ratio is positive throughout spring and summer, implying that the specific total loss rate of Phaeocystis is higher than that of diatoms (γ PA total >γ D total , see Eq. 6), which favors the accumulation of diatom biomass relative to that of Phaeocystis (Fig. 5e-h). In fact, the total loss rate of Phaeocystis is on average 17%/38% (60-90 • S) and 18%/40% (Ross Sea) higher than that of diatoms in spring/summer in ROMS-BEC ( Fig. 5g & h). In the model, the relative total loss ratio is only negative in early fall in both subareas (γ D total >γ PA total , Fig. 5e & f), but in this season, the difference 470 between diatoms and Phaeocystis in their specific total loss rates is overall smaller than in the other seasons (9% and 3% between 60-90 • S and in the Ross Sea, respectively, Fig. 5g & h). The simulated seasonality of the total loss ratio is the result of the interplay between losses through grazing, aggregation, and non-grazing mortality of each phytoplankton type in ROMS-BEC (Eq. 6, colors in Fig. 5e-h). Of all three loss pathways, differences in aggregation losses are largest between Phaeocystis and diatoms both between 60-90 • S (up to 200% higher aggregation losses for Phaeocystis in summer, yellow in Fig. 5e & g) 475 and in the Ross Sea (up to 250% higher in summer, Fig. 5f & h). In comparison, differences between Phaeocystis and diatoms in grazing (up to 16% and 14% between 60-90 • S and in the Ross Sea, respectively) and mortality losses (50% everywhere) are considerably smaller (see blue and red areas in Fig. 5e-h, respectively), suggesting that aggregation losses predominantly contribute to the simulated differences in the total loss rates between Phaeocystis and diatoms. In summary, between 60-90 • S, the simulated growth advantage of Phaeocystis early in the season (facilitated by advantages 480 in the temperature limitation of their growth) are not large enough to outweigh the disadvantages in iron limitation of their growth and in the biomass losses they experience. As a result, in spring and summer, Phaeocystis do not accumulate substantial biomass relative to (or even dominate over) diatoms in this subarea in ROMS-BEC. In the Ross Sea, however, the simulated growth advantages of Phaeocystis (resulting from advantages in the light and temperature limitation of their growth) are large enough to outweigh the disadvantages in iron limitation and specific biomass loss rates, allowing them to dominate over 485 diatoms early in the growing season in our model and explaining the simulated succession from Phaeocystis to diatoms close to the Antarctic continent (see also section 3.2). Ultimately, this simulated spatio-temporal variability in the relative importance of Phaeocystis and diatoms has implications for both SO and global carbon cycling, which we will assess subsequently.

Quantifying the importance of Phaeocystis for the cycling of carbon
Phaeocystis is an important member of the SO phytoplankton community in our model, particularly south of 60 • S, where it 490 contributes 46% and 40% to total annual NPP and POC formation, respectively. (Table 3 & Fig. 6). Even when considering the entire region south of 30 • S, the contribution of Phaeocystis to NPP (15%) and POC production (16%) is sizeable. The simulated spatial differences in phytoplankton community structure have direct implications for the fate of organic carbon upon biomass loss, and Fig. 6 illustrates the different pathways of POC formation for each PFT in ROMS-BEC. Overall, in our model, the p ratio, i.e., the fraction of NPP that is transformed to sinking POC (Laufkötter et al., 2016), is higher at high latitudes south 495 of 60 • S (45%) than the domain average (37%, Fig. 6). This is a direct result of the higher fraction of large phytoplankton  types, i.e., Phaeocystis and diatoms, in the phytoplankton community in this area (Fig. 2), facilitating more carbon export relative to NPP in the model. In fact, our model results suggest that these two large phytoplankton types contribute more to POC formation than to total biomass (compare yellow and green boxes in Fig. 6). Integrated annually, diatoms contribute most of all PFTs to POC formation in our model (60% and 49% between 30-90 • S and 60-90 • S, respectively, Fig. 6). For both 500 diatoms and Phaeocystis, grazing by zooplankton (i.e., the formation of fecal pellets) is the most important pathway of POC production in ROMS-BEC (black arrows in Fig. 6, 13%/52% and 29%/37% of total POC production for Phaeocystis/diatoms between 30-90 • S and 60-90 • S, respectively). Yet, at high latitudes (60-90 • S), aggregation of Phaeocystis biomass contributes significantly to POC formation (20% of total POC production, 9% of NPP, grey arrows in Fig. 6b). Given that the loss of biomass via a given pathway is a function of the local biomass concentrations of each PFT at any given point in time (see January between 30-90 • S (blue bars in Fig. 7e) and even dominantly contributes to POC formation between November and January at high SO latitudes (up to 65%, blue bars in Fig. 7f). Altogether, this implies that both spatial and temporal variations in SO phytoplankton community structure critically impact the fate of carbon beyond the upper ocean.

Biogeochemical implications of SO Phaeocystis biogeography
Our simulated contribution of Phaeocystis to SO NPP and net POC formation is higher than that found in the previous modeling study by Wang and Moore (2011), particularly at higher latitudes: They diagnosed a contribution of 13% and 23% to NPP south of 40 • S and 60 • S, respectively, compared to our estimates of 15% for south of 30 • S and 46% for the region south 525 of 60 • S. We interpret the difference to stem primarily from differences in parameter choices of the PFTs between the two models. As an example, the half-saturation constant of iron of Phaeocystis is 125% larger than that of diatoms in the model by Wang and Moore (2011, 0.18 µmol m −3 for Phaeocystis as compared to 0.08 µmol m −3 for diatoms), but only 33% larger in ours (Table 1). The resulting smaller difference in the growth limitation by iron in ROMS-BEC leads to more competitive Phaeocystis relative to diatoms in the iron-depleted SO and can at least partially explain the higher relative importance of 530 Phaeocystis in our model. In fact, differences in model parameters between Phaeocystis and diatoms in ROMS-BEC can alter the simulated contribution of Phaeocystis to total NPP from 5-32% and 17-63% between 30-90 • S and 60-90 • S, respectively (see section A1, as well as section 2.2 and Table 2). This illustrates how single model parameters sensitively impact the competitive success of Phaeocystis in the SO. Still, the simulated community structure in the Baseline simulation with ROMS-BEC is supported by available observations (see section 3.1), giving us confidence in our estimates.

535
The simulated contribution of Phaeocystis to POC export in ROMS-BEC (16% and 40% south of 30 • S and 60 • S) is in broad agreement with the previous estimate from Wang and Moore (2011, 19% and 30% south of 40 • S and 60 • S, respectively). This is despite the differences in phytoplankton community structure between the two models (see above) and demonstrates our on-going limited quantitative understanding of the fate of biomass losses (see also Laufkötter et al., 2016). Across the parameter sensitivity runs in ROMS-BEC (section 2.2), the contribution of Phaeocystis to POC production and export varies from 4-540 23% and 13-59% south of 30 • S and 60 • S, respectively. In addition to this uncertainty resulting only from the growth and loss parameters of Phaeocystis in the model (Table 2), further uncertainty arises from parameters describing the partitioning of biomass losses amongst dissolved and particulate carbon species, which we did not assess in this study. Acknowledging that the exact numbers are highly sensitive to parameter choices in the model, our analysis reveals how the pathways of POC production, in particular the relative importance of fecal pellets from zooplankton and aggregated phytoplankton cells, are 545 impacted by the simulated spatio-temporal variability in phytoplankton community structure throughout the year (Fig. 7). In this regard, the simulated strong temporal coupling between POC fluxes and biomass distributions in ROMS-BEC is a direct result of the model formulations describing particle sinking (Lima et al., 2014). This coupling is supported by observations, e.g., from the Ross Sea, where the POC flux from the upper ocean has been found to be closely linked to biomass levels in the overlying surface layer (with aggregates being an important vector for POC export when Phaeocystis dominated the 550 community, Asper and Smith, 1999). Yet, the coupling in our model is potentially too strong in other areas, where reprocessing of POC by zooplankton in the upper ocean or lateral advection of POC could decouple the seasonal evolution of phytoplankton biomass and POC export (e.g. Lam and Bishop, 2007;Stange et al., 2017), the effect of which we can currently not assess.
Given the possibly large importance of different POC production pathways for carbon and nutrient cycling through their impact on the remineralization depth of organic matter, these processes should be better constrained in the future, in order to further 555 quantify the imprint of spatio-temporal variations in the relative importance of Phaeocystis for the high-latitude cycling of carbon. Despite these uncertainties, Phaeocystis is clearly a substantial player on the global scale. Comparing the modelsimulated integrated NPP and POC export across the entire model domain with data-based estimates of global NPP and POC export suggests that SO Phaeocystis alone contribute about 5% to globally integrated NPP (58±7 Pg C yr −1 , Buitenhuis et al., 2013a), and about the same percentage to global POC export (9.1±0.2 Pg C yr −1 , DeVries and Weber, 2017).

560
Besides its impact on the carbon cycle, Phaeocystis is the major contributor to the marine sulphur cycle in the SO through its production of DMSP (Keller et al., 1989;Liss et al., 1994;Stefels et al., 2007). Integrating the modeled Phaeocystis biomass loss rates via zooplankton grazing and non-grazing mortality and making assumptions regarding the corresponding DMS(P) release (see scetion 2.3.1), our estimated annual DMS production by Phaeocystis in ROMS-BEC amounts to 3.3-11.5 Tg S and 1.8-6.4 Tg S south of 30 • S and 60 • S, respectively. Consequently, assuming that all of this DMS production quickly escapes 565 to the atmosphere, our estimates correspond to 11.6-40.1% (30-90 • S) and 6.5-22.7% (60-90 • S) of the global flux of DMS to the atmosphere previously estimated by Lana et al. (2011, 28.1 Tg S yr −1 ). Our estimate is an upper bound, however, as not all DMS produced in seawater is readily released to the atmosphere. In fact, a fraction is likely broken down by bacteria, by photolysis, or is mixed down in the water column (see e.g. Simó and Pedrós-Alló, 1999;Stefels et al., 2007). Still, given that other phytoplankton types also produce DMS(P) (Keller et al., 1989;Stefels et al., 2007), the ROMS-BEC-based contribution of SO Phaeocystis alone (3.3-11.5 Tg S yr −1 ) to the global flux of DMS to the atmosphere is in agreement with the flux suggested in Lana et al. (2011, 8.1 Tg S yr −1 south of 30 • S, i.e., 29% of their global estimate), and the substantial contribution of SO Phaeocystis underpins its major role for the global cycling of sulphur.

Drivers of SO Phaeocystis biogeography and their competition with diatoms
In ROMS-BEC, the interplay of iron availability with temperature (60-90 • S) and light levels (Ross Sea), respectively, largely controls the competitive fitness of Phaeocystis relative to diatoms in the high-latitude SO. Yet, differences in the simulated biomass loss rates between the two PFTs (in particular via aggregation) need to be considered in order to explain why peak Phaeocystis biomass levels precede those of diatoms only close to the Antarctic continent in the model. In the literature, the spatial distribution of Phaeocystis and diatoms and the temporal succession from Phaeocystis to diatoms is almost exclusively discussed in terms of light and iron availability (see e.g. Arrigo et al., 1999;Smith et al., 2014). In this context, regions/times of 580 low light and/or high mixed layer depth are typically associated with high Phaeocystis abundance (Alvain et al., 2008;Smith et al., 2014), explaining their bloom in spring, whereas iron availability has been suggested to largely control the magnitude of the summer diatom bloom (Peloquin and Smith, 2007;Smith et al., 2011). This is in agreement with the simulated dynamics and parameters chosen in ROMS-BEC, in which the difference in light limitation between growth of Phaeocystis and diatoms facilitates early Phaeocystis blooms in the Ross Sea. Yet, it has to be noted that advantages in temperature limitation contribute 585 to the growth advantage of Phaeocystis in the high-latitude SO in ROMS-BEC as well. Currently, this growth advantage of Phaeocystis at temperatures <4 • C is possibly underestimated in the model, as diatom growth at low temperatures is currently overestimated when comparing to available laboratory measurements (Fig. A1a). Nevertheless, in agreement with Peloquin and Smith (2007) and Smith et al. (2011), when diatoms reach peak chlorophyll levels in summer in our model, the simulated difference in iron limitation between the two PFTs is largest across the high-latitude SO ( Fig. 5a & b), suggesting that any 590 change in summer iron availability will indeed strongly impact peak diatom and hence total chlorophyll levels in ROMS-BEC.
Still, an important limitation in the assessment of the role of iron in controlling the relative importance of Phaeocystis in the high-latitude phytoplankton community is the assumption of a constant k Fe of Phaeocystis in the model (0.2 µmol m −3 , Table   1). In laboratory experiments, the affinity of Phaeocystis for iron has been shown to be sensitive to light (Garcia et al., 2009), which is not accounted for in the Baseline simulation of ROMS-BEC. In order to assess the possible effect of a varying k Fe on 595 the competition between Phaeocystis and diatoms, we fit a polynomial function to describe the k Fe of Phaeocystis as a function of the light level (VARYING_kFE simulation in Table 2, Fig. A1b, Garcia et al., 2009). Acknowledging the substantial uncertainty in the fit, our model simulates k Fe <0.2 µmol m −3 and even k Fe <0.15 µmol m −3 (corresponding to the k Fe of diatoms in the Baseline simulation) only at highest light intensities in summer and mostly close to the surface, and 0.2 µmol m −3 < k Fe ≤ 0.26 µmol m −3 elsewhere as a result of low light levels ( Fig. S9a & b). While the contribution of Phaeocystis to total 600 NPP is only affected to a lesser extent as a consequence (37% and 13% south of 60 • S and 30 • S, respectively, instead of 46% and 15% in the Baseline simulation), the simulated phytoplankton seasonality is impacted substantially. The maximum chlorophyll levels of diatoms occur earlier than those of Phaeocystis in many more places of the SO compared to the Baseline simulation, both in coastal areas and in the open ocean ( Fig. S9c & d). Thus, in order to include light-iron interactions in future modeling efforts with Phaeocystis and to assess their impact on the competition of Phaeocystis with diatoms throughout 605 the SO, additional measurements are needed for how k Fe varies as a function of the surrounding light level. Taken together, given the likely lower k Fe of Phaeocystis at high-latitude light levels in summer than what we currently assume in ROMS-BEC ( Fig. A1b and Garcia et al., 2009) and given the likely underestimation of their growth advantage in temperature, we probably currently underestimate the competitive advantage in growth of Phaeocystis relative to diatoms in the model. However, such a potential underestimation in growth advantage does not automatically mean that the contribution of Phaeocystis to the phyto-610 plankton community is underestimated as well. This is because of the important role of biomass loss processes to explain why Phaeocystis do not outcompete diatoms everywhere in the high latitudes in ROMS-BEC (Fig. 5). Furthermore, the simulated phytoplankton community structure is in good agreement with available observations (Fig. 2).
Loss processes, such as aggregation and grazing, clearly matter for the competitive advantage of one PFT over another, but these loss processes are generally not well quantified and often not studied with sufficient detail. For example, while the 615 modeling study by Le Quéré et al. (2016) demonstrates the importance of such top-down control for total SO phytoplankton biomass concentrations, an analysis of the impact on phytoplankton community structure is lacking. In fact, in the literature, only few studies discuss the role of top-down factors for the relative importance of Phaeocystis and diatoms in the high-latitude SO (Granéli et al., 1993;van Hilst and Smith, 2002). In agreement with our results, aggregation has been suggested to be an important process facilitating high POC export when Phaeocystis biomass is high (Asper and Smith, 1999), but to what extent 620 this process significantly contributes to the observed relative importance of Phaeocystis and diatoms throughout the year in the high-latitude SO remains largely unknown.
Here, our findings suggest an important role for biomass loss processes in controlling the relative importance of Phaeocystis and diatoms in ROMS-BEC, but very little quantitative information exists to constrain model parameters (see section 2.1) or to validate the simulated aggregation, grazing, or non-grazing mortality loss rates of Phaeocystis and diatoms over time.

625
Additionally, as discussed in Nissen et al. (2018), the single zooplankton grazer is a major limitation of ROMS-BEC. To what extent accounting implicitly for grazing by higher tropic levels in the non-grazing mortality term makes up for not including more zooplankton PFTs remains unclear. Nevertheless, by changing the overall coupling between phytoplankton and zooplankton and through the distinct grazing preferences of the different zooplankton types, the addition of larger zooplankton grazers would likely change the simulated temporal evolution of Phaeocystis and diatom biomass in the model (Le Quéré et al.,630 2016). Therefore, the above mentioned uncertainties should be addressed by future in situ or laboratory measurements in order to better constrain the simulated biomass loss processes, as our findings suggest these to be necessary to explain the seasonal evolution of the relative importance of Phaeocystis and diatoms in the high-latitude SO.

Limitations & Caveats
Our results may be affected by several shortcomings regarding the parameterization of Phaeocystis, in particular the represen-635 tation of its life cycle, the fate of its biomass losses, and its nutrient uptake stoichiometry. We considered here only colonial Phaeocystis, thereby implicitly assuming that a seed population of solitary cells is always available for colony formation. Not including an explicit parameterization for single cells and hence life cycle transitions might substantially impact both the seasonal Phaeocystis biomass evolution and the competition with diatoms, as solitary cells have been proposed to require less iron (Veldhuis et al., 1991) and are possibly subject to higher loss rates due to e.g. zooplankton grazing compared to colonies 640 (Smith et al., 2003;Nejstgaard et al., 2007). The transition from solitary to colonial cells is a function of the seed population and light and nutrient levels (Verity, 2000), and transition models have been applied in SO marine ecosystem models (e.g. Popova et al., 2007;Kaufman et al., 2017). Yet, implementing morphotype transitions of Phaeocystis into a basin-wide SO model such as ROMS-BEC is severely hindered by data availability. At the moment, 390 Phaeocystis biomass observations are included in the MAREDAT data base south of 30 • S, and the distinction between solitary and colonial cells is often dif-645 ficult (Vogt et al., 2012), impeding the basin-wide model evaluation of both Phaeocystis life stages. In addition, colonies of Phaeocystis are surrounded by a gelatinous matrix, which contains nutrients and carbon , leading to an underestimation of modeled Phaeocystis carbon biomass estimates if not accounting for this mucus (Vogt et al., 2012). In ROMS-BEC, this underestimation is likely small, as <20% of the total Phaeocystis biomass is reportedly encorporated into the mucus in the SO ( Fig. 9 in Vogt et al., 2012). Nevertheless, through its function as a nutrient storage, the mucus promotes the accumulation of Phaeocystis biomass relative to other phytoplankton types when the latter become limited by low nutrient availability. While the gelatinous matrix is additionally thought to prevent grazing, the literature on grazing losses of Phaeocystis colonies is non-conclusive . This is possibly a result of the large range of sizes of both Phaeocystis and the respective grazers, with smaller zooplankton typically grazing less on Phaeocystis colonies than larger zooplankton (see reviews by Schoemann et al., 2005;Nejstgaard et al., 2007). As discussed above, the fate of biomass losses of Phaeocystis 655 is still poorly constrained (this applies to all model PFTs, see also Laufkötter et al., 2016). Currently, ROMS-BEC treats POC from all formation pathways equal, i.e., once produced, there is no differentiation between POC originating from diatoms or Phaeocystis or from grazing or aggregation. In reality, Phaeocystis aggregates might be recycled more readily than those from diatoms. This could reconcile our model results, i.e., the substantial simulated contribution of Phaeocystis to POC export at 100 m, with observations which suggest that the contribution of Phaeocystis to the POC flux across 200 m is small (<5%, 660 Gowing et al., 2001;Accornero et al., 2003;Reigstad and Wassmann, 2007). Ultimately, the C:P and N:P nutrient uptake ratios by Phaeocystis and diatoms are higher (147±26.7 and 19.2±0.61) and lower (94.3±20.1 and 9.67±0.33), respectively (Arrigo et al., 1999(Arrigo et al., , 2000, than those originally suggested by Redfield and currently used in ROMS-BEC (117:16:1 for C:N:P uptake by Phaeocystis and diatoms, Anderson and Sarmiento, 1994). Consequently, this suggests that not accounting for the non-Redfield ratios in nutrient uptake by these PFTs leads to an over(under)estimation of carbon fixation per unit of P and 665 hence POC export where/when Phaeocystis (diatoms) dominate the phytoplankton community.

Conclusions
In this modeling study, we present a thorough assessment of the factors controlling the relative importance of SO Phaeocystis and diatoms throughout the year and quantify the implications of the spatio-temporal variability in phytoplankton community structure for POC export. In ROMS-BEC, Phaeocystis colonies are an important member of the SO phytoplankton community, 670 contributing 15% (16%) to total annual NPP (POC export) south of 30 • S. Moreover, their contribution is threefold higher south of 60 • S in our model. Given that our results imply a contribution of approximately 5% of SO Phaeocystis colonies to total global NPP and POC export, respectively, we recommend the inclusion of an explicit representation of Phaeocystis in ecosystem models of the SO. This will allow for a more realistic representation of the SO phytoplankton community structure, in particular the relative importance of silicifying diatoms and non-silicifying phytoplankton, which we here find to significantly impact the simulated high-latitude nutrient distributions. A follow-up study with both regional SO and global marine ecosystem models should more closely assess what the impact of this simulated change in the relative concentrations of silicic acid and nitrate in the high-latitude SO is on subantarctic and low latitude phytoplankton dynamics.
On a basin-scale, we find that the competition of Phaeocystis and diatoms is controlled by seasonal differences in temperature and iron availability, but that variations in light levels are critical on a local scale, e.g. facilitating early blooms by Phaeocystis 680 in the Ross Sea, in agreement with previous studies. Yet, our model suggests that the relative importance of Phaeocystis and diatoms over a complete annual cycle is ultimately determined by differences in their biomass loss rates (such as zooplankton grazing and aggregation), which in turn impacts the formation of sinking particles and hence carbon transfer to depth. Despite knowing of the importance of top-down factors for global phytoplankton biomass distributions (Behrenfeld, 2014) and for the formation of sinking particles (e.g. Steinberg and Landry, 2017), model parameters describing the fate of carbon after its 685 fixation during photosynthesis are still surprisingly uncertain (Laufkötter et al., 2016), complicating the assessment of the role of biomass loss processes in regulating global biogeochemical cycles.
Environmental conditions in the SO have changed considerably in the last million years (see e.g. Martínez-García et al., 2014), as well as during the past decades (Constable et al., 2014), and are projected to change further during this century (IPCC, 2014). These changes will impact the competitive fitness of Phaeocystis and diatoms (see e.g. Hancock et al., 2018;690 Boyd, 2019) and hence affect the entire phytoplankton community with likely repercussions for the entire food web (Smetacek et al., 2004). Consequently, based on our results, future laboratory and modeling studies should assess how uncertainties in marine ecosystem models surrounding e.g. the parameterization of the life cycle of Phaeocystis and the fate of biomass losses impact the simulated relative importance of this phytoplankton type and carbon transfer to depth at high SO latitudes. Thereby, such studies will allow us to better constrain how potential future changes in the high-latitude phytoplankton community 695 structure impact global biogeochemical cycles. We assess the sensitivity of the simulated annual mean Phaeocystis biogeography to parameter choices by performing a set of sensitivity experiments (runs 1-6 in Table 2). Overall, the simulated surface Phaeocystis biomass concentrations change by ±50% for each of the experiments in the high-latitude SO (Fig. A2). Between 60-90 • S and in the Ross Sea, the largest increase in Phaeocystis biomass concentrations is simulated for AGGREGATION (+112% and +96%, respectively, Fig. A2b & c), whereas the strongest decline is simulated for ALPHA PI (-76% and -87%, respectively, Fig. A2b & c). As a response to changes 705 Figure A1. a) Growth rates of Phaeocystis antarctica colonies as a function of temperature (conditions of nutrients and light are non-limiting) in laboratory data (grey triangles, see compilation by Schoemann et al., 2005) and as used in ROMS-BEC (black line, see Eq. 1). Green circles and the green line show the temperature-limited growth rate of diatoms in laboratory data (see compilation by Le Quéré et al., 2016) and as used in ROMS-BEC, respectively (see also Table 1). b) Half-saturation constant of Fe (kFe) of Phaeocystis as a function of light intensity I (W m −2 ) in laboratory data (red circles) and the polynomial fit (k PA Fe (I) = 2.776 · 10 −5 · (I + 20) 2 -0.00683 · (I + 20) + 0.46) without (black) and with (dashed grey, as used in ROMS-BEC in simulation VARYING_kFe, see Table 2) the correction at low and high light intensities to restrict kFe to the range measured in the laboratory experiments by Garcia et al. (2009). The green lines correspond to the half-saturation constants used for Phaeocystis (solid), diatoms (dashed), and coccolithophores (dotted) in the Baseline simulation in this study (see Table 1).
In both panels, the blue lines correspond to the simulated annual range in a) sea surface temperature [ • C] and b) light intensity [W m −2 ] between 50-60 • S, 60-70 • S, and 70-80 • S, respectively.
in Phaeocystis parameters, diatom biomass changes overall more than that of SP on a basin scale, suggesting Phaeocystis is indeed mostly competing with diatoms for resources in the high-latitude SO. Between 60-90 • S, the magnitude of change is similar for the experiments TEMPERATURE (-73%), ALPHA PI (-76%), and IRON (+70%), while in the Ross Sea, the response in IRON is substantially smaller (+17%) than that for the other two experiments (-82% and -87% for TEMPERATURE and ALPHA PI , respectively; Fig. A2b & c). This supports our findings from section 3.4, namely that the difference in light 710 sensitivity between Phaeocystis and diatoms is more important in coastal areas than on a basin scale in controlling the relative importance of Phaeocystis for total phytoplankton biomass.
Author contributions. MV and CN conceived the study. CN set up the model simulations, performed the analysis, and wrote the paper. MV contributed to the interpretation of the results and the writing of the paper. Competing interests. The authors declare that they have no conflict of interest.

715
Acknowledgements. We acknowledge all the scientists who contributed phytoplankton and zooplankton cell count data to the MAREDAT initiative and William Balch, Helen Smith, Mariem Saavedra-Pellitero, Gustaaf Hallegraeff, José-Abel Flores, and Alex Poulton for providing additional cell count data. Furthermore, GlobColour data (http://globcolour.info) used in this study has been developed, validated, and distributed by ACRI-ST, France. We would like to thank Nicolas Gruber, Matthias Münnich, and Domitille Louchard for valuable discussions and Damian Loher for technical support. Additionally, we would like to thank Gianna Ferrari for the analysis of early ROMS-BEC simulations 720 with Phaeocystis. This research was financially supported by the Swiss Federal Institute of Technology Zürich (ETH Zürich) and the Swiss National Science Foundation (project SOGate, grant no. 200021_153452). The simulations were performed at the HPC cluster of ETH Zürich, Euler, which is located in the Swiss Supercomputing Center (CSCS) in Lugano and operated by ETH ITS Scientific IT Services in Zürich. Model output is available upon request to the corresponding author, Cara Nissen (cara.nissen@usys.ethz.ch).