Sea ice melt drives vertical pCO₂ variability modulating air–sea gas exchange

Strong spatial and temporal gradients in salinity, temperature, and carbonate chemistry in Arctic coastal surface waters complicate the estimation of air-sea carbon dioxide (CO₂) exchange, particularly during sea ice breakup. The present study evaluates the applicability of the widely used bulk flux model under such conditions. This approach assumes homogeneous surface conditions and no vertical pCO₂ gradients in the bulk seawater. However, our observations in a stratified Arctic fjord reveal pronounced vertical variability in pCO₂ within the upper water column, including non-linear gradients near the air-sea interface. This results in, widely varying flux estimates depending upon the depth of the pCO₂ measurement used to establish air-sea disequilibrium. Importantly, similarly structured nonlinear pCO₂ profiles were observed across distinct fjord systems and years, suggesting that this vertical heterogeneity may be a characteristic feature of Arctic stratification during sea ice breakup. We recommend incorporating both micrometeorological techniques and high-resolution vertical profiling in Arctic fjords to improve flux estimates of CO₂ in this rapidly changing region.

1 Introduction

High latitude coastal oceans are strong sinks for atmospheric carbon dioxide (CO₂), absorbing more CO₂ per unit area than lower latitude regions (Dai et al., 2022; Roobaert et al., 2019). This strong uptake results from both the high solubility of gases in cold water and the intense biological activity typical of these regions. However, climate change is rapidly transforming this carbon sink. The Arctic is warming more than twice as fast as the global average, and sea ice extent has been shrinking by over 13 % per decade (Perovich et al., 2020). The loss of sea ice increases CO₂ uptake by exposing larger areas of open water for longer periods, which can further stimulate biological productivity (Arrigo and van Dijken, 2015; Bates and Mathis, 2009; Perovich et al., 2020). However, at the same time, melting sea ice freshens the surface layer and strengthens stratification, limiting vertical mixing with deeper water. Freshwater from melting sea ice and terrestrial run-off creates pronounced gradients in physical properties such as salinity and temperature, as well as chemical properties like dissolved inorganic carbon (DIC) and total alkalinity (TA) (e.g. Henson et al., 2025a). As a result, the partial pressure of CO₂ (pCO₂) can vary markedly with depth under melt conditions (Miller et al., 2019).

This vertical variability in pCO₂ poses a challenge for air-sea CO₂ flux estimation. The transfer of gases between the atmosphere and ocean depends on the difference in concentration between the two as well as the efficiency of the transfer process. Therefore, the bulk flux of CO₂ across the air-sea interface is commonly described as the product of the gas transfer velocity, k (m s¹), CO₂ solubility s ($\mathrm{mol}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{atm}}^{-\mathrm{1}}$), and the partial pressure gradient (µatm) across the air-sea interface (Wanninkhof et al., 2009):

$$\begin{array}{}\text{(1)}& F=ks(p{\mathrm{CO}}_{{\mathrm{2}}_{\text{sea}}}-p{\mathrm{CO}}_{{\mathrm{2}}_{\text{air}}})\end{array}$$

While widely applied, this formulation simplifies a complex process influenced by surfactants on the water surface, bubble-mediated gas exchange, and turbulence (Wanninkhof et al., 2009). Furthermore, surface water heterogeneity, driven by sea ice melt and freshwater runoff from land, complicates the physical and chemical processes governing air-sea CO₂ exchange. As a result, simplified parameterizations commonly used in global carbon flux estimates may be inadequate in these settings.

In most studies, pCO₂ is measured several meters below the surface, assuming vertical homogeneity under well-mixed conditions (Jørgensen et al., 2020). However, in stratified waters, where temperature, salinity, and pH can vary with depth, this assumption may lead to substantial errors in flux estimates (Ahmed et al., 2020; Dong et al., 2021; Miller et al., 2019; Watts et al., 2022). In this study, we distinguish between the air-sea interface, the near-surface freshwater lens (approximately the upper 0–0.5 m), and subsurface bulk waters (>1 m depth), as these layers may exhibit substantially different physical and carbonate chemistry conditions during sea ice melt. Previous work in Young Sound demonstrated that snow and sea ice melt can substantially reduce under-ice pCO₂ through dilution and enhanced stratification during the melt season (Verdugo et al., 2025). However, these observations were conducted at a fixed depth below the ice, leaving unresolved how pCO₂ variability evolves vertically during sea ice breakup.

Although Arctic surface waters are often undersaturated with respect to atmospheric CO₂ levels and act as CO₂ sinks (e.g. Burgers et al., 2017; Dai et al., 2022; Henson et al., 2024; Laruelle et al., 2014; Roobaert et al., 2019), such assessments typically rely on sparse data collected from 0.5–5 m depth during limited periods. Dong et al. (2021) illustrate that high latitude fluxes of CO₂ calculated using the bulk method (based on measurements sampled at 6 m depth) differ significantly from those measured using direct eddy covariance in sea ice melt regions.

Gas transfer velocity (k) is often parameterized as a function of wind speed. However, the true driver is mixing in the surface waters, which governs k. Fick's first law of diffusion, which underlies Eq. (1), assumes a linear concentration gradient within the diffusive sublayer (Fig. 1) and steady-state conditions (Garbe et al., 2014). Jørgensen et al. (2020) argued that, due to seawater's high buffer capacity, chemical gradients do not significantly affect CO₂ equilibration, supporting the use of measurements at 3–4 m depth. However, this conclusion relies on the assumption of horizontal and vertical homogeneity and neglects the effects of shallow surface stratification, particularly when alkalinity dilution is involved.

Figure 1 Schematic illustrating the interface between the air and the water in conjunction with pCO₂ concentration gradients. In Eq. (1), the concentration gradient is assumed to occur in the diffusive layer between the air and water, and the concentrations are assumed to be vertically constant in the turbulent layers (Adapted from Liss and Slater, 1974; Wanninkhof et al., 2009).

In Arctic spring, the upper ocean is often strongly stratified due to freshwater input from glacier melt, snowmelt, river runoff, and sea ice meltwater (Ahmed et al., 2020; Granskog et al., 2011; Meire et al., 2017; Miller et al., 2019). These inputs can extend vertical CO₂ gradients beyond the diffusive sublayer, complicating flux estimates during ice breakup and early open-water periods. Several studies have demonstrated strong vertical heterogeneity in pCO₂ in Arctic coastal waters, with implications for air-sea flux calculations (Ahmed et al., 2020; Dong et al., 2021; Miller et al., 2019).

Surface freshening from ice melt and runoff strongly influences carbonate chemistry in Arctic coastal waters, which can either suppress or enhance oceanic CO₂ uptake. For example, Burgers et al. (2017) reported large horizontal variability in surface pCO₂ (144–364 µatm) linked to riverine input in the Eastern Canadian Arctic. Similarly, Sejr et al. (2011) observed strong surface pCO₂ gradients associated with salinity and temperature in Young Sound, and later documented a long-term decline in surface salinity (Sejr et al., 2017). Freshwater-induced stratification has also been shown to create vertical gradients in pCO₂ and pH with important implications for flux calculations (Miller et al., 2019). Finally, Bates et al. (2014) demonstrated that sea ice meltwater and melt ponds exhibit extreme variability in pCO₂ (<10 to >1500 µatm) and pH (6.1 to >10.8), highlighting the complex chemical landscape of ice-influenced waters. Together, these studies underscore the high spatial and temporal variability of carbonate chemistry in freshened waters across the Arctic.

To project future CO₂ uptake or outgassing in the Arctic, we must better understand the physical and chemical drivers of near-surface carbonate variability. In this study, we investigate the vertical and temporal variations in pCO₂ in a stratified Arctic fjord during sea ice breakup. By examining water-column pCO₂ profiles across the transition from ice-covered to open water, we evaluate usage of the bulk flux model under Arctic seasonal transitions.

2 Study Site and Measurement Methods

2.1 Study Site

This study was conducted in Young Sound, a high Arctic fjord system located near the Daneborg Research Station in Northeast Greenland (Fig. 2). The fjord system comprises the Tyrolerfjord (inner fjord) and Young Sound (outer fjord), extending approximately 90 km from Tyroler River to the Greenland Sea. A sill at about 45 m depth separates Young Sound from the open ocean. Young Sound is 2 to 7 km wide, with an average depth of 100 m (maximum 350 m), and a total surface area of ∼ 390 km². Tidal amplitudes range from 0.8 to 1.5 m, with mean current velocities of approximately 2 cm s⁻¹ (Rysgaard et al., 2003). Freshwater inputs are primarily derived from Greenland Ice Sheet runoff, local glaciers, precipitation, and snowmelt from adjacent ice-free terrain. The drainage basin of the Tyrolerfjord-Young Sound system spans 2846 km², of which 33 % is glaciated.

Figure 2 Map of Greenland and the sampling area at the coast of Young Sound in Northeast Greenland. The red circle indicates the location of the micrometeorological measurement tower at research station Daneborg while the marine sampling site (Standard Station in the Greenland Ecosystem Monitoring program) is indicated as a blue circle (74.310, −20.300). Three Copernicus Sentinel true-color images of the fjord on 12, 22, and 31 July illustrate the transition between sea ice cover and open water.

Sampling was conducted from 12 to 31 July 2017. Sampling occurred during and immediately after a period of sea ice breakup. On 15 July, ice coverage was approximately 30 %, decreasing to less than 10 % by 16 July. Water sampling was conducted both from an inflatable boat and via sea ice leads, all in close proximity to the Greenland Ecosystem Monitoring (GEM) program's standard station (Fig. 2).

2.2 pCO₂ Measurements Using the HydroC Sensor

Surface water pCO₂ was measured with a CONTROS^® HydroC CO₂ sensor, which utilizes a membrane equilibrator coupled with a non-dispersive infrared detector (Henson et al., 2025b). The instrument is equipped with a built-in water pump that provides flow rate of 35 mL s⁻¹ across the membrane. At each sampling depth, the sensor was allowed to equilibrate for 10 to 20 min, and values were recorded once stable for at least 2 min. The sensor operates over a range of 200–1000 µatm and temperatures of −2 to 35 °C. Annual calibration has been conducted using a certified 400 ± 2 % ppm CO₂ gas that was traceable to WMO standards. The sensor showed remarkable stability (397–401 ppm), supporting a measurement uncertainty of ±2 µatm.

2.3 pCO₂ Estimation from TA and DIC

In addition to direct measurements, pCO₂ was calculated from total alkalinity (TA) and dissolved inorganic carbon (DIC) using the Seacarb package (Gattuso et al., 2021) in R. Due to the low salinity and cold temperatures characteristic of Arctic coastal waters, no universally accepted set of equilibrium constants (K1 and K2) exists. For consistency with previous studies in the region (Henson et al., 2023), we used the refitted constants from Lueker et al. (2000). The selection of equilibrium constants introduces assumptions regarding seawater composition. Raimondi et al. (2019) showed that different constants can lead to discrepancies between measured and calculated pCO₂ values, ranging from −3.1 to −35.8 µatm, with Lueker et al. (2000) demonstrating the best internal consistency under polar conditions. Still, Sulpis et al. (2020) found that the calculation of pCO₂ from DIC and TA can lead to uncertainty up to 15 % under cold conditions, which is far greater than when pCO₂ is measured directly.

2.4 Sea Ice TA and DIC Sampling

TA and DIC in sea ice were assessed using three ice cores. Each core was sectioned into 5–10 cm segments and sealed in gas-tight NEN/PE bags with sampling valves (Hansen et al., 2000). Samples were transported in thermally insulated boxes to a nearby field laboratory. Cold (1 °C) deionized water of known mass and carbonate composition (10–30 mL) was added to each bag, which was then resealed after removing air and weighted.

The samples were melted in the dark over ∼ 48 h. Meltwater was transferred to 12 mL Exetainer vials (Labco, UK) pre-dosed with 20 µL of saturated HgCl₂ solution (5 % $w/v$) to prevent microbial alteration. DIC was measured by on Apollo SciTech^®'s AS-C3 analyzer while TA was determined via potentiometric titration on an Apollo SciTech AS-ALK2 total alkalinity titrator (Haraldsson et al., 1997).

2.5 Physical Parameters

Vertical profiles of conductivity, temperature, and depth (CTD) were obtained using a Seabird^® SBE19plus CTD. On 16 July 2017, additional surface conductivity measurements were taken using a Thermo Orion-Star^® instrument with an Orion 013610MD conductivity cell. Surface water temperatures were independently measured with a Testo^® thermometer.

2.6 Historical Data

For contextual comparison, pCO₂ time series data (Greenland Ecosystem Monitoring, 2020) from the Greenland Ecosystem Monitoring program are also included in the analysis. pCO₂ data from 2007–2023 was measured using the same HydroC CO₂ sensor in August each year.

2.7 Eddy Covariance

Sensible and latent heat fluxes were estimated using micrometeorological instrumentation mounted on a 3 m mast positioned approximately 0.5 m from the waterline. Three-dimensional wind vectors were recorded using a METEK^® uSonic-Scientific sonic anemometer. To enhance reliability, we applied complementary analysis techniques for flux estimation: (1) the standard eddy covariance (EC) method using EddyPro software (Version 7.0.6, LI-COR Inc., 2019); (2) the ogive optimization method (OGM) (Sievers et al., 2015). Among these, the OGM was deemed most robust due to its ability identify and filter out low-frequency noise, sensor dampening, and large-scale turbulent motions that can bias flux measurements. These issues often introduce large relative bias associated with flux measurement over Arctic marine surfaces (Sievers et al., 2015). OGM's superior ability to isolate relevant turbulent scales and reduce contamination from mesoscale variability is based on the accumulation and modelling of each cospectra over each 20 min averaging period (Figs. S1 and S2 in the Supplement). Uncertainty in sensible and latent heat fluxes was estimated directly from the OGM procedure. The reported values correspond to the standard error associated with the fitted ogive tail and reflect random uncertainty in flux integration.

3 Data and results

Vertical profiles of surface water pCO₂ were measured using the CONTROS^® HydroC CO₂ sensor across three distinct periods in July 2017 (Fig. 3a–c). Each observational period corresponded to different sea ice conditions: before, during and after sea ice breakup (Fig. 2). These high-resolution profiles revealed substantial vertical variability within the upper 2 to 3 m of the water column. Under ice-covered conditions, pCO₂ measurements were taken through an open melt pond. At this time, elevated CO₂ concentrations were observed at the very surface (0.1 m), followed by a sharp decrease to approximately 1 m depth, coinciding with the ice-water interface. Below this depth, pCO₂ increased again, though remained well below atmospheric concentrations (Fig. 3a).

Figure 3 Measured Young Sound pCO₂ profiles (a) prior to sea ice breakup (measured through open melt pond), (b) during sea ice breakup and (c) after sea ice break up measured through CO₂ equilibration and calculation from carbonate chemistry parameters (DIC & TA).

During the period of sea ice breakup, when ice coverage ranged from approximately 30 % to 10 %, the vertical distribution of pCO₂ exhibited a similar structure. Concentrations were highest near the surface, declined to a local minimum at 1 to 2 m, and then stabilized below 3 m (Fig. 3b). Following the complete breakup of sea ice, pCO₂ showed a more gradual decrease from the surface down to about 3 m, beneath which concentrations remained relatively constant (Fig. 3c). Across all three observational periods, a shallow surface layer approximately 5 m thick was identified, within which most of the pCO₂ variability occurred. Below this depth, pCO₂ remained relatively constant.

These vertical structures are consistent with strong physical stratification, likely driven by freshwater input from glacial melt and surface heating. Temperature and salinity profiles collected concurrently support the presence of sharp vertical gradients in the upper water column, with salinity ranging from 1.4 to 29.6 PSU and temperature from −0.4 to 6.2 °C. These physical profiles, shown in Fig. 4, confirm that vertical mixing was strongly suppressed during the observational period. The small tides present in Young Sound, combined with the continued input of freshwater from first sea ice and then glacial melt, lead to a fjord system with distinct surface stratification throughout the months of July and August.

Figure 4 Measured Young Sound profiles of under-ice water and open water salinity and temperature.

Measurements from a different fjord (Henson et al., 2025b) in East Greenland on 4 June 2025, revealed strikingly similar vertical pCO₂ heterogeneity (Fig. 5). Elevated pCO₂ at 0.1 m decreased to a minimum around 1–1.5 m before increasing again and stabilizing near 3 m depth. Extreme stratification in the upper few meters caused pCO₂ levels in each profile to vary by more than 100 µatm between the surface and 1 m. This repeated observation of similarly structured nonlinear pCO₂ profiles, measured 8 years later and in a different fjord system, suggests that this heterogeneity is not an isolated phenomenon. Instead, these C-shaped vertical pCO₂ distributions may represent a characteristic response of strongly stratified Arctic surface waters during sea ice breakup. Indeed, Arctic surface stratification during sea ice breakup induces chemical changes that may influence the way we estimate air-sea exchange of CO₂.

Figure 5 Measured pCO₂ (a) and salinity (b) profiles at 4 locations in Tasiilaq Bay. Profiles were measured on 4 June 2025 during the period of sea ice breakup following the method in Sejr et al. (2011).

4 Discussion

Air-sea CO₂ fluxes in Arctic coastal areas are generally estimated using bulk parameterization models (Henson et al., 2024; Meire et al., 2015; Roobaert et al., 2019; Sejr et al., 2011). These models rely on several key assumptions, including unstratified surface conditions, a linear pCO₂ gradient within the diffusive boundary layer, and a vertically uniform pCO₂ profile within the mixed layer. Our observations challenge the applicability of these assumptions in Arctic coastal waters in several important ways. We observe distinct vertical heterogeneity in both physical and chemical parameters, which leads to distinct differences in flux estimates based on the depth of the measurements used.

4.1 Stratified conditions in the marine CO₂ system

Vertical pCO₂ profiles collected during July 2017 revealed pronounced non-linear behavior in the upper 3 to 5 m of the water column (Fig. 3). This directly contradicts the assumption that the ΔpCO₂ accurately represents the difference between the atmosphere and the “well-mixed bulk fluid” below the diffusive layer (Wanninkhof et al., 2009). Under ice-covered conditions, the lowest pCO₂ values (∼ 150 ppm) were consistently observed just beneath the sea ice, with concentrations increasing with depth and stabilizing around 5 m (Fig. 3a). During the ice breakup stage, a similar pattern emerged, although the minimum pCO₂ was higher (∼ 250 ppm).

More recent measurements from Tasiilaq Bay in June 2025 demonstrate very similar vertical pCO₂ profiles. Indeed, 4 high-resolution profiles with measurements every 0.25 m reveal the same C-shaped pCO₂ variation. Like in Young Sound, the most elevated pCO₂ levels were observed near the surface, and pCO₂ minima occurred near 1–2 m depth before increasing and becoming stable. This repeated observation in a different fjord system, but during the period of sea ice breakup indicates this vertical variability may be representative during stratified Arctic conditions.

Several interacting processes influence surface water chemistry during ice breakup. Low surface water pCO₂ values reflect the influence of low-salinity meltwater from snow and sea ice or glacial meltwater found in freshened Arctic waters (Geilfus et al., 2015; Henson et al., 2025a). However, surface water chemistry during the ice breakup period is further complicated by processes such as ikaite (CaCO₃⋅6H₂O) dissolution (Miller et al., 2011; Rysgaard et al., 2013; Søgaard et al., 2013) and high under-ice primary production (Søgaard et al., 2021; Verdugo et al., 2025). Additionally, snowmelt, characterized by low pH and ionic strength (de Caritat et al., 2005), may further alter carbonate system dynamics in the upper water column.

Two mechanisms may explain the nonlinear C-shaped trend in pCO₂ observed in the top few meters. First, as demonstrated by Henson et al. (2025a) mixing between glacial meltwater and seawater can result in nonlinear behavior in pCO₂, even when DIC and TA mix conservatively. In such cases, initial freshwater dilution leads to dramatically reduced pCO₂, but at very low salinities, the diminished buffering capacity can cause acidification to occur and pCO₂ to increase again. In fact, Henson et al. (2025a) present a U-shaped pCO₂ curve along a salinity gradient that appears extraordinarily similar to the C-shaped pCO₂ curve we observe with depth. Indeed, the salinity gradient created in their mixing experiments can be observed vertically in the highly stratified surface waters during sea ice breakup (Figs. 3a and 5a). Although, Henson et al. (2025a) focused on the influence of glacial meltwater, our results suggest similar processes could occur in systems influenced by sea ice and snowmelt.

Both glacial meltwater and sea ice have low DIC concentrations and act to dilute the inorganic carbon of the surface ocean (Fig. S4). However, changes in alkalinity can also impact the buffering capacity of the water mixture, leading to nonlinear effects. If the meltwater has a lower TA : DIC ratio than seawater, due to the absence of ikaite, acidification and a shift in carbonate equilibria at very low salinities could lead to higher pCO₂ values at the surface. During July 2017, Young Sound showed both diluted DIC and TA levels in upper few meters, suggesting pH change during sea ice break up could occur more easily (Fig. S4). Indeed, calculated pH profiles indicated variable surface conditions between periods of sea ice cover and sea ice breakup (Fig. S5). In this very fresh surface layer, diminished pH may elevate pCO₂ relative to waters around 1 m depth, where freshwater-seawater mixing ratios are more moderate and seawater buffering leads to very low CO₂ concentrations.

A second, but less likely, explanation involves atmospheric equilibration of sea ice melt ponds before draining into open leads. The relatively elevated pCO₂ observed at ∼ 0.1 m depth could reflect such partial equilibration. While chamber-based studies (e.g. Geilfus et al., 2012, 2015; Nomura et al., 2010; Semiletov et al., 2004) have demonstrated both uptake and efflux of CO₂ in melt ponds, equilibrium times between melt-pond water and atmosphere depend upon pond depth, wind speed, and carbonate chemistry. For example, a 0.1 m deep pond under low wind conditions (∼ 2 m s⁻¹) may reach atmospheric equilibrium in 1–4 d. However, in our case, pCO₂ values calculated from TA and DIC in melt ponds did not indicate equilibrium with the atmosphere, making this explanation less likely than the freshwater mixing mechanism. Furthermore, Verdugo et al. (2025) demonstrate that melt pond drainage in Young Sound was responsible for reducing pCO₂ levels. Though this is conclusion relies on measurements from 2.5 m depth. Nevertheless, atmospheric equilibration may play a role after the sea ice barrier is removed. Elevated pCO₂ levels at the surface (0.1 m) post sea ice breakup may result from the combination of the chemical changes described above, heating from solar radiation, and from atmospheric CO₂ uptake (partial equilibration) in the limited volume of this freshwater lens.

As melt progresses and sea ice recedes, riverine input and vertical mixing become more influential. Yet even after ice breakup, surface waters often remain fresh due to glacial meltwater runoff, and the resulting low salinities help maintain stratification. In August 2017, vertical structure remained pronounced even after sea ice breakup, with elevated pCO₂ at 0.1 m which stabilized below ∼ 3 m. In other words, near-surface conditions remained decoupled from deeper waters. This persistent shallow layer, characterized by low salinity, higher temperature, and elevated pCO₂, suppresses gas exchange with the colder, more undersaturated water below, consistent with observations by Dong et al. (2021). In such environments, bulk flux models that assume homogeneity and linear gradients are likely to yield biased or inaccurate estimates.

To place these 2017 measurements in historical context, we examined long-term surface water pCO₂ data collected at 0.5–1 m depth by the Greenland Ecosystem Monitoring (GEM) program between 2007 and 2023. These data, measured using consistent protocols, are presented in Fig. 6 alongside our open-water profiles. Over the 17-year record, August pCO₂ concentrations at ∼ 1 m depth had ranged from 220 to 408 µatm and had consistently remained below atmospheric levels. This apparent stability has contributed to the perception of sustained CO₂ uptake throughout the summer season. However, the high-resolution vertical profiles obtained during the 2017 field campaign add nuance to this assumption. Elevated pCO₂ levels confined to the uppermost meter of the water column may go undetected in standard monitoring approaches that rely on fixed-depth sampling. These results suggest that the dynamic changes during sea ice melt can induce episodes of slowed uptake, air-sea equilibrium, or CO₂ outgassing. Consequently, existing sampling protocols may underestimate surface variability and bias flux estimates, especially in stratified conditions where near-surface chemistry is decoupled from subsurface layers.

Figure 6 Measured Young Sound pCO₂ profile after ice break up in 2017 compared with historical variation in pCO₂ at 1 m depth in the same location.

4.2 Influence of surface stratification on fluxes

To assess how bulk models function when estimating CO₂ fluxes in an Arctic fjord influenced by sea ice and snow melt, we calculated fluxes using seawater pCO₂ measurements from multiple depths and two gas transfer velocity parameterizations. Specifically, we computed fluxes throughout July using pCO₂ measured at 0.1, 1, 2, and 4 m. To estimate the surface (interface) pCO₂ at 0 m, we adjusted the 1 m pCO₂ measurements to a derived skin temperature (Table 1), estimated from sensible heat fluxes (Fig. S6) following the parameterization of Smedman et al. (2007). Accounting for this skin layer correction is critical, as Woolf et al. (2016) demonstrated that neglecting the thermal skin and relying only on bulk sea surface temperature can introduce significant errors in flux estimates.

Table 1 CO₂ fluxes calculated based on pCO₂ measured at the different depth. The fluxes are calculated using the bulk model of Ho et al. (2006) and Nightingale et al. (2000). We have used locally measured wind speeds for the calculations to match flux measurements captured by eddy covariance.

^a Denotes skin temperatures derived from heat fluxes. ^b Denotes pCO₂ values estimated from measurements at 1 m depth and adjusted to derived skin temperatures.

The resulting calculations (Table 1) show that estimated CO₂ fluxes vary significantly depending on the depth of the pCO₂ measurement. Notably, fluxes derived from 0.1 m differ markedly from those based on deeper values. Since many studies rely on pCO₂ measured at a fixed depth (often at 1 m or at a ship's seawater intake below 5 m), these results underscore the potential for misrepresentation of flux magnitude and direction due to vertical heterogeneity in surface water chemistry.

Measurements from both Young Sound and Tasiilaq demonstrate that during sea ice breakup, pCO₂ levels are most elevated at the surface. This may be linked to acidification of the most freshened 0.5 m and a shift in the marine carbonate system, or partial equilibration due to air-sea gas transfer. If this acidified freshwater lens warms, for instance, due to solar radiation, pCO₂ may rise leading to oversaturation relative to atmospheric concentrations. Indeed, when pCO₂ measurements on 31 July were corrected for skin temperature, to estimate pCO₂ at the boundary layer, they suggested a transition from undersaturation to oversaturation (Table 1). While we did not directly observe this oversaturation in the vertical profiles, this likely reflects the inability to sample at the sea surface. Nevertheless, the occurrence of stratification-related vertical pCO₂ heterogeneity, with levels most elevated at the surface, will slow air-sea gas transfer compared to when conditions are well mixed. Meanwhile, warming or acidification of the thin surface layer may periodically induce a reversal of flux direction, as seen in some micrometeorological studies in Arctic coastal environments during sea ice breakup (e.g. Butterworth et al., 2025).

Overall, these findings echo those of Miller et al. (2019), who reported pronounced spatial heterogeneity in Arctic coastal pCO₂ and large differences in estimated fluxes depending on the sampling depth. The broader implications of this heterogeneity for seasonal or regional flux estimates remain unclear. However, if fluxes are upscaled from sparse, single-point measurements (e.g., once per month, as in Laruelle et al., 2014), substantial errors may result due to unrecognized spatial and temporal variability. Thus, our results emphasize the need for continuous, high-resolution observations of air-sea CO₂ fluxes, particularly in Arctic coastal systems affected by stratification and meltwater input. These observations will be essential for refining flux parameterizations, reducing uncertainty in carbon budget estimates, and improving the representation of Arctic shelf systems in global carbon models.

5 Conclusions

During the summer thaw, carbon chemistry and pCO₂ dynamics in Arctic coastal surface waters are significantly altered by the combined effects of snow and sea ice melt, terrestrial runoff, and biological activity. These influences lead to substantial variability in surface temperature, pH, dissolved inorganic carbon (DIC), and total alkalinity (TA), ultimately disrupting carbonate system equilibrium in the upper water column. As a result, estimating air-sea CO₂ fluxes using traditional bulk models becomes highly uncertain during this period.

The sea ice breakup period, typically lasting 2–4 weeks, represents a particularly dynamic and complex phase in the annual cycle. Despite its brevity, this phase may have a disproportionate influence on total summer CO₂ uptake, given that open-water conditions in high Arctic fjords are limited to only 80–120 d yr⁻¹ (Sejr et al., 2011).

Air-sea gas exchange rates depend not only on the pCO₂ gradient between the atmosphere and surface water, but also on rapid, nonlinear changes in surface water chemistry driven by the composition and volume of meltwater and runoff. Accurate flux estimation will therefore require knowledge of the depth at which surface water pCO₂ becomes vertically homogeneous, combined with gas exchange parameterizations tailored to highly stratified and ice-affected conditions. Profiling pCO₂ in the upper water column is therefore essential to identify this depth and to constrain surface flux estimates reliably.

Several eddy covariance studies in other arctic environments report variable uptake and efflux of CO₂ during the sea ice breakup period (e.g. Butterworth et al., 2025). Building on these findings, our data provides the first step in understanding potential drivers behind this variability. Though, proper quantification of the mechanisms driving nonlinear pCO₂ profiles and the resulting uncertainty of flux estimates will require observations spanning the air-sea boundary. The current lack of these datasets underscores the need for studies that integrate continuous, direct CO₂ flux measurements with detailed observations of surface water carbonate chemistry, atmospheric forcing, skin temperature, and turbulence at the air-ice-water interface.

Such integrated measurements are critical to improving our understanding of the drivers and net effect of sea ice melt-driven changes on CO₂ fluxes in Arctic coastal systems. Ultimately, this knowledge is essential to accurately quantify the seasonal and regional uptake of atmospheric CO₂ in the rapidly changing Arctic.