Sea-ice derived meltwater stratification slows the biological carbon pump: results from continuous observations
Table S1 lists the data used in this paper, the instruments that it is based on, the data repositories, and in which figures the data are used.Global data setsBathymetryBathymetric data was taken from the International Bathymetric Chart of the Arctic Ocean (IBCAO 30 sec V3)64 available at https://www.ngdc.noaa.gov/mgg/bathymetry/arctic/grids/version3_0/.Sea ice concentrationWe use data derived from the Advanced Microwave Scanning Radiometer sensor AMSR-2 for the years 2013–18 processed in accordance with65 and downloaded from https://seaice.uni-bremen.de/sea-ice-concentration-amsr-eamsr2/66. At each grid point the sum of days during all April/May/June of 2013–2018 when the sea ice concentration at the grid point was >20% was divided by the total number of days with data in those months to obtain the percentage of days with ice concentration >20% (Fig. 1). For separate 7-day periods in April/May/June 2017 and 2018 the mean ice concentration over those 7 days was calculated and the 20% contour of this mean was plotted separately for each of those 7-day periods. For each mooring and each day, the ice concentration at the grid cell closest to the mooring was calculated (Fig. 4a and S1a), and if the ice concentration at the mooring was below 20%, the shortest distance to grid cells where the ice concentration exceeded 20% was calculated (Fig. 4a and S1a). If the ice concentration at the mooring exceeded 20%, the shortest distance to grid cells where the ice concentration was below 20% was calculated and the distance was defined as negative.Sea ice velocity and sea ice area exportIce area flux estimates in Fig. 2a are calculated using CERSAT (Center for Satellite Exploitation and Research, France) motion estimates together with CERSAT ice concentration information67. Fluxes are estimated along a zonal gate positioned at 82°N between 12°W and 20°E and a meridional gate at 20°E between 80.5°N and 82°N (Fig. 1) for the period 1994–2020 (January–May). The ice area flux at the gate is the integral of the product between the meridional and zonal ice drift and ice concentration. For a more detailed description we refer to ref. 68. Arctic-wide sea ice velocity anomalies (Fig. 2b, c) were computed from the OSI-405-c motion product provided by the Ocean and Sea Ice Satellite Application 635 Facility (OSISAF)69.Satellite chlorophyllSurface chlorophyll concentrations measured with the Sentinel 3 A OLCI (Ocean and Land Colour Instrument) were downloaded from https://earth.esa.int/web/sentinel/sentinel-data-access. The 8-day satellite data were averaged for the time series over grid points within boxes of 60 km by 60 km around the moorings.Atmospheric reanalysisERA-Interim reanalysis70 data at the surface on a 0.25° latitude by 0.25° longitude grid at 12 hourly resolution was downloaded from https://apps.ecmwf.int/datasets/data/interim-full-daily/levtype=sfc/. Incoming shortwave radiation (ssr) and outgoing longwave radiation (str), sensible heat flux (sshf), and latent heat flux (slhf) were extracted and averaged to daily values.Physical numerical modelsFESOMIn this study, we used model data from the Finite-Element Sea ice-Ocean Model (FESOM) version 1.471. FESOM is a sea ice-ocean model that solves the hydrostatic primitive equations for the ocean and comprises a finite element sea ice component. It uses triangular surface meshes for spatial discretization, allowing for a refined mesh in regions of interest, while keeping a coarser mesh elsewhere. In the model configuration used here, a mesh resolution of nominally 1° was applied in the global oceans. The mesh was refined to 25 km north of 40°N, and to 4.5 km in the Nordic Seas and Arctic Ocean. In the wider Fram Strait (20°W-20°E/76°N-82°30′N), the mesh was further refined to 1 km. In this region, the simulation can be considered as eddy-resolving, as the local internal Rossby radius of deformation is about 2–6 km72,73. In the vertical, the model used 47 z-levels with a resolution of 10 m in the upper 100 m, and coarser resolution with depth (with a resolution of ~100 m at 800 m depth). For bottom topography, the RTopo-2 data set was used74. The model simulation covers the period 2010–2018 and has daily model output. It was forced with atmospheric reanalysis data from Era-Interim70, and was initialized with model fields from the simulation described in ref. 75. River runoff (except for Greenland) was taken from the JRA-55 data set76, and Greenland ice-sheet runoff was taken from ref. 77. Tides were not taken into account in this simulation. Here we studied the model data of 2016 to 2018 in Fram Strait for comparison with our observations.1-dimensional mixed layer depth modelThe PWP78 1-dimensional mixed layer model simulates the response of the ocean to surface fluxes. It ignores horizontal gradients and horizontal advection. This allows to judge whether certain surface flux conditions can on their own explain observed conditions. We ran the PWP model (as implemented for Matlab by http://www.po.gso.uri.edu/rafos/research/pwp/) with four different scenarios (Fig. S6: P17-M17, P17-M18, P18-M17, P18-M18) where: P17: An idealized initial profile based on the observed profiles (Fig. 3) representing the conditions in 2017: constant temperature of 2 °C in the vertical, linear salinity gradient from 30.5 at the surface to 35 at 50 m and another linear salinity gradient from 35 at 50 m to 35.1 at 200 m. P18: An idealized initial profile based on the observed profiles (Fig. 3) representing the conditions in 2018: Same as P17 except that the surface to 50 m salinity gradient is from 34.8 to 35. M17: A time series of the the meteorological forcing (10 m wind velocity, heat fluxes, and evaporation minus precipitation) from the ERA-Interim reanalysis (Fig. 4b) at the grid point closest to mooring HG-IV for the period 15-May-2017 to 01-Aug-2017. M18: Same as M17 but for the period 15-May-2018 to 01-Aug-2018. M17 and M18 are provided in Supplementary Data 1.Shipboard CTD dataShipboard CTD casts of a standard dual sensor Seabird 911+ CTD-rosette were occupied in spatial and temporal vicinity to the moored observations (Tab. S2) on three cruises: PS107 in 2017 (https://doi.org/10.1594/PANGAEA.894189), PS114 in 2018 (https://doi.org/10.1594/PANGAEA.898694) of RV Polarstern, and JR17005 in 2018 (https://doi.org/10.5285/84988765-5fc2-5bba-e053-6c86abc05d53) of RRS James Clark Ross. The data were processed according to standard routine79. Additionally, we use underway CTD data from an OceanScience underway CTD collected during PS107 in 2017 (https://doi.org/10.1594/PANGAEA.886146) and processed according to ref. 21.Mooring dataThe mooring data discussed in this paper is from two mooring clusters in the central and eastern Fram Strait (named “HG-IV” at ~79°N 4°20’E and “F4” at ~79°N 7°E) where moorings were located as close to each other as possible (the horizontal separation was equal to the water depth) in order to enable more measurements than could be fit physically onto a single mooring. Tab. S2/S3 list the deployment and recovery details of the moorings including the exact latitudes/longitudes as well as the individual instruments on the moorings. Note that all data shown in this paper from ~30 m depth and the temperature/salinity/oxygen data from ~55 m is from the HG-IV-S-* and F4-S-* moorings, while all other data is from the HG-IV-FEVI-* and F4-* moorings. The AZFP data is from F5-17 located roughly half way between the two clusters. All sensor based mooring raw data (except for the ASL AZFP data) is available at ref. 80.It is known that conversion factors for biogeochemical sensors (e.g., chlorophyll fluorescence) change over the seasons, depths, and regions81,82. In order to make as few assumptions as possible, we used the following approach: we could have determined the conversion factors from the instance when the ship was there with the CTD-rosette, but these conversion factors might not be appropriate for the majority of the time series. Hence, simply using the manufacturers’ calibrations, as we do here, introduces fewer uncertainties. Where we have different estimates of the same parameter, we present them together and demonstrate that they agree qualitatively and also mostly quantitatively (e.g., Fig. 5b). In particular the timing of events is robust.At some locations, the target variables were not measured the whole time or the measurements failed, hence we present what is available. The vertical location of the instruments (Fig. 4c and S1c) varied substantially (intermittently up to 200 m) as a result of mooring blow downs caused by strong intermittent ocean currents. Time series have not been corrected for this vertical motion, but data are not used during blow downs in order not to bias the time series interpretation by temporal changes introduced by instruments traversing through vertical property gradients.Physical sensor measurementsThe physical sensors (for pressure, temperature, conductivity, and oxygen) were pre-cruise manufacturer calibrated and processed similar to ref. 83; the processed data is also available at ref. 80.Mixed layer depth (MLD)Since there are no autonomous vertically profiling measurements available, we can only determine the minimum value of the mixed layer depth. At each hourly time step, the potential density difference (Δσ) between the uppermost (~30 m) temperature/salinity recorder and the underlying temperature/salinity recorders is calculated. The 0.5th percentile of each Δσ time series is added to the Δσ time series for the different deployments. This fixes slight offsets in the temperature and/or conductivity calibrations which result in too negative or too positive density differences. The minimum estimate of the mixed layer depth at hourly resolution is then determined as the depth of the deepest instrument where Δσ 0.05 kg m−3 for all depths at a time step, then the minimum mixed layer depth can only be determined as 0 for that time step. Daily values of the MLD were defined as the depth at which three hourly realizations of MLD were shallower within a 24 h time span and at which the remaining 21 MLD realizations were deeper. This biases the daily MLD estimate towards situations where phytoplankton is kept in the surface ocean rather than also being mixed down for some amount of time.Stratification estimated between 30 m and 55 mBased on the temperature and salinity time series observed at ~30 m and ~55 m, we estimate the buoyancy frequency as ({N}^{2}=frac{-g}{{rho }_{0}}frac{Delta rho }{Delta z}) where g is the acceleration due to gravity, Δσ is the potential density difference over the vertical distance of Δz = 25 m, and ρ0 is the average density. The contributions to stratification due to temperature (N2T) and salinity (N2S) are estimated as ({N}_{T}^{2}=g*alpha frac{Delta T}{Delta z}) and ({N}_{S}^{2}=-g*beta frac{Delta S}{Delta z}), respectively, where ΔT/ΔS are the temperature/salinity differences and α/β are the thermal expansion/haline contraction coefficients estimated from the average temperature/salinity at the two measurement depths.Apparent oxygen utilization (AOU)Oxygen concentration from the microcats was calculated using the pre-cruise manufacturer calibrations. AOU was calculated as the atmospherically equilibrated oxygen concentration (calculated from measured pressure, temperature, and salinity with sw_satO2 from the Seawater toolbox available at http://www.cmar.csiro.au/datacentre/ext_docs/seawater.htm) minus the measured oxygen concentration.LightPolar night/polar dayThe length of day (hours per 24 h that the sun is above the horizon) was calculated from the sunrise equation as implemented for Matlab by https://de.mathworks.com/matlabcentral/fileexchange/55509-sunrise-sunset.Photosynthetically available radiation (PAR)The WetLabs Eco PAR measured PAR for 5 (in 2016–2017) or 10 (in 2017–2018) individual measurements 1 s apart from each other before it slept for 1 h before repeating the measurement cycle. These 5 or 10 individual measurements are averaged linearly to obtain hourly values at ~30 m depth (Fig. 5a blue). Values below the detection limit are set to a constant of 10−1.32 μmol m−2 s−1. Hourly values are linearly averaged to daily values (Fig. 5a black). The incoming solar shortwave radiation varies as a function of season and latitude as well as cloud cover as represented in the ERA-Interim reanalysis (parameters ssr). Its unit of W m−2 is converted to PAR assuming a constant spectral distribution as 1 W m−2 = 2.1 μmol m−2 s−184. In order to compare the PAR measured at a depth of approximately 30 m to the surface values, we approximate a spectrally averaged diffuse attenuation coefficient for PAR in clear water using the values of85 as kd = 0.02 m−1 and apply it to calculate a constant exponential extinction applied to the reanalysis surface values (Fig. 5a yellow). The average PAR available (PARavailable) to phytoplankton being moved around in the clear water mixed layer of depth MLD was calculated as the depth averaged vertical integral of the clear water extinguished PAR at the surface (PARsurf from the shortwave radiation of ERA-Interim): ({{PAR}}_{{available}}=frac{1}{{MLD}}*{int }_{z=0}^{z={MLD}}{{PAR}}_{{surf}}*{e}^{-{k}_{d}z}{dz}) (Fig. 5a red).Chlorophyll concentration and optical backscatteringChlorophyll fluorescenceThe WetLabs ECO Triplet measures fluorescence at a “chlorophyll wavelength” and at a “CDOM wavelength” as well as optical scattering at 700 nm. The conversion from fluorescence to chlorophyll a concentration (in μg l−1) follows a manufacturer determined conversion determined for a mono-culture of phytoplankton (Thalassiosira weissflogii), which typically overestimates the chlorophyll concentration. Hence, we applied the community-established calibration bias of 2 for the WetLabs ECO-series fluorometer to these in situ fluorometric chlorophyll values81. This conversion factor may be different in ocean waters of Fram Strait, but it still gives reasonable agreement with independent estimates.Optical backscatteringThe EcoTriplet measured 8 individual measurements 1 s apart from each other before it slept for 1 h before repeating the measurement cycle. For the chlorophyll fluorescence, the individual measurements are averaged to hourly values. For the scattering, times when individual 1-second measurements exceed 0.002 m−1 sr−1 are indicative of strong optical backscattering not due to small particles in the water column, but rather to larger potentially aggregated particles. The times of strong backscattering are marked individually (Fig. 5b red).NutrientsNitrate (SUNA sensor)Prior to deployment (11 and 15 days for sensors deployed at HG-IV and F4, respectively), the reference spectrum of the sensors were updated as per manufacturer specifications. We first let the sensors cool down for 24 h at 0 °C in a temperature controlled laboratory. Next, the reference spectrum update was achieved by measuring Milli-Q water (i.e., no nitrate present). To verify if this update was successful, solutions with three different nitrate concentrations (3, 7, and 14 μmol l−1) were then measured, with the output being monitored live (expected to be within ±2 μmol l−1 of each concentration). A measuring time of 20 s yields stable results and was thus applied during the deployments with an interval of 6 h. Upon recovery, SUNA data were processed using the SeaBird UCI software package version 1.2.1. Here, temperature and salinity data were used to remove the spectrum of bromide and compensate for temperature dependent absorption using an algorithm developed by ref. 86. This step yields the spectrum of nitrate only, at a precision of ±0.3 μmol l−1. The sensor is characterized by a drift of 0.3 μmol l−1 per hour lamp time. Given the deployment settings, a total operational time of about 8 h was accumulated. Therefore, a linear drift correction of 2.4 μmol l−1 (365 days)−1 was applied. Up to this point, however, accuracy remains at 2 μmol l−1 as per manufacturer specifications. Therefore, an offset correction is then applied based on the in situ concentrations observed at the beginning of the deployment as well as with the RAS (see below) where available, with outliers excluded.Inorganic nutrients from Remote Access Samplers (RAS)McLane RAS were programmed to draw two 500 ml samples (1 h apart, starting at noon) approximately every other week. Samples within the RAS were collected in sterile plastic bags and fixed with 700 μl of 50% mercuric chloride solution. Upon recovery, two samples from a given sampling date were combined to yield a volume of 1 l, required for bacterial and phytoplankton genetic analyses (see below), and a 50-ml aliquot destined for the measurement of dissolved inorganic nutrients. Aliquots for nutrient analysis were collected in PE bottles, which were then stored frozen (−20 °C) until analysis on land. Analyses for inorganic nutrients were carried out using a QuAAtro Seal Analytical segmented continuous flow autoanalyser following standard colorimetric techniques. The accuracy of the analysis was evaluated through the measurement of KANSO LTD Japan Certified Reference Materials and corrections were applied accordingly. Finally, we evaluated pressure, temperature, and salinity data from the CTD (SBE37-SMP-ODO) attached to the RAS to determine whether the two samples taken one hour apart on a given date drew water from the same depth and with consistent properties.Carbonate system
pCO
2 and pH
The calibration of SAMI pH and SAMI CO2 sensors was carried out by the manufacturer, approximately 2 months prior to deployment. The calibration certificates specify accuracy and precision of ±0.003/±0.001 pH units and ±3/ More
