Widespread surface water pCO₂ undersaturation during ice-melt season in an Arctic continental shelf sea (Hudson Bay, Canada)

Global oceans slow the effects of climate change by absorbing anthropogenic carbon dioxide (CO₂) emitted by the burning of fossil fuels, land-use change, and cement production. However, the invasion of this anthropogenic CO₂ into the ocean has caused a decline in surface ocean pH leading to ocean acidification (Orr et al., 2005; Gattuso and Hansson, 2011). Friedlingstein et al. (2019) quantified the major components of the global carbon budget (atmosphere, land, and ocean reservoirs) over the previous decade (2009–2018) and estimated that the global ocean absorbed about 23% of the total anthropogenic CO₂ emissions. However, this global budget neglected the carbon cycle in coastal ocean waters largely because the spatial and temporal heterogeneity of biogeochemical processes in these waters make it difficult to accurately scale up results from observational studies. Despite great efforts over the last decade to assess the coastal contribution to the global carbon budget (e.g., Chen and Borges, 2009; Bauer et al., 2013; Regnier et al., 2013; Laruelle et al., 2014; Bourgeois et al., 2016; Roobaert et al., 2019), the CO₂ budgeting of these waters is still subject to uncertainty and lacks consensus. The coastal ocean has been impacted significantly by human activities (hydroelectric damming, eutrophication, overfishing, etc.) and is now considered one of the most sensitive parts of the marine environment to climate change (e.g., Arctic Climate Impact Assessment, 2005; Liu et al., 2010; Liao et al., 2015; Breitburg et al., 2018). Furthermore, the coastal shelves play a disproportionately significant role in the global carbon cycle in comparison with its surface area (about 7%–10% of the oceanic surface area) by linking the oceanic, terrestrial, and atmospheric carbon reservoirs (Wollast, 1998; Liu et al., 2010). Hence, studying the carbon cycle in these waters is essential to understand future shifts in response to climate change.

The coastal ocean consists of diverse but closely connected ecosystems that include estuaries, rivers, tidal wetlands, and the continental shelf. The continental (or coastal) shelves are shallow and have been the focus of many studies due to their key role in the global biogeochemical cycling of carbon and nutrients (e.g., Chen et al., 2013; Bourgeois et al., 2016) as well as their importance to local livelihoods. For instance, the ocean’s biological pump operates more efficiently in coastal shelves relative to the open ocean as a result of high terrestrial inputs and efficient use of nutrients (Gattuso et al., 1998; Chavez et al., 1999). The high biological production in these waters usually leads to atmospheric CO₂ drawdown and subsequent CO₂ export to the subsurface waters via the biological pump (Ducklow et al., 2001). The eventual outflow of CO₂-enriched deep water from coastal shelf regions to the subsurface waters of the adjacent deep oceans constitutes “the continental shelf pump,” which is thought to contribute significantly to the global ocean’s atmospheric CO₂ uptake (Tsunogai et al., 1999). Furthermore, the production, export, and degradation of organic matter in these waters are several times higher than in the open ocean (Wollast, 1998; Borges et al., 2005). The coastal shelves contribute 10%–30% of global marine primary production, 80% of organic carbon, and 30%–50% of inorganic carbon burial in sediments, and about 50% of the organic carbon supplied to the deep open ocean (Wollast, 1998; Liu et al., 2010). As a result, the sea–air CO₂ exchanges are expected to be more intense in these waters than in the open ocean.

Coastal shelves show a latitudinal pattern in CO₂ fluxes, with the mid-to-high-latitude (30°–90°) shelves acting as oceanic sinks of CO₂, whereas the low latitude shelves (equator to 30°) act as sources of CO₂ to the atmosphere or are nearly neutral (Borges et al., 2005; Cai et al., 2006). This latitudinal pattern can be attributed to a complex interaction among temperature (which alters the solubility of CO₂), the amount of biological activity, the magnitude of upwelling and mixing from deeper waters, and the export of terrestrial organic and inorganic carbon from rivers (Takahashi et al., 2009; Dai et al., 2012). Adding further complexity, the processes controlling CO₂ dynamics vary with proximal distance from shore (Wollast, 1998; Bauer et al., 2013). Inner shelf waters close to land tend to be sources of CO₂ due to respiration of terrestrial organic carbon and lateral transport of CO₂-enriched waters from the adjacent inshore system, whereas the mid-to-outer shelf waters generally act as a sink of CO₂ due to increased primary production as a result of nutrient availability from upwelling and mixing across the shelf break (Walsh, 1991; Jiang et al., 2008). Moreover, coastal shelves show a temporal variability in sea–air CO₂ fluxes as a result of variability in riverine input, phytoplankton blooms, ice formation and melting, and sea surface temperature (e.g., Takahashi et al., 2009; Rysgaard et al., 2012; Ahmed et al., 2019). For example, Borges et al. (2005) documented seasonal changes of surface seawater CO₂ partial pressure (pCO₂) as high as 400 μatm in the coastal shelves as a result of riverine input and other coastal processes.

In high-latitude coastal seas, changes in sea-ice coverage (Notz and Stroeve, 2016), water temperature (Shepherd, 2016), riverine input (Carmack et al., 2016), glacial melt input (Pierre et al., 2019), and primary production (Arrigo and van Dijken, 2015) will have a significant impact on the distribution and magnitude of pCO₂ and thus the strength of atmospheric CO₂ uptake. Given these potential changes, as well as the spatial and temporal variability in pCO₂, studying the spatiotemporal patterns of pCO₂ independently in the Arctic and sub-Arctic coastal shelves is important to improve our knowledge of the main drivers affecting pCO₂ and to increase the reliability of sea–air CO₂ flux budgets.

The Hudson Bay system (including James Bay, Foxe Strait, Ungava Bay, and Hudson Strait) is a sub-Arctic/Arctic continental shelf sea that receives nearly one-third of Canada’s river discharge from more than 42 rivers (Déry et al., 2011). Hudson Bay and its watershed have been recognized as highly sensitive to climate variability relative to other Arctic shelves (Arctic Climate Impact Assessment, 2005; Markus et al., 2009). For instance, several studies have documented changes in timing and volume of river discharge (Déry et al., 2011, 2018), the timing of sea-ice formation and melting (Saucier et al., 2004; Hochheim et al., 2011), and biological production (Sakshaug, 2004; Sibert et al., 2011). Permafrost melting due to climate change may also profoundly affect the quality, quantity, and composition of organic matter in Hudson Bay (Gough and Leung, 2002; Mallory et al., 2010; Godin et al., 2017). These potential changes could have a drastic impact on sea–air CO₂ fluxes in Hudson Bay.

Although included in several global coastal ocean CO₂ budgets (e.g., Chen and Borges, 2009; Chen et al., 2013), Hudson Bay has not been included in other global and Arctic-specific budgets (Cai et al., 2006; Bates and Mathis, 2009; Laruelle et al., 2014). This limited inclusion is due to a lack of observations on spatial and temporal scales. Only two studies investigating sea–air CO₂ fluxes in Hudson Bay have been published (Else et al., 2008a, 2008b). Using pCO₂ measurements at discrete stations collected in September and October, Else et al. (2008a) found that the bay acts as a minor sink in October, while Else et al. (2008b) applied a remote sensing approach to determine that Hudson Bay acts as a weak CO₂ source during August and September. The limited temporal and spatial extent of those observations means that we have no reliable estimate of sea–air CO₂ exchange in Hudson Bay over the complete open water period. Therefore, our baseline estimates in this large inland sea need to be extended, especially into the spring and early summer seasons, when the ecosystem is likely subject to rapid and pronounced physical and biological changes (e.g., Else et al., 2012).

Here, we report the first spring and early summer measurements of pCO₂ and estimates of sea–air CO₂ fluxes in Hudson Bay and Hudson Strait. The main objectives of this study were to (1) examine the spatial and temporal patterns of pCO₂; (2) identify the key processes driving the variability in the observed pCO₂; and (3) estimate the sea–air CO₂ fluxes during the 2018 spring and early summer seasons in Hudson Bay and Hudson Strait.

Hudson Bay is a large shallow inland sea located in central Canada, with an average depth of 125 m and a maximum depth of 250 m, while adjacent areas of Foxe and Hudson Straits can reach depths of 400 m (Prinsenberg, 1986; Saucier et al., 2004). Hudson Bay transitions from complete ice cover in the winter to open water in summer, which adds about 140–160 cm of freshwater when sea-ice melts from May to July (Prinsenberg, 1988). Freshwater from river runoff contributes about 85 cm, distributed around the perimeter of the bay (Prinsenberg, 1986), and largely constrained near the coast (Granskog et al., 2007). Peak river runoff occurs in May, whereas the minimum occurs in March, with a notable change in recent decades toward higher winter discharge for hydroelectric power generation (Déry et al., 2011). This massive seasonal freshwater input to the bay leads to a strong stratification in the upper surface layer (Granskog et al., 2011; Burt et al., 2016; Ahmed et al., 2020), which impacts the mixing of heat and nutrients between the surface and deep waters (Drinkwater and Jones, 1987; Prinsenberg, 1988). Water flowing out of the Canadian Arctic Archipelago enters Hudson Bay at its northwestern end through Foxe Strait, while the northern portion of Hudson Strait allows water inflow of both Atlantic and Arctic origin from the Labrador Sea. Water exits Hudson Bay via the southern portion of Hudson Strait (Ingram and Prinsenberg, 1998). The general circulation in the bay is cyclonic (counterclockwise), which can be attributed to wind-driven transport and estuarine circulation buoyancy resulting from river runoff and ice melt/formation (Prinsenberg, 1986; Saucier et al., 2004). In this study, we focus on Hudson Bay (approximately 0.8 × 10⁶ km²) and Hudson Strait (approximately 0.2 × 10⁶ km²), where we conducted the majority of our field measurements.

We collected our field measurements during a 6-week survey from May 25 to July 13, 2018, onboard the research icebreaker CCGS Amundsen (Figure 1). This study was part of the BaySys project (see other papers in this special feature of Elementa), a comprehensive study examining the relative influences of climate change and hydroelectric regulation on Hudson Bay’s physical, chemical, and biological systems. The field program was the first-ever research cruise in Hudson Bay in the spring and early summer, at a time when sea-ice cover was still partially intact, and sea-ice melt and river runoff are near their peaks.

Surface water pCO₂ was measured during the cruise using an underway pCO₂ system (General Oceanics model 8050; Pierrot et al., 2009). This automated pCO₂ system continuously sampled water through a shower-type equilibrator system at a nominal flow rate of 2.5 L min^–1 from a through-hull inlet located approximately at 7 m beneath the water surface. The precision of the underway pCO₂ system is 0.1 μatm, with an overall accuracy estimated at 2 μatm (Pierrot et al., 2009). The system was calibrated every 8 h using four certified standard gases (ultra-high purity N₂ as a zero gas and three CO₂/air mixtures at concentrations of 295.8, 448.6, and 551.7 ppm) traceable to World Meteorological Organization standards. The system was very stable during our sampling time, with drift typically less than 0.2 μatm between calibrations, contributing to uncertainty in our measurements less than 0.5%. We calculated the seawater pCO₂ using the mole fraction of CO₂ (xCO₂) based on Dickson et al. (2007) and following the recommendation of Pierrot et al. (2009) for pCO₂ calculations and xCO₂ corrections using the measured standards.

Using discrete water samples, Ahmed et al. (2020) found that the sea-ice melt created a shallow, stratified freshwater layer that played a major role in diluting the surface pCO₂ across Hudson Bay. To correct for vertical pCO₂ gradients resulting from freshwater stratification during our study, we applied a linear regression that predicts the pCO₂ gradient as a function of observed pCO₂ (fully described in Ahmed et al., 2020) to about 45% of the total underway measurements (those measured 1–5 weeks after ice breakup with wind speed less than 6 m s^–1). For pCO₂ in areas measured less than 1 week after ice breakup, stratification was considered to be too strong for confident corrections, and those records (about 10% of the total measurements) were discarded. For pCO₂ in areas measured more than 5 weeks after ice breakup or areas 1–5 weeks after ice breakup with wind speeds greater than 6 m s^–1, no correction was applied (45% of the measurements) as the near-surface vertical pCO₂ and salinity gradients were not likely significant under those conditions (Ahmed et al., 2020).

The underway pCO₂ system has flow-through temperature, salinity, and oxygen sensors (Idronaut model Ocean Seven 315), which were factory-calibrated before deployment. However, we opted to use the temperature and salinity recorded by a separate thermosalinograph (TSG) that draws water from the same intake line but is installed closer to the ship inlet compared to the underway pCO₂ system (distance difference between two instruments in the engine room was about 3 m). A chlorophyll a (Chl a) fluorescence sensor (WestStar fluorometer with an accuracy of 0.03 µg L^–1) and a fluorescing dissolved organic matter (FDOM) sensor (WETLabs ECOFLD; 370 nm excitation, 460 nm emission; Guéguen et al., 2015) were also installed on the TSG line. Both the FDOM and Chl a sensors were also factory-calibrated before deployment. Measurements of surface nitrate (NO₃) concentrations and surface photosynthetic active radiation (SPAR) across Hudson Bay were collected using the Seabird Satlantic MBARI-ISUS V3 nitrate and the Biospherical Instruments SPAR (QCR 2200) sensors, respectively, which were mounted on the ship’s rosette. Additional details about sensors, calibrations, and data processing of the TSG and rosette data are reported in Amundsen Science Data Collection (2019).

Despite the proximity of the pCO₂ system to the inlet (<5 m), we observed a slight increase in the recorded water temperature in the equilibrator relative to the on-station measurements made by the ship’s conductivity-temperature-depth (CTD) sensor on the rosette. Therefore, we corrected the pCO₂ measurements for this temperature increase following the procedure described by Takahashi et al. (1993). Using the maximum uncertainty in the in situ seawater temperature measurements through the pCO₂ temperature correction, we estimate that the total uncertainty in the final pCO₂ measurements is less than 2%. We processed the underway data to remove any measurements collected when the system was undergoing cleaning, calibration, or restricted water flows from the intake line (i.e., water flow <1.5 L m^–1) during icebreaking operations. Further details about the data processing of the underway pCO₂ system are described in Ahmed et al. (2019).

Atmospheric CO₂ mixing ratio was measured every minute using a closed-path infrared gas analyzer (Li-COR model LI-7000). The analyzer was installed in a container on the foredeck and connected to a gas sampling system that drew air from an inlet located on the ship’s foredeck meteorological tower at a height of 17 m above the sea surface. The LI-7000 was calibrated twice a day with two certified gas standards traceable to World Meteorological Organization standards, and the mixing ratio measurements were converted into pCO_2atm following Dickson et al. (2007). Wind speed and air temperature were measured using a conventional propeller anemometer (RM Young Co. model 15106MA) and HMP45C212 temperature sensor, respectively. These sensors were mounted on the meteorological tower at a height of approximately 16 m above sea level. Wind speed was corrected to a height of 10 m assuming a log-linear wind profile and a neutral surface layer (Stull, 1988). We discarded the data when instruments on the tower were cleaned and/or calibrated or when the wind was blowing from the stern of the vessel.

For logistical reasons, the foredeck tower and associated sensors were removed 4 weeks into the cruise (on July 2, 2018). As a result, we opted to utilize wind velocities at 10 m height and surface air temperature provided by the National Centers for Environmental Prediction North American Regional Reanalysis (NARR) to have complete measurements along the ship track. The NARR wind and air temperature data are available eight times daily (i.e., every 3 h) with a horizontal grid spacing of about 32 km × 32 km (Mesinger et al., 2006). We clipped and resampled the NARR data using a bilinear algorithm to 1-km spatial resolution and then matched it with the ship track based on spatial location and closest estimate to our sampling time. The NARR and in situ measurements show good agreement over the period when the in situ sensors were deployed, especially for the wind speed data (Table S1). The atmospheric pCO₂ concentrations from May 25 to July 2, 2018, ranged from 400 to 419 μatm, with a mean of 408 ± 2.8 (1σ) μatm (n = 17,308).

Daily sea-ice observations along the ship track were obtained from the Advanced Microwave Scanning Radiometer 2 (AMSR2) passive microwave ice product (approximately 3-km spatial resolution) provided by the University of Bremen. The AMSR2 characteristics and the ARTIST Sea Ice algorithm used for calculating the daily ice concentrations are described in detail by Spreen et al. (2008). We used the spatial locations of our ship track to index the daily ice maps during our study time and then retrieve the sea-ice concentration from the closest temporally located observation. We compared the AMSR2 ice data with weekly ice charts from the Canadian Ice Service Digital Archive (CISDA) to ensure that the AMSR2 product did not underestimate ice concentrations along the ship track due to melt ponds on ice and/or weather conditions. The CISDA data are based on integration and analysis of data from satellite imagery (i.e., synthetic aperture radar and optical), weather and oceanographic information, and visual observations from ship and aircraft. Past studies have demonstrated the reliability of CIS data to monitor sea-ice cover change in the Canadian Arctic (Tivy et al., 2011) and its application to sea–air CO₂ exchange studies (Ahmed and Else, 2019). Figure S1 shows the sea-ice observations obtained from AMSR2, CISDA, and ship observations recorded by the coast guard officers on the ship. Although the AMSR2 product showed relatively better agreement with the ship observations compared to the CISDA (Figure S1), we noticed that the AMSR2 underestimated ice concentrations in a few locations in Hudson Bay, especially along the southwest and east coasts of Hudson Bay (June 27–30, 2018), where we observed more melt ponds on the ice. For these reasons and to make sure we accounted for any biases in sea-ice concentrations, we decided to use both AMSR2 and CISDA products as upper and lower estimates, respectively, when calculating the sea–air CO₂ fluxes in Hudson Bay.

Figure 2 shows the spatial variations of sea-ice concentration, salinity, temperature, and pCO₂ throughout the study area. At the start of our cruise, we observed a large polynya on the northwest side of Hudson Bay, with higher ice concentrations (>50%) encountered mainly in the south of Hudson Strait and the central and eastern side of the bay (Figures 2a and S2). The northwest polynya (hereafter NW polynya) in Hudson Bay has been observed occasionally in winter and spring (e.g., Gough et al., 2004; Landy et al., 2017), resulting from strong offshore westerly winds that advect the existing ice cover eastward and enhance local ice production (Saucier et al., 2004). However, sea-ice melt was actively underway throughout our study period, especially close to the ice edge and coastal areas. The rapid change in the sea-ice cover is evident from maps of sea-ice concentrations throughout the bay at the start, middle, and end of our cruise (Figure S2). At the end of our cruise, sea ice was mainly limited to the eastern side of the bay (Figure S2).

The spatial patterns of salinity and temperature were associated with the spatial patterns of sea-ice concentrations (Figure 2). Salinity ranged from minima of <20 close to the sea-ice edges and the river-dominated coastal region to a high value of 33 in offshore ice-free regions (Figure 2b and Table 1). We observed higher surface seawater temperatures (>8 °C) close to the Churchill and Nelson rivers, and a low value of surface temperature (about –1.7 °C) close to the sea-ice edge (Figure 2c and Table 1). Typically, colder (<1 °C) and less saline waters (salinity between 10 and 25) were encountered in the southern and eastern side of the bay near the ice edge, while relatively warmer (>4 °C) more saline (salinity from 15 to 29) waters were observed in the river-dominated coastal zone (Figure 2a and c). We also observed high salinity (>30) associated with cold water (<0 °C) in Hudson Strait and the southwest side of Baffin Island. The observed variability in salinity matches the general cyclonic (counterclockwise) circulation in Hudson Bay (Prinsenberg, 1986; Ingram and Prinsenberg, 1998), as the marine water entering Hudson Bay from Hudson Strait and Foxe Basin is usually diluted by sea-ice melt and river runoff before exiting the bay to the Labrador Sea via Hudson Strait.

We calculated the area-weighted average and standard deviation of these observations based on a 1 km × 1 km grid to minimize the effect of sampling biases on our calculation of the bay-wide flux resulting from long occupations of some station locations. This binning resulting in decreasing the number of ship observations on the waterside from about 19,500 to 7,817 and on the atmospheric side from about 17,500 to 5,580. Hereafter, all the statistics reported in this study are area-weighted and based on n = 7,817 for the sea and n = 5,580 for the atmosphere. A summary of the range, average, and variability of the observed underway measurements is available in Table 1. The pCO₂ ranged between 125 and 650 μatm (Figure 2d, Table 1), with the lowest measurements observed close to the ice edge on the east side of the bay, the NW side of the bay, and in Hudson Strait (125–280 μatm), while the highest pCO₂ values were mainly observed in ice-covered waters (380–550 μatm) and close to the Churchill and Nelson Estuaries (390–465 μatm). We calculated the average of pCO₂ measurements across Hudson Bay to be 317 ± 61 μatm, which is significantly lower (p = 0) than the average of atmospheric pCO₂ (408 ± 3 μatm), indicating that overall Hudson Bay was undersaturated in pCO₂ relative to the atmosphere.

The measured pCO₂ showed no clear evidence of temporal changes during the study period, although it was clearly correlated in various locations with variations in salinity and temperature (Figure 3a, b). This result suggests that the observed variability was linked to the ship’s spatial proximity to sea-ice or riverine input rather than the timing of sampling. However, we observed a relatively narrow range of salinity and temperature values (cooler, more saline water) during the first 2 weeks (May 30 to June 15), with higher variability and frequency of warmer, less saline water thereafter. This pattern is consistent with the heavy sampling in the NW polynya at the beginning of the cruise and increasing freshwater abundance as the ship entered the southern coastal domain where most of the runoff enters the bay. We also observed an increase in pCO₂ along the SW coast, associated with a decrease in salinity, an increase in temperature, and low sea-ice concentrations. On the east coast, we observed a decrease in pCO₂ values associated with low temperature and salinity. Low pCO₂ values were also observed in Hudson Strait associated with relatively higher temperatures and salinity. Additionally, the pCO₂ values showed some correlation with variability in sea-ice concentrations (Figure 3c). For instance, we observed that the high spikes of pCO₂ were usually associated with a higher sea-ice concentration (>80%), whereas the lowest pCO₂ observations were mainly associated with lower sea-ice concentration (<40%; Figure 3c).

The influence of seawater temperature on pCO₂ is governed by the thermodynamic equilibrium of the inorganic carbon system (Takahashi et al., 1993). For example, if pCO₂ is only controlled by temperature change with no other physical or biogeochemical processes playing a role, the pCO₂ in water will double for every 16 °C increase in temperature (Takahashi et al., 1993). The weighted-area average of seawater temperature observed during our sampling period was 0.8 °C ± 1.7 °C, with values ranging from –1.8 °C to 11.1 °C (Table 1). In order to assess the relative importance of the competing effects of thermal (i.e., temperature) and nonthermal (i.e., biological activity and/or physical water mixing) effects on changes in pCO₂ values, we calculated the temperature-normalized pCO₂ (pCO_2nonthermal) and the effects of temperature on the pCO₂ (pCO_2thermal) using the decomposition used by Takahashi et al. (2002). More details about the thermal and nonthermal pCO₂ calculations are described in Text S1.

In Figure 4, we show the temporal variability of pCO_2nonthermal and pCO_2thermal. The high temporal variability of pCO_2nonthermal compared to pCO_2thermal in the first 3 weeks of our sampling reflects the relatively stable temperatures and significant contribution of nonthermal processes (i.e., biological activity, water mixing, and sea–air CO₂ exchange) to the observed variations in pCO₂. For example, the high pCO_2nonthermal (>450 μatm) in the NW region could be attributed to the inflow of pCO₂-rich water from rivers, mixing with subsurface waters, and/or respiratory CO₂ production.

Starting at about the 4th week of our study period (June 27), during the shipping transit along the SW coast of the bay, pCO_2thermal became more variable (Figure 4). This increased variability is not surprising as we measured higher water temperatures during this time (Figure 3b), which were mostly linked to an increase in the air temperature with time (Figure S3a) and/or the presence of high-temperature riverine waters from Nelson and Churchill Rivers (Figure 2c). The ratio of thermal to nonthermal effects on pCO₂ (ΔpCO_2ratio) shown in Figure 4 illustrates that the surface water pCO₂ in Hudson Bay and Hudson Strait was controlled mainly by nonthermal processes at the start of the sampling period, with the thermal effects more dominant at the end of our study.

As pCO_2nonthermal varies substantially through the entire study, we show the spatial variability of pCO_2nonthermal in Figure 5a. The average of pCO_2nonthermal is 319 ± 66 μatm, with a somewhat smaller range (124–548 μatm) than the observed pCO₂ values (Table 1). The largest differences between pCO_2nonthermal and observed pCO₂ were close to river mouths. For example, the differences between pCO_2nonthermal values and the observed pCO₂ measurements in the Nelson and Churchill Estuaries on average were about 55 μatm (Figures 5a and 2d), indicating that mixing with relatively warm water introduced by Nelson and Churchill Rivers was contributing to the observed high pCO₂.

The distributions of oxygen saturation, Chl a, and FDOM (Figure 5b-d) indicated that for the most part, biological processes and terrestrial organic matter were not a major driver of pCO₂ distribution in Hudson Bay at the time of sampling. This interpretation is supported by the low Pearson’s correlation coefficients (r < .2; p < .05) among pCO₂, Chl a, and FDOM across the study area (see Section 4.3). However, there were some notable, localized exceptions. For example, the FDOM signal was relatively high close to the coast, especially in the SW region of Hudson Bay (Figures 5b and S5), associated with high terrestrial input in this region (Guéguen et al., 2011), while low FDOM signals (<700 counts) were observed in the offshore ice-free water and the NW polynya, likely due to photobleaching in the clearer waters (Dainard et al., 2015). The highest FDOM values (up to 9,000 counts) in our study area were measured east of the Nelson Estuary close to sediment-laden sea ice. Similar to our observations, Granskog et al. (2007) and Guéguen et al. (2011) found high chromophoric dissolved organic matter and FDOM values near the coast in Hudson Bay due to terrestrial input by rivers. Overall, we did not observe any consistent relationship between FDOM and pCO₂ in our study area, indicating that FDOM played a minor role in pCO₂ variability within Hudson Bay during the study period.

The average of Chl a during our sampling time was 1.1 ± 1.2 µg L^–1, with minimum and maximum observations of 0 µg L^–1 and approximately 10 µg L^–1, respectively (Table 1). O₂ saturation averaged 111 ± 9% with a range varying from 90% to 199% across the study area (Table 1). Typically, the response of primary production to ice melt during the spring varied both spatially and temporally as a result of local differences in nutrient supply (e.g., Tremblay et al., 2015). We measured relatively high surface Chl a concentration (>4 µg L^–1; Figure 5d) in Hudson Strait (especially in our transit at the end of the cruise) and in many locations in Hudson Bay (particularly close to the Nelson and Churchill Rivers and in the NW Polynya), associated with oxygen supersaturation (>120%; Figure 5c) and depleted nitrate (NO₃) concentrations (<1 mmol m^–3; Figure S4), but variable pCO₂ values (Figure 2d). Because dissolved O₂ equilibrates much more quickly with the atmosphere than dissolved CO₂, we cannot use O₂ saturation to rule out primary production. However, we can use the observed O₂ supersaturation as evidence of primary production. For example, the strong correlation among low pCO₂ (< 350 μatm), high Chl a, and O₂ supersaturation measured in the NW polynya and Hudson Strait, where the measurements were made before significant gas exchange (i.e., the sea ice was actively melting; Figure 2a), suggest that the O₂ supersaturation is related to ice-associated blooms and that primary production played a significant role in lowering the pCO₂ values in these locations.

By contrast, in the Nelson Estuary, we observed a high Chl a concentration associated with high pCO₂ values (>440 μatm), suggesting that the relatively high productivity in this area was not large enough to overcome other processes acting to increase pCO₂. This is not surprising as the Nelson River has the highest annual discharge to Hudson Bay (Déry et al., 2011), with large loads of organic compounds and nutrients (e.g., Granskog et al., 2007), which potentially can sustain both the high biological production and the respiration that generates high pCO₂ values. Indeed, Demarty et al. (2009) observed the annual average of pCO₂ in the Kettle generating station reservoir located on the Nelson River (about 150 km upstream of Hudson Bay) to vary between 550 and 650 μatm, showing the highly supersaturated nature of this river that might not be completely offset with high biological productivity. Furthermore, we observed variable Chl a concentrations close to the ice edge and under the ice, with values varying from less than 1–4 µg L^–1, and associated with O₂ saturation varying from saturated to supersaturated values.

Similar to our high Chl a observations near the coast, Roff and Legendre (1986) observed higher primary production in these waters compared to offshore waters in Hudson Bay and attributed differences to nutrient replenishment of the upper productive layer through coastal resuspension processes (Griffiths et al., 1981), internal waves (Ingram et al., 1989), or freshwater plume dynamics (Monti et al., 1996). As well, Kuzyk et al. (2010) found that river runoff supports new primary production around the margin of the bay, where river flow is constrained. In contrast, the low surface Chl a (<1 µg L^–1) in many locations in Hudson Bay (Figure 5d), suggests either low productivity or we sampled either before or after the spring bloom. Recently, Matthes et al. (2020) observed that 30% of the annual primary production occurs during the ice melt season (May to June) with stratification likely playing a significant role in delaying the bloom onset into the summer months. That said, our findings of low biological production in many areas of Hudson Bay are consistent with several previous studies that described Hudson Bay as an oligotrophic basin characterized by relatively weak biological production due to a short ice-free season and strong stratification (e.g., Legendre et al., 1996; Ferland et al., 2011; Sibert et al., 2011; Lapoussiere et al., 2013).

Although we cannot quantify whether deep primary production affects our surface pCO₂ measurements, it seems unlikely to have an impact given the depth of the deep Chl a maximum (which varied from 20 to 40 m based on similarities in CTD profiles across the study area) and its position beneath the strong and persistent mixed layer that exists in the region during spring (Ahmed et al., 2020). Similarly, Ferland et al. (2011) observed deep Chl a maxima at 40–60 m during summer and noted a negative correlation between primary production and stratification strength in the upper water column across Hudson Bay. We therefore conclude that biological production had a limited effect on the observed pCO₂ values during our sampling period, except in the NW polynya and Hudson Strait where primary production may have contributed to the observed pCO₂ undersaturation.

To study the impact of freshwater from sea-ice melt and river runoff on pCO₂ variability in Hudson Bay, we plotted temperature versus salinity visualized with pCO₂ values throughout the study area (Figure 6). Past studies have shown that sea-ice meltwater is initially undersaturated in pCO₂ (e.g., Geilfus et al., 2015) with low seawater temperature and salinity, whereas Arctic river waters tend to be relatively warmer and often supersaturated in pCO₂ as a result of degrading terrestrial organic carbon (e.g., Semiletov et al., 2011; Denfeld et al., 2013).

We measured low surface salinity seawater (i.e., <28) associated with colder water (<1 °C) and low pCO₂ (<350 μatm) in many locations close to the ice edge in Hudson Bay (Figures 6 and 2b and d), indicating that sea-ice melt contributed to pCO₂ undersaturation in these regions. We also measured high temperature (>5 °C) associated with relatively low salinity observations (<29) and high pCO₂ (>420 μatm) in some of the coastal regions in Hudson Bay (especially near the Nelson and Churchill Estuaries; Figures 6 and 2b and d), suggesting that river water contributed to pCO₂ supersaturation in these regions. Similarly, Capelle et al. (2020) noted that the inputs of terrestrial organic matter into Hudson Bay were a significant driver of CO₂ evasion in the coastal surface waters. Interestingly, the observations from the Churchill Estuary were warmer (0.5–11.2 °C) and had relatively lower pCO₂ (averaged 350 μatm with a range of 230–480 μatm) compared to the temperature (0.9–8.3 °C) and pCO₂ (averaged 416 μatm with a range of 370–470 μatm) observations recorded in the Nelson Estuary (Figure 6). These differences between both estuaries could be attributed to differences in watershed chemistry, associated with the different physiographic regions, which likely influences the biology and chemistry of the two rivers (Rosenberg et al., 2005); river discharge (Déry et al., 2016), which is expected to influence carbon loading and marine biogeochemistry; and/or atmospheric conditions (e.g., air temperature and solar radiation) prior to sampling. For example, water and air temperatures in the Churchill Estuary were about 4 °C and 9 °C warmer, respectively, than in the Nelson Estuary during sampling (Figures 6 and S3a). Also, the Churchill River has higher concentrations of dissolved organic carbon than the Nelson River (Mundy et al., 2010), which is consistent with the Churchill watershed having a lower carbonate rock content and a higher permafrost abundance than the Nelson watershed.

Results from past studies in Hudson Bay are consistent with our observations, with low dissolved inorganic carbon associated with sea-ice melt (Burt et al., 2016), which typically results in low pCO₂ values (e.g., Ahmed et al., 2020). Additionally, our observations of high pCO₂ values close to river-dominated nearshore waters match the measurements of Else et al. (2008a). Overall, our results suggest that Hudson Bay was undersaturated in pCO₂ relative to the atmosphere during our study time as a result of the dominance of low-pCO₂ ice-melt waters, with few observations of supersaturated pCO₂.

To more closely examine the effect of ice melt and riverine input on pCO₂ in Hudson Bay, we visualized the ship transects from ice-covered waters to coastal areas in the NW polynya and the Nelson and Churchill Estuaries (Figure 7). This plot shows typically high pCO₂ (>410 μatm) in the coastal domain and low pCO₂ (<330 μatm) observed mostly in the offshore waters. Figure 7a confirms that the Nelson and Churchill Rivers were responsible for the supersaturated pCO₂ during our sampling time, but their impact diminishes rapidly offshore toward the ice-melt-influenced waters. Similar to the transects in Nelson and Churchill Estuaries, off Chesterfield Inlet, we observed relatively high pCO₂ values that decreased rapidly as we moved offshore across the NW polynya (Figure 7b). The undersaturated pCO₂ observations in offshore waters of the NW polynya are probably a result of ice melt and/or biological productivity promoted by the delivery of upwelled and/or riverine nutrients. Many of these processes have been observed in past studies of Arctic polynyas (e.g., Yager et al., 1995; Else et al., 2012).

We often observed high pCO₂ (>415 μatm) associated with relatively high salinity (>31), cold (< –0.5 °C) waters when the ship was navigating through the heavy ice cover (Figure 7). The observed high pCO₂ measurements were also associated with high nitrate concentration (Figure S4), suggesting that the under-ice phytoplankton bloom had not started yet when we were sampling at these locations. These high under-ice pCO₂ observations are likely remnants of biogeochemical and physical processes during the winter season. For example, community respiration has been observed to exceed photosynthesis during the winter as a result of low incoming solar radiation attenuated by sea ice and the overlying snowpack (Shadwick et al., 2011; Else et al., 2012). Another potential source of the observed high pCO₂ might be the expulsion of CO₂-enriched brines into the water column during sea-ice formation and the retention of alkalinity by calcium carbonate crystals within the sea ice (Rysgaard et al., 2007; Dieckmann et al., 2008). Overall, Figure 7 shows both CO₂ supersaturation and undersaturation in the marginal ice zone, suggesting very high spatial variability in the distributions of sea-ice melt at these locations.

The observed relatively high salinity (>32) and low temperature (<0 °C) on the east side of Southampton Island in Hudson Bay (Figure 7b) and the southwest side of Baffin Island indicate that upwelling contributed to high pCO₂ values in these regions (>420 μatm). Upwelling is a wind-driven process that brings the deep, cold, and usually high pCO₂ water to the surface, often resulting in oceanic CO₂ outgassing. Several past studies (e.g., Prinsenberg, 1986; Kuzyk et al., 2010) have observed upwelling events in the NW region of Hudson Bay during the summer and fall seasons, suggesting that upwelling in this area may be a long-lasting event rather than an episodic phenomenon. Therefore, although the pCO₂ supersaturation we observed in association with upwelling during our sampling was limited in spatial extent, it may persist through much of the open water season and possibly occur in other locations within Hudson Bay.

To explain statistically the key variables affecting pCO₂ in Hudson Bay, we performed a linear regression analysis between pCO₂ as the dependent variable and salinity, temperature, Chl a, and FDOM as independent or explanatory variables. We limited our spatial analysis to the ice-free regions (i.e., ice concentrations less than 10%) in Hudson Bay where the majority of our field measurements were collected. As well, we excluded the pCO₂ values impacted by the upwelling event in the NW region to ensure that our regression relationship is not biased by this spatially limited phenomenon.

We found that salinity described the major variability in pCO₂ (Pearson’s product-moment: r = –.71; p < .05), followed by temperature (r = .65; p < .05; Figure 8 and Table 2). Although Chl a and FDOM measurements showed a significant linear relationship with pCO₂, the Pearson’s correlation coefficients were low (Table 2), showing that their impact on the pCO₂ variability is minor but statistically significant.

The regression relationship between salinity and pCO₂, statistically significant at a 95% confidence level, was

A root mean square error (RMSE) of 36 μatm and relatively strong coefficient of determination (R² = 0.51) indicate that this algorithm is explaining more than half of the variability in pCO₂, which is a good result, given that this algorithm operates over a range of about 270 μatm (i.e., pCO₂ ranged from 198 to 468 μatm in the ice-free regions). To examine the significance of other variables on pCO₂, we added temperature, FDOM, and Chl a, sequentially, and selected the best model based on the lowest RMSE, Akaike information criterion, and Bayesian information criterion. The resulting multivariate relationship combined salinity, temperature, FDOM, and Chl a, with RMSE equal to 35 μatm and R² of .55:

To estimate CO₂ fluxes, we used a bulk flux equation:

where F_CO2 is the CO₂ flux (mmol CO₂ m^–2 day^–1), α is the solubility of CO₂ in seawater (mol m^–3 atm^–1) as a function of temperature and salinity (Weiss, 1974), and ΔpCO₂ is the difference between the measured surface seawater pCO₂ and the average atmospheric pCO₂ (in μatm). We calculated sea–air CO₂ flux scaled by sea-ice concentration (C_i, reported as a percentage) following the studies of Butterworth and Miller (2016) and Prytherch et al. (2017). The gas transfer velocity k (cm h^–1) was calculated according to Wanninkhof (2014):

where $U10n2$ is the square of wind speed (m s^–1) corrected to a height of 10 m, and Sc is the Schmidt number (unitless). Negative sea–air CO₂ fluxes represent atmospheric CO₂ uptake (oceanic sink), whereas positive fluxes represent oceanic CO₂ outgassing to the atmosphere (oceanic source).

As the shipboard wind speed data were not available after July 2, we estimated the CO₂ fluxes in Hudson Bay using both shipboard wind speed data and NARR wind data (Figure 9). The NARR wind speeds were usually underestimated compared to the in situ measurements (Figure 9a), with the highest wind speeds (>11 m s^–1) during our study time observed in the NW polynya and the Nelson Estuary (Figure 10a). The shipboard wind speeds ranged up to 15 m s^–1 with an average of 5.8 ± 2.9 m s^–1 (Table S1 and Figure 9a), whereas the NARR wind data showed a similar range with a lower spatial average of 4.4 ± 2.9 m s^–1 (Table S1 and Figure 10a). Consequently, CO₂ fluxes calculated using the shipboard wind data had a higher magnitude than those calculated using the NARR wind data (Figure 9b). The average CO₂ fluxes using the NARR and the ship wind data were –5.8 and –7.8 mmol CO₂ m^–2 day^–1, respectively (Table S1). Therefore, we consider the CO₂ fluxes calculated based on the NARR wind data as lower-bound estimates of the oceanic CO₂ sink in Hudson Bay.

The calculated ΔpCO₂ varied from about –283 μatm to 242 μatm, with an average of –91.5 ± 60.9 μatm (i.e., the surface waters were generally undersaturated, with respect to the atmosphere). During our study period in Hudson Bay, we observed supersaturated values of ΔpCO₂ (>10 μatm) mainly close to the Nelson and Churchill Estuaries and along the southeast of Southampton Island (Figure S3b). Whereas, we observed strongly undersaturated values of ΔpCO₂ (< –100 μatm) mostly in Hudson Strait and along the east side of the bay and in the NW polynya.

The high variation in the CO₂ fluxes shown in Figure 10b is attributed mostly to the variability in seawater pCO₂ resulting from the regional and dynamic variability of ice melt and riverine input (as discussed in Section 4.2). However, the high CO₂ uptake observed in the NW polynya was associated with high wind speeds (Figure 10a) over ice-free waters (Figure 2a). To examine this spatial variability of sea–air CO₂ fluxes more closely, we visualized the ship transects in the Nelson and Churchill Estuaries and the NW polynya (Figure S6). Overall, the strong CO₂ outgassing observed in these locations (Figure S6) was correlated with the supersaturated pCO₂ (Figure 2d), resulting in positive ΔpCO₂ values (Figure S3b) calculated for these ice-free areas.

Using all of the observations collected during our study, we found that Hudson Bay overall acted as a moderate to weak net sink for atmospheric CO₂ with an average CO₂ uptake of 5.1 (±9.3) mmol CO₂ m^–2 day^–1 (Table 1). We estimated the uncertainty in the CO₂ flux to be 38% (±1.9) mmol CO₂ m^–2 day^–1 based on the cumulative uncertainties derived from the field observations, the gas exchange parameterization, and sea-ice concentration (see Text S2 for additional details on the uncertainty calculations). The moderate average sink of CO₂ that we calculated is a result of the dominance of uptake versus outgassing over the study area. Using the 1-km spatial bin analysis, 90% of the sampled area acted as a sink for CO₂, while only 10% of the area acted as a source of CO₂. Varying the source of sea-ice data did not seem to have a significant effect on these results; using the CISDA ice data in Equation 3, we calculated the average CO₂ uptake to be 4.8 (±9.0) mmol CO₂ m^–2 day^–1.

While reporting the instantaneous CO₂ fluxes based on ship observations, understanding the hourly to daily trends and patterns in spring and early summer seasons is important (Figures 9b and 10b); otherwise, estimates might not be representative of the overall source-sink status of the bay. The problem with extrapolating ship observations to calculate CO₂ fluxes over space and time is that the processes influencing seawater pCO₂ are variable, hence biasing estimates to the sampling times and locations. For example, the biological activity and the magnitude of river inputs almost certainly change over time and space beyond our sampling survey. Although this spatial–temporal variability should be kept in mind, we expect that our resulting errors will be relatively low, given that Hudson Bay is known as an oligotrophic basin with low biological activity (e.g., Ferland et al., 2011) and that river discharge to Hudson Bay has relatively consistent seasonal hydrographs (Déry et al., 2016; 2018). Assuming relatively consistent or slow change of freshwater input and biological production, we can extrapolate the average CO₂ fluxes over space and time to estimate the sea–air CO₂ fluxes in the open water areas in Hudson Bay. This approach of extrapolating ship CO₂ fluxes has been adopted by many studies in the Arctic regions to calculate the regional and annual CO₂ fluxes (e.g., Bates and Mathis, 2009, and references therein).

Considering the sea-ice coverage during our sampling time, we estimated the area of open water in each month separately (i.e., May, June, and July in 2018) using the CISDA ice data. Using the average sink of –5.1 mmol CO₂ m^–2 day^–1 and the calculated open water areas, we estimated the total CO₂ fluxes for Hudson Bay in May, June, and July to be –0.3, –1.1, and –1.9 Tg C (10¹²g C), respectively, for a total net sea–air CO₂ flux of –3.3 (± 1.2) Tg C.

Compared to other Arctic shelf regions, Hudson Bay during the spring and early summer seasons acts as a relatively moderate to weak CO₂ sink (Table 3). Our calculated average uptake rates are similar to rates estimated for the Laptev and East Siberian Seas (Pipko et al., 2017), which are shelf seas similarly affected by strong river runoff. The calculated rates are also close to estimates for the Canadian Arctic Archipelago, which like Hudson Bay overall experiences low primary production (Ahmed and Else, 2019). In contrast, the Chukchi Sea (which has high rates of primary and net community production during the ice-free period) experiences much stronger CO₂ uptake (Bates, 2006). The high CO₂ uptake observed in the Kara sea (a shelf sea that receives a lot of riverine input) compared to Hudson Bay was more likely due to sampling close to sea-ice meltwater under high wind speeds (Pipko et al., 2011).

Our estimates of CO₂ flux through May, June, and July complement past research in Hudson Bay that calculated CO₂ fluxes in late summer and autumn (Else et al., 2008a, 2008b). Else et al. (2008a) found that warm sea surface temperature in late summer (August) causes weak outgassing, which eventually reverts to a weak CO₂ sink as the surface waters cool through the fall. Combined with the moderate uptake we observed during spring and early summer, a fairly complete picture emerges that is consistent with the typical seasonal CO₂ cycle for high-latitude shelf seas described by Bates and Mathis (2009, and references therein). However, the temporal separation of the two studies (13 years) limits our confidence in any calculation of the total sea–air CO₂ flux throughout Hudson Bay over the entire open water season solely from the observational results.

Instead, Else et al. (2008b) found that the variability in seawater temperature is a good predictor for the pCO₂ variability during the summer and fall seasons in Hudson Bay. Therefore, to compute the total CO₂ fluxes during the entire 2018 open water season in Hudson Bay, we can use satellite-derived sea surface temperature and the average pCO₂ during our study (approximately 317 μatm), corrected for temperature following the approach of Takahashi et al. (2002) to estimate pCO₂ during the 2018 summer and fall seasons. To estimate the average seawater temperature from August to October in 2018, we used the monthly Level-3 MODIS Aqua (https://oceancolor.gsfc.nasa.gov/), clipped the images to the study area, and resampled them to 1-km spatial resolution. The calculated average water temperatures for August, September, and October months were 5.6 °C, 4.3 °C, and 1.9 °C, respectively. The resulting estimated pCO₂ values for August, September, and October were about 402, 380, and 344 μatm, respectively. We calculated k (gas transfer velocity) using the NARR wind data and Equation 4 and estimated the sea–air CO₂ fluxes using Equation 3. The estimated CO₂ fluxes during the 2018 August, September, and October were –0.6, –1.2, and –2.2 Tg C, respectively. Adding these values to those we estimated for the spring and early summer seasons, we obtain a rough estimate of the total open water CO₂ uptake in Hudson Bay of –7.2 Tg C. While a complete annual budget should account for the sea–air CO₂ exchange during winter, through-ice fluxes are orders of magnitude smaller than open water fluxes and are not expected to have a large impact on the annual budget (Butterworth and Else, 2018). At an areal extent of about 1 × 10⁶ km², Hudson Bay and Hudson Strait make up approximately 8% of the Arctic Ocean and its associated shelf seas (approximately 13 × 10⁶ km²). Here, we estimated a net sea–air CO₂ flux of –7.2 Tg C during the open water season, suggesting that Hudson Bay and Hudson Strait are likely to account for 4–6% of the Arctic Ocean carbon sink, which previous estimates have shown is on the order of 118–180 Tg C year^–1 (e.g., Manizza et al., 2019, and references therein). As a result, Hudson Bay is considered a slightly weaker sink than might be expected based on its surface area likely due to the influence of high-pCO₂ river water along the coastal regions. Yet, this contribution remains significant at the pan-Arctic scale and underscores the importance of including Hudson Bay in carbon budgets for the Arctic Ocean.

This study presents the first spring and early summer measurements of pCO₂ and estimates of sea–air CO₂ fluxes in Hudson Bay and Hudson Strait and thus fills up a major gap in our understanding of sea–air CO₂ fluxes in the Arctic/Sub-Arctic shelf seas. During our study, sea surface pCO₂ exhibited large spatial variability, with localized pCO₂ supersaturation in areas influenced by rivers and upwelling, as well as under the sea ice, and more generalized pCO₂ undersaturation in ice-melt-influenced waters. Overall, Hudson Bay was mostly undersaturated in pCO₂. Temperature variations, and thereby their effect on pCO₂, were small at the beginning of our study, when water mixing and biology activity were the primary drivers of pCO₂ variability. Along the southern coast, where most of the relatively warm riverine waters enter the bay, the thermal effect on pCO₂ became more important. The contribution of biological production to the variability in pCO₂ was minor but nonnegligible, especially in the NW polynya and Hudson Strait. In the open water areas, pCO₂ showed the highest correlations with salinity and temperature, with minor but statistically significant correlations with Chl a and FDOM. Our observations in the Nelson and Churchill Estuaries highlight the need to better understand variability in the biogeochemical makeup of rivers entering Hudson Bay and their impacts on the receiving marine system. Hudson Bay as a whole acted as a CO₂ sink with an average flux of about –5 mmol CO₂ m^–2 day^–1 (–3.3 Tg C) during the study, with an estimated uncertainty of 38% that corresponds to about ±1.9 mmol CO₂ m^–2 day^–1. The calculated sink in this study is significantly larger than previous estimates for Hudson Bay in late summer and autumn, underlining the significance of properly accounting for seasonal variability in Arctic shelf seas. Nevertheless, on average, the bay was a weaker CO₂ sink than Arctic shelves due to the strong influence of high-pCO₂ river water. We estimated the total CO₂ uptake in Hudson Bay during open water season to be –7.2 Tg C, contributing about 4%–6% of the annual CO₂ flux estimated for the Arctic Ocean. Future work will attempt to use remotely sensed data and machine learning techniques to extend our observations and generate a more comprehensive estimate of regional sea–air CO₂ fluxes in Hudson Bay.

The underway data presented in this article were uploaded to the Surface Ocean CO₂ Atlas (SOCAT) and will be available to the public in SOCAT Version 2021.

The supplemental files for this article can be found as follows:

Figure S1. Sea-ice concentration from daily AMSR2 and weekly CISDA sea-ice products and ship-recorded observations. Sea-ice concentration (%) from daily AMSR2 products (blue), weekly CISDA sea-ice products (black), and in situ observations recorded by the ship officers (red circles). AMSR2 = Advanced Microwave Scanning Radiometer 2; CISDA = Canadian Ice Service Digital Archive.

Figure S2. Regional maps of changes in sea-ice cover over time. Sea-ice concentration (%, color-coded) at approximately the start (May 28), middle (June 18), and end (July 16) of our cruise in Hudson Bay. Sea-ice concentrations were provided by the Canadian Ice Service. Dates are year-month-day.

Figure S3. Air temperature and sea–air CO₂ gradient in Hudson Bay. Spatial variability of (a) NARR air temperature and (b) the measured sea–air gradient of CO₂ partial pressure (ΔpCO₂) during our study time (May 25 to July 13, 2018).

Figure S4. Nitrate and SPAR variability in Hudson Bay. Surface distribution of (a) nitrate (NO₃) concentrations and (b) surface photosynthetic active radiation (SPAR) across Hudson Bay during our study. The dots indicate locations of CTD casts, and the white area represents sea-ice concentration greater than 90%, as of July 9, 2018, based on weekly ice charts provided by the Canadian Ice Service. SPAR = surface photosynthetic active radiation; CTD = conductivity-temperature-depth.

Figure S5. pCO₂ versus FDOM. Scatterplot of pCO₂ versus fluorescent dissolved organic matter (FDOM) in Hudson Bay during our study, with the measurements from the SW coast in red. pCO₂ = CO₂ partial pressure.

Figure S6. CO₂ fluxes in the Nelson Estuary and NW polynya. Spatial variability of sea–air CO₂ fluxes in (a) the Nelson and Churchill Estuaries with the sea-ice cover as of July 2, 2018, and (b) the NW polynya with the sea-ice cover as of June 16, 2018. NW polynya = northwest polynya.

Table S1. Shipboard and NARR match-up measurements of wind speed, air temperature, and sea–air CO₂ fluxes from May 31 to July 2, 2018.^a NARR = National Centers for Environmental Prediction North American Regional Reanalysis.

Text S1. Thermal and nonthermal pCO₂ calculations. pCO₂ = CO₂ partial pressure.

Text S2. Uncertainties of CO₂ flux estimation.

The authors thank the Captain and the crew of the Canadian research icebreaker CCGS Amundsen for their hard work and help during the BaySys cruise in 2018. They thank Rachel Mandryk and Yekaterina (Kate) Yezhova for their help during sampling; Tonya Burgers and Sebastian Luque for logistics and help during ship mobilization. They also thank Dr Elizabeth Jones and the second anonymous reviewer for their valuable comments and suggestions to improve the article. The ship thermosalinograph and conductivity-temperature-depth data used in this study were collected by systems onboard the Amundsen and made available by the Amundsen Science Program.

This study was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) Collaborative Research and Development (CRD) Program, Manitoba Hydro, and Fisheries and Oceans Canada.

The authors have no competing interests to declare. LAM serves as Elementa Associate Editor. She was not involved in the review of this article.

Contributed to the conception and design: MMMA, BGTE, TP.

Contributed to the acquisition of data: MMMA, DWC, BB, TP.

Contributed to the analysis and making the figures: MMMA.

Contributed to the interpretation of results: All authors.

Drafted the article: MMMA.

Revised the article: All authors.

Approved the submitted version for publication: All authors.

How to cite this article: Ahmed, MMM, Else, BGT, Butterworth, B, Capelle, DW, Guéguen, C, Miller, LA, Meilleur, C, Papakyriakou, T. 2021. Widespread surface water pCO₂ undersaturation during ice-melt season in an Arctic continental shelf sea (Hudson Bay, Canada). Elementa: Science of Anthropocene 9(1). DOI: https://doi.org/10.1525/elementa.2020.00130

Domain Editor-in-Chief: Jody W. Deming, University of Washington, Seattle, WA, USA

Associate Editor: Jean-Éric Tremblay, Department of Biology, Université Laval, Quebec, Canada

Knowledge Domain: Ocean Science

Part of an Elementa Special Feature: BaySys