Southern Ocean BGC-Argo Detect Under Ice Phytoplankton Growth Before Sea Ice Retreat

. The seasonality of sea ice in the Southern Ocean has profound effects on the life cycle (phenology) of phytoplankton residing under the ice. The current literature investigating this relationship is primarily based on remote sensing, which often lacks data for half the year or more. One prominent hypothesis holds that following ice retreat in spring, buoyant melt waters enhance irradiance levels, triggering a bloom which follows the ice edge. However, an analysis of BGC-Argo data sampling under Antarctic sea ice suggests that this is not necessarily the case. Rather than precipitating rapid accumulation, we show that 5 melt waters enhance growth in an already highly active phytoplankton population. Blooms observed in the wake of the receding ice edge can then be understood as the emergence of a growth process that started earlier under sea ice. Indeed, we estimate that growth initiation occurs, on average, 4-5 weeks before ice retreat, typically starting in August and September. Novel techniques using on-board data to detect the timing of ice melt were used. Furthermore, such growth is shown to occur under conditions of substantial ice cover (>90% satellite ice concentration) and deep mixed layers (>100 m), conditions previously thought 10 to be inimical to growth. This led to the development of several 0D model experiments in which we sought to investigate the mechanisms responsible for such early growth. The results of theses experiments suggest that a combination of higher light transfer (penetration) through sea ice and extreme low light adaptation by phytoplankton can account for the observed phenology.

. Distribution of great-circle distances of under ice profiles to the estimated satellite sea ice edge (latitude of 15% sea ice concentration contour). Negative values indicate that the profile is poleward of the ice edge. salinity closest to the surface, is computed and compared to the temperature measured at the same depth. It is important to note that the depth of these near-surface measurements will vary from ∼20 -25 m in winter, to ∼0 -5 m in summer. This is because of the on-board ice avoidance algorithm described above: in winter the temperature threshold is generally exceeded, and so sampling ceases ∼20 m from the surface, while in other months this condition is generally not met and so floats are able to 90 sample much closer to the surface. Therefore, since the on-board ice avoidance algorithm is intentionally conservative, it may assume there is ice present when in fact melting has already occurred.
While the fact that winter profiles generally only sample up to 20 m may seem unimportant, it actually has significant bearing on the under ice detection algorithm used in this study. This is because within the upper 20 m in winter, water is generally above its freezing point (if salinity is taken into account). Therefore, in order to delineate under ice from open ocean profiles, 95 one has to assume some degree of cooling from the last measurement in the profile to the surface. Since the realized degree of cooling over this winter surface layer cannot be observed, we tested several values (the corresponding effect they had on a key result of the paper is shown in Figure 5). The orange and green curves in this Figure depict the change in the probability density function when increasing and decreasing the cooling threshold by 20%, respectively (the details of what is depicted in the figure is discussed below in the "Growth Initiation" section). The blue curve and associated histogram depicts the chosen 100 4 https://doi.org/10.5194/bg-2020-257 Preprint. Discussion started: 17 July 2020 c Author(s) 2020. CC BY 4.0 License. value of 0.1°C used in this study (i.e. we assume a decrease of 0.1°C from ∼20 m to the surface). We would note that the essential features of the distribution remain unchanged in this sensitivity test.
In addition to the above testing, two further checks were performed to assess the validity of using an assumed rate of cooling to detect under ice profiles. The first approach is shown in Figure 1, which plots the distribution of distances of under ice profiles to the satellite ice edge. Here the ice edge is defined by the 15% sea ice concentration contour, following previous 105 satellite-based studies (e.g. Stroeve et al. (2016)) We found that the vast majority of profiles where located 100 km or more south of the ice edge, with 13.6% being north of the edge. It is important to point out here that while sampling under ice, floats do not communicate their location, since they are prevented from surfacing. A simple linear interpolation is used to estimate the location of the under ice profiles (based on the relative time stamp difference), with an approximate maximum error of 100 km as reported by Riser et al. (2018). It is precisely because of this uncertainty that we chose to use on-board data to detect 110 under ice profiles (as well as to detect melting), as opposed to flagging profiles as under ice based on their relative position to the satellite ice edge. The distribution shown in Figure 1 is included to illustrate that there is broad agreement between the two methods, although the use of on-board data should be more accurate given the uncertainty of the float location.
The second approach used to assess the under ice detection method involved visual inspection of time series of mean mixed layer temperature and salinity like those shown in Figure 2A. This consisted of comparing the timing of the transition from 115 under ice to open ocean (depicted by the black vertical line in Figure 2A), with the associated changes in temperature and salinity. Both raw (dashed) and filtered (solid) time series are shown, with a first-order, low-pass Butterworth digital filter employed with a cut-off frequency of 0.1 Hz. By inspecting a subset of floats sampling under ice, we found good visual agreement between the computed timing of the transition (black vertical line) and the corresponding tendency of the curves toward freshening and warming of the surface ocean. The best agreement is achieved by assuming a relative temperature 120 difference (between the last winter measurement and the surface) 0.1°C as described above, which is why this value was chosen over other candidates. A sample of time series for floats other than that shown in Figure 2A is provided in supplementary Figures S1-S4.

Melt Detection
Once a transition from ice cover to open ocean has been established, our algorithm then verifies that these changes are associ-125 ated with melting. This is done by computing time derivatives of surface temperature and salinity at the time of transition (data are taken from measurements closest to the surface). In order to be classified as a melt event, the temperature derivative must be positive (i.e increasing temperature) with a negative salinity derivative that is persistent for 1 month following transition.
An example of a such a melt event is shown in Figure 2A, where salinity (blue lines) decreases gradually for ∼1 month prior to transition, while temperature (red lines) begins to steadily increase after remaining consistently below freezing. At least 130 three consecutive under ice profiles (equivalent to ∼1 month since profiles are at 10-day frequency) are needed to detect a melt event. In cases where multiple transitions occur in one season, the transition with the strongest signal (i.e steepest time derivative) of warming and freshening is chosen. This enables us to filter out transitions which occur as a result of advection or high frequency warming associated with synoptic variability. Apart from the three criteria discussed above (transition from under ice to open ocean, positive temperature derivative, neg-135 ative salinity derivative), additional inspection of time series of stratification depth (our chosen metric for assessing vertical mixing, termed N d ) was performed. This depth is defined as the point at which the Brunt-Väisälä frequency reaches its maximum value in the upper water column, implying a region of maximum resistance to mixing (Gill, 1982). Furthermore, this measure of the depth of the mixed layer has been shown to be more ecologically relevant in the Southern Ocean than other more traditional methods involving density/temperature thresholds (Carvalho et al., 2017). As is discussed in Section 1, the release of melt waters tends to stratify the surface ocean, and so N d should rapidly decrease following the detected melt event. In Figure 2B  Our main metric for assessing the relationship between melting and phenology is termed growth initiation (GI). It is defined here as the point at which the time derivative of mean mixed layer chl-a exceeds the median time derivative computed for the growth period in question. These time derivatives, here taken as a proxy for growth rates, are only computed over the period of positive growth. This period is determined from a filtered time series of mean mixed layer chl-a used to remove variability 150 at the 10-day sampling frequency (the actual value of the median is computed from the raw signal). A first-order, low-pass Butterworth digital filter is employed with a cut-off frequency of 0.1 Hz. An example of the resulting filtered time series is shown in Figure 2B and compared to the original raw signal. Also shown in the figure by the black vertical line is the timing of GI for this particular season. The distance between the 2 black vertical lines in panels A and B in Figure 2 then denotes the timing difference between melting and GI as shown in Figure 5 for all float data. growth rates (as opposed to an absolute threshold value) to be more appropriate, since it avoids any biases in the median which may be created by long periods of close to zero chl-a concentration under ice (followed by a rapid increase).

Model Set-up
A biogeochemical box model is employed in this study to investigate the drivers of under ice growth. The model is based 165 on the Biogeochemical Flux Model (BFM) framework, for which documentation can be found in Vichi et al. (2015). Our particular configuration is a "0.5D" box model where all the major components of the marine biogeochemical system are simulated, namely, phytoplankton, zooplankton, organic and inorganic matter, nutrients and bacterioplankton. The model is termed "0.5D" due to the fact that the depth of the box is able to vary to simulate the effect of vertical mixing. In this case the only effect taken into account is the attenuation of light with increasing mixed layer depth.

170
Since our study is process-oriented, we chose to simplify the model as much as possible, while still retaining the major features of interest. Accordingly, only 1 phytoplankton (diatoms) and 2 zooplankton groups (omnivorous mesozooplankton and heterotrophic nanoflagellates) are simulated. In terms of nutrients, all those standard in biogeochemical models are included, as well as silicate and iron. Initial nutrient conditions where chosen to be representative of the Southern Ocean south of ∼60°S, with non-limiting concentrations of nitrate (31.8 mmol/m 3 ), phosphate (2 mmol/m 3 ) and silicate (40 mmol/m 3 ). An initial 175 dissolved iron concentration of 0.3 µmol/m 3 is applied to all experiments, which gave the most realistic magnitude of summer growth when compared to float data.
The model is forced daily with solar radiation, satellite sea ice concentration, float temperature and salinity, as well as mixed layer depth derived from float data (refer to Section 2.1 for data sources). Light available at the surface is scaled by the sea ice concentration by simply multiplying the incident radiation by the percentage of open ocean derived from remote sensing data.

Experiment Design
Three core experiments were conducted in 4 study regions, with each run having a spin-up time of 10 years to allow for adjustment to a repeating annual cycle (although in most cases adjustment took only a few years). In Table 1 we provide an overview of the available float data in each study region. For every complete time series of float observations we performed the set of three core experiments. First, 2 sets of experiments were run to test the effect of sea ice cover on phytoplankton 185 phenology: a run with no ice (OPEN) and a run with imposed satellite sea ice concentration (ICE). A third experiment sought to test the combined effect sea ice cover and increased low light adaptation by phytoplankton had on phenology (LLA). This was achieved by increasing the initial slope of the Photosynthesis-irradiance curve by a factor of 10, thus enhancing photosynthetic efficiency at light levels close to zero. This value is equivalent to what is commonly used for sea ice algae (Tedesco et al., 2010).

Results
The results presented here fall under two general themes. In the first section we will test the melt water hypothesis outlined in Section 1, by comparing the timing of growth initiation (GI) with that of sea ice retreat. Following this we will present results from a set of simple model experiments in an attempt to explain for the observed phenology. In these experiments we 8 https://doi.org/10.5194/bg-2020-257 Preprint. Discussion started: 17 July 2020 c Author(s) 2020. CC BY 4.0 License. investigate the role sea ice cover and phytoplankton low-light adaptation play in controlling winter/spring growth. By placing 195 the experiments in 4 distinct study regions with different physical conditions, we also utilize the spatial and temporal variability available in the float dataset to derive results of wider regional applicability.
9 https://doi.org/10.5194/bg-2020-257 Preprint. Discussion started: 17 July 2020 c Author(s) 2020. CC BY 4.0 License. In Figure 4A we show more explicitly the relationship between growth initiation timing and latitude. Here we find a statistically significant correlation (p = 0.003) of -0.44, implying that 19% of the variance in GI may be explained by variability in latitude alone. Conversely, the relationship between GI and the timing of melt water release is insignificant at the 5% level 210 (p = 0.08) with a lower correlation of 0.27 ( Figure 4B). Furthermore, almost all events fall below the 1:1 line in Figure 4B, revealing that GI tends to precede the release of melt waters.

Observed Under Ice Growth
Consequently, in Figure 5 we plot the distribution of the difference in timing between GI and melting. For the majority of the observed events, GI occurs well before the release of melt waters (the mean timing difference is 7 weeks). Furthermore, for 35% of the events, GI is observed more than 35 days before melting, with a further 25% preceding melting by 25-35 days.

215
Only 10%, or 4 events, occur either at the same time as or after sea ice retreat.
We would note that our definition of GI (detailed in section 2.2) is likely to be more conservative than methods employing a threshold value (i.e likely to delay growth initiation). In addition, GI was also computed using the mean mixed layer particulate organic carbon as opposed to chl-a, which resulted in a timing difference distribution similar to that shown in Figure 5 (albeit with a smaller time difference, see Supplementary Figure S6). In terms of vertical mixing, average stratification depth at GI is 220 129m. Table 2 highlights some of the salient properties of the data set investigated in this study, as well as summarizing the major findings discussed above.
While the results discussed up to this point incorporate data from all available under ice floats, in Figure 6       Three core experiments were conducted for each region, consisting of first running with no sea ice forcing (OPEN), then with satellite derived ice concentration (ICE), and finally with the low light efficiency of phytoplankton enhanced by a factor of 10 (LLA; sea ice forcing is also kept for these runs). Within each of the four study regions, this set of experiments is conducted in the same manner as in the float data (see section 2.2), although there was no need for filtering. While the LLA set of experiments generally performs best at reproducing GI, there are notable exceptions in each of the 4 study regions. In the W60 region the observed GI occurs between early September and mid-October, with OPEN experiments having growth too early and LLA experiments too late. Moving further south to W65, we see that only LLA is able to capture the observed variability in GI, but that in some cases OPEN provides the best fit to data. Continuing south and west, both B70 and R75 contain cases 275 where GI is best described by ICE simulations. In the following section we will bring together both the observational and modelling results discussed thus far, thereby shedding light on the possible mechanisms leading to under ice growth in the Antarctic winter and spring.

280
The central question of the present study relates to what conditions are necessary to trigger phytoplankton growth in the Antarctic SSIZ. As has been outlined in Section 1, a popular hypothesis holds that the release of buoyant melt waters following sea ice retreat shoals the mixed layer, relieving light limitation and triggering rapid growth. In contrast to previous studies (e.g. Smith and Nelson (1985); Smith and Comiso (2008); Sokolov (2008); Taylor et al. (2013)) relying on satellite data or models, we were able to thoroughly test this hypothesis by utilizing a unique in-situ dataset of under ice profiles from bio-ARGO floats.

285
In particular, we were able to test two predictions of the hypothesis; first, that at least part of the variability in the timing of growth initiation (GI) may be explained by the timing of sea ice melt, and second, that GI should either be synchronous with or occur after the release of melt waters.
Based on the data analysed here, we do not find evidence which convincingly supports either claim. Figures 5 and 6, we clearly demonstrate that phytoplankton are able sustain growth long before significant freshening of the surface ocean. It is 290 important to reiterate here that GI is based on the rate of growth exceeding the median rate, and so the tendency of GI to precede melting (as illustrated by the timing differences between these events shown in Figure 5) suggests that the rate of growth is already well above average prior to ice retreat. This explains why GI and melting are not correlated in time ( Figure   4B); the release of melt waters does not appear to relieve light and/or nutrient limitation and so variability in melt timing cannot account for variability in GI. GI is instead correlated more strongly with latitude ( Figure 4A), suggesting that phytoplankton 295 are responding to changing incident light conditions rather than fresh water fluxes. To be clear, the latitudes plotted in Figure   4A are computed based on the approximate location of the float at GI, which in almost all cases corresponds to an under ice condition. Therefore, the correlation found in this figure implies that light may be non-limiting under Antarctic sea ice (at least in the conditions sampled by the floats), provided it is late enough in the season for there to be sufficient light available at the surface.

300
Also noteworthy is the extent to which growth occurs prior to melting, with ∼60% of events preceding melting by a month or more. As is discussed in 3.1, our float dataset samples in a wide variety of environmental conditions, which exhibit very different sea ice and vertical mixing regimes. This suggests that the results presented here are fairly representative of the SSIZ In our regional box model experiments we explore both physical and biological factors. The fact winter and spring phenology is brought closer to observations when low light efficiency is enhanced by an order of magnitude (to a value typical of sea ice algae) certainly suggests a role for phytoplankton adaptation (Figure 7 -LLA experiments). However, the interpretation is 340 complicated somewhat by the fact that under certain conditions phenology may be best described by simulations with no ice (OPEN) or with ice but standard physiology (ICE).
For example, in the Weddell Sea ( Figure 7A and B) early growth in August is best captured by OPEN experiments, but subsequent spring growth rates (October-November) more closely align with LLA simulations. The inference here would be that in this region sea ice is unconsolidated and highly permeable to light, allowing growth to initiate as soon as incident 345 radiation is sufficient. This corresponds well with the correlation between GI and latitude shown in Figure 4A. This is despite the apparently near 100% sea ice concentration suggested by satellite data (see Figure 6 and Supplementary Figures S5-S7).
Indeed, at these latitudes we may actually be in the MIZ, which would explain the higher light permeability. Yet, this is not to say that sea ice has no effect, later in the season growth rates are slowed by its presence, explaining why LLA experiments perform better here. These findings generally agree with previous studies which point to light (as opposed to dissolved iron) Further south in the Bellingshausen and Ross seas, sea ice is expected to be more consolidated in winter and spring, and so phenology is better captured by LLA simulations. However, in two cases the timing of GI most closely matches ICE experiments (see Figure 8, regions B70 and R75). This may be accounted for by especially thick snow and ice layers in those 355 cases, which led to delayed growth. This highlights the importance of the particularities of ice morphological features and their effect on the light environment, something which does not seem to be captured by satellite sea ice concentration.
Thus, it is both the character of ice and snow overhead, and the physiological response to severe light limitation that may address the question raised at the start of this section. A crucial point here is that 100% sea ice cover (in the winter Antarctic sea ice) as seen from satellite does not necessarily imply a completely consolidated ice surface (Vichi et al., 2019). While the 360 ocean may indeed be completely covered, the ice itself may be unconsolidated, being primarily composed of pancakes loosely connected by frazil or brash ice. Such a condition is common in the Southern Ocean, and is maintained by wind and wave action far from the ice edge. Waves are known to propagate several hundred kilometres into the ice, effectively preventing the formation of pack ice-like conditions Meylan et al., 2014). Wind forcing is also known to be highly effective in causing ice break-up and motion, with intense synoptic events in the Weddell and Eastern Indian oceans occurring 365 frequently (Vichi et al., 2019;Uotila et al., 2000). Such events, along with interactions with the westerly wind belt, drive the formation of gaps within the MIZ, as well as within pack ice. Therefore, the highly dynamic nature of Antarctic sea ice may lead to a general enhancement of light availability in the underlying ocean. The presence of even a tiny amount of light may be expected to induce acclimation in primary producers that are adapted to low light, thereby explaining why the model configurations presented here which take this into account produce a more realistic phenology.

370
This study has characterised under ice phytoplankton phenology using a unique dataset of Bio-ARGO profiles, complemented by a set of process-oriented biogeochemical model experiments. We have shown that rather than acting as a trigger as postulated in previous studies, the release of melt waters enhances growth in an already highly active phytoplankton population. This may explain the decline in phytoplankton stocks observed by Veth et al. (1992) in melt water lenses of the north-western Weddell 375 Sea. That is, the decline (in a still highly stratified surface ocean) may be accounted for by the natural reduction occurring in a bloom that already started prior to the melting. Such unexpected early growth (under presumed severe light limitation) may be accounted for by a combination of low light adaptation by phytoplankton and sea ice permeability with respect to light. We argue that such permeability is related to wind and wave forcing, which together preserve an unconsolidated ice morphology that is not captured by current satellite sea ice concentration algorithms.

380
However, our investigation has not been exhaustive of all possible mechanisms leading to under ice growth. Future research directions could include an examination of potential discrepancies between the timing of shoaling of the mixed layer and that of active turbulent mixing (e.g. Carranza et al. (2018); Sutherland et al. (2014)). An earlier reduction in mixing would increase ambient light, and help explain the observed under ice growth. Other ecological factors could also be explored, such as potential interactions between pelagic and sympagic communities, which are known to be highly efficient at low light intensities 385 (Tedesco and Vichi (2014) and citations therein). Nevertheless, the findings presented here have important implications for our understanding how the biogeochemistry of the region may change in the future. With possible earlier sea ice retreat, as well as a generally thinner and more dynamic ice in some regions (including the Arctic), we may expect even earlier growth then reported here, which would likely alter the seasonal air-sea carbon flux and thus the biological carbon pump.
Author contributions. Mark Hague conducted the float data analysis, performed the model experiments and wrote the manuscript. Marcello