Inorganic Carbon and Water Masses in the Irminger Sea since 1991

. The subpolar region in the North Atlantic is a major sink for anthropogenic carbon. While the storage rates show large interannual variability related to atmospheric forcing, less is known about variability in the natural Dissolved Inorganic Carbon (DIC) and the combined impact of variations in the two components on the total DIC inventories. Here, data from 15 cruises in the Irminger Sea covering the 24-year period between 1991 and 2015 were used to determine changes in total 5 DIC and its natural and anthropogenic components. Based on the results of an extended Optimum Multiparameter Analysis (eOMP), the inventory changes are discussed in relation to the distribution and evolution of the main water masses. The inventory of DIC increased by 1.43 (cid:6) 0.17 mol m (cid:0) 2 yr (cid:0) 1 over the period, mainly driven by the increase in anthropogenic carbon (1.84 (cid:6) 0.16 mol m (cid:0) 2 yr (cid:0) 1 ), but partially offset by a loss of natural DIC (-0.57 (cid:6) 0.22 mol m (cid:0) 2 yr (cid:0) 1 ). Changes in the carbon storage rate can be driven by concentration changes in the water column, for example due to ageing of water masses, or by 10 changes in the distribution of water masses with different concentrations, either by local formation or advection. A decomposition of the trends into their main drivers showed that variations of natural DIC inventories are mainly driven by changes in the layer thickness of the main water masses, while anthropogenic carbon is most affected by concentration changes. The storage rates of anthropogenic carbon are sensitive to data selection, while changes in DIC inventory show a robust signal on short timescales, associated with the strength of convection.


Introduction
Since the industrial revolution, atmospheric CO 2 levels have been increasing almost exponentially as a result of human activities such as fossil fuel burning, cement production, and land use changes. The global ocean has acted as strong sink for this anthropogenic CO 2  and is currently taking up approximately 25 % of the annual emissions (Le Quéré et al., 2016). While the ocean has capacity to store almost all of the anthropogenic CO 2 released to the atmosphere, the emissions 20 currently outpace the oceanic absorption rates (Sabine and Tanhua, 2010). This is because the transport of anthropogenic CO 2 from the atmosphere into the ocean interior is limited by the rate of vertical exchange between the surface and deep ocean (Sarmiento and Gruber, 2002). Warming of the ocean will decrease this rate as a consequence of the increased strat-ification, and Earth System Models predict a decline in oceanic anthropogenic CO 2 uptake efficiency over the 21st century (Friedlingstein et al., 2006;Schwinger et al., 2014). Warming of the ocean will also affect CO 2 solubility, primary production and other factors governing the distribution and inventory of natural carbon in the ocean (Arora et al., 2013;Schwinger et al., 2014). It is important to constrain the magnitude of these feedbacks for policy planning, but current estimates vary significantly among models. Observational based quantitative and qualitative insight in carbon cycle climate interactions are important for 5 the further improvement of projections of the future ocean carbon cycle.
Over the past three decades, ocean CO 2 chemistry data have been collected on a regular basis in the world's oceans. As the observational record grows, direct evidence of the climate sensitivity of the marine carbon cycle emerges. For example, the Southern Ocean carbon sink exhibits clear variations in response to atmospheric circulation patterns; the sink was 10 weakening from the early 1980s to the early 2000s (Le Quéré et al., 2007), but has strengthened in the more recent decades (Landschützer et al., 2015). In the subarctic western North Pacific, measurements from 1992 to 2008 at the two time series stations KNOT and K2 reveal decadal trends in total dissolved inorganic carbon (DIC), related to alkalinity driven reductions in CO 2 outgassing (Wakita et al., 2010). In the Mediterranean Sea, changes in the large-scale circulation result in variability of the anthropogenic CO 2 concentration (Touratier and Goyet, 2009). Within the subpolar North Atlantic, high-quality carbon 15 data have been collected almost every second year since the early 1990s , enabling the determination of subdecadal variability. This shows relationships between the anthropogenic CO 2 storage rate and the extent and intensity of ventilation processes, primarily driven by the North Atlantic Oscillation (NAO) Wanninkhof et al., 2010;Fröb et al., 2016;Woosley et al., 2016). 20 The subpolar North Atlantic is a key region for the storage and transport of CO 2 in the global ocean .
While a response of anthropogenic CO 2 storage to atmospheric forcing has been determined as mentioned above, less is known about variations in natural DIC, any relations to atmospheric forcing and relevance for total DIC inventories (Tanhua and Keeling, 2012). Here, we analyse changes in total DIC and its natural and anthropogenic components in the central subpolar North Atlantic, the Irminger Sea, in relation to the distribution and evolution of water masses over a 24-year period from 1991 to 2015, 25 covering three periods of variable convective activity (Fröb et al., 2016).

Hydrographic setting
The Irminger Sea, a central sea in the subpolar North Atlantic (Fig. 1), is a climatically sensitive area with strong hydrographic contrasts. The subpolar North Atlantic circulation pattern has been extensively presented in the literature; here the description follows Lavender et al. (2005) and Våge et al. (2011). In the upper ocean, the East Greenland Current carries cold and fresh 30 water of Arctic origin southwards, in the west, close to the shelf of Greenland. In the east, the Irminger Current carries warm and salty water northwards along the Reykjanes Ridge. The salinity and temperature signature of these warm water masses is affected by the strength and shape of the subpolar gyre (Häkkinen and Rhines , 2004;Hátún et al., 2005). South of the Denmark the Earth's surface (Amante and Eakins, 2009) using a nearest-neighbour interpolation. The ratio between the bottom depth of this bathymetry and the reported cruise station bottom depth was multiplied with the pressure coordinates of each station. This normalization step mainly affected the adjusted FOUREX data, while for the other cruises the normalization changed sampling depths by less than 20 meters. 5 The accuracy of the GLODAPv2 data product is better than 0.005 in salinity, 1 % in O 2 , 2 % in nitrate (NO 3 ), 2 % in silicate (SiO 2 ), 2 % in phosphate (PO 4 ), 4 µmol kg −1 in DIC and 6 µmol kg −1 in total alkalinity (A T ) . For 29AH20120623, the overall accuracy of NO 3 , PO 4 and SiO 2 was 1 %; the accuracy of DIC was 2 µmol kg −1 and for A T it was 4 µmol kg −1 (Ríos et al., 2015;García-Ibáñez et al., 2016). For the SNACS cruise in 2015, pressure, conductivity, temperature and dissolved O 2 were directly measured with a Seabird 911+ CTD profiler. At every station, water samples were obtained at 12 10 depths using Niskin bottles, and used to calibrate the CTD measurements following the Global Ocean Ship-based Hydrographic Investigations Program (GO-SHIP) calibration procedure (Hood et al., 2010). The accuracy of bottle salinities, analysed with a salinometer, was ±0.003. The accuracy of O 2 concentration measured with Winkler titration using a potassium iodate solution as a standard was 0.2 µmol kg −1 . The precision was better than 2 % in PO 4 , 1 % in SiO 2 and 1 % in NO 3 as evaluated using samples drawn from sets of Niskin bottles tripped at the same depth. DIC and A T was measured according to Dickson et al. 15 (2007), with an accuracy of 2 µmol kg −1 for both (Fröb et al., 2016).
The seawater CO 2 chemistry can be fully described if at least two of the four variables DIC, A T , CO 2 partial pressure or pH are known. The measured variables at each of the 15 cruises are listed in Table 1. For six cruises, A T and pH were measured, therefore DIC was calculated for these, using the dissociation constants of Lueker et al. (2000). For three cruises, only DIC was 20 measured. For these, A T was approximated using the salinity-alkalinity relationship for the North Atlantic of Lee et al. (2006). This relationship is defined for the surface ocean only, therefore its validity for the deep Irminger Sea was tested (Appendix A).
The mean difference between approximated and the measured A T data available was less than 5 µmol kg −1 ; this is better than the target accuracy of A T of the GLODAPv2 data product. No bias with depth or position was evident.

25
The total DIC concentration is partitioned into its natural and anthropogenic components (DIC = DIC nat + C ant ). The C ant concentration was estimated with the φC • T method (see Sect. 4.1). The DIC nat concentration is the difference between DIC and C ant . For all cruises, the column inventories were estimated for DIC, DIC nat and C ant . The inventories are sensitive to depth, therefore column inventories were only estimated for the part of the transect covered by all 15 cruises, between 40.5 • W and 31.5 • W. The column inventory is the concentration profile integrated over the entire water column (Tanhua and Keeling, Here, Inv s is the column inventory of any species s, c s its concentration, ϱ the density at in situ temperature and pressure and z the depth of the water column. The storage rate is the slope of a linear least-squares regression over the mean column inventories with time. The standard error of the slope is the error of the storage rate. Changes in inventories can be caused by changes in the distribution of water masses with different species concentrations or by changes in species concentration within the water masses. The distribution of water masses was determined using an extended Optimum MultiParameter analysis 5 (eOMP, see Sect. 4.2). While the hydrographic parameters that describe a set of source water types (SWTs) used for the eOMP analysis are assumed to be time independent, the concentrations within each water mass of species such as DIC or C ant vary over time and can therefore only be resolved by evoking the concept of water mass mixing averaged concentration, i.e.
archetypal concentration (Álvarez-Salgado et al., 2013) (see Sect. 4.2). Finally, the inventory changes can be decomposed into contributions from changes in the archetypal concentration of the source water types and from changes in layer thickness of 10 each water mass assuming linearity: Here, ∂Inv ∂c dc dt is the mean layer thicknesses with variable archetypal SWT concentrations, while ∂Inv ∂z dz dt can be calculated as the mean archetypal SWT concentrations multiplied by the layer thickness changes over a specific time period. Hence, the two drivers for the observed inventory variability of total DIC and its natural and anthropogenic components can be identified.

Anthropogenic CO 2 calculation
The φC • T method was applied to all cruises in the Irminger Sea to estimate C ant concentrations (Pérez et al., 2008;Vázquez-Rodríguez et al., 2009, 2012. The φC • T method is a back-calculation method that follows the same principles as the ∆C * method of Gruber et al. (1996). In the φC • T method, C ant is quantified as the difference between the preformed DIC at the time t and at preindustrial times (π): C ant = DIC •,t − DIC •,π . DIC •,t is calculated by correcting the measured DIC for changes due to remineralisation of 20 organic matter and CaCO 3 dissolution, while DIC •,π is quantified as the sum of the saturated DIC concentration with respect to the preindustrial atmosphere and the air-sea CO 2 disequilibrium (∆C dis ). The approach involves the following basic features: The subsurface layer (100-200 m) preserves conditions during water mass formation and is therefore taken as a reference.
∆C dis is parameterized based on subsurface data using a short-cut approach to calculate Cant. The set of parameterizations for preformed alkalinity (A • T ) and ∆C dis obtained from the subsurface data are applied directly to waters with temperatures larger 25 than 5 • C. For waters below the 5 • C isotherm, A • T and ∆C dis estimates are based on an eOMP analysis, which was successfully used in previous studies (Pérez et al., 2008;Vázquez-Rodríguez et al., 2009, 2012. This eOMP determines in each sampling point the fraction of 6 water masses that ventilate the global ocean, taking different formation histories and water mass origins into account. Each water mass has assigned values for A • T and ∆C dis and together with the obtained fractions, A • T and ∆C dis are calculated. Note that this eOMP analysis is only used to determine A • T and ∆C dis , which are used in the C ant calculation. The major advantage of the φC • T method over other back-calculation methods is that is does not rely on measurements of age tracers, such as chlorofluorocarbons (CFCs). Further, the parameterized A • T is corrected for effects of CaCO 3 dissolution changes and the sea surface temperature increase since preindustrial times and any spatial and temporal variability of ∆C dis is taken into account. Overall, the uncertainty of φC • T derived C ant has been reported to be 5 µmol kg −1 (Pérez et al., 2008;Vázquez-Rodríguez et al., 2009). The Optimum MultiParameter (OMP) analysis (Tomczak and Large, 1989) is used to estimate the contribution of water masses, which are represented through SWTs, to each water parcel along the Irminger Sea sections. The OMP analysis assumes that all hydrographic parameters describing the water masses are affected by the same mixing processes. For each sampling point the contribution of the various water masses is quantified from an over-determined system of linear mixing equations, which is 10 solved in a non-negative least square sense assuming that the parameters are linearly independent: where G is the SWT matrix containing their properties, x the relative contributions of each SWT to the sample, d the observed data and Res the residual, which is minimized in the calculation. The OMP was further developed into the extended OMP (eOMP) analysis, by Karstensen and Tomczak (1998); in the present analysis an eOMP has been adopted. This accounts 15 for the non-conservative behaviour of O 2 and nutrients by using Redfield ratios. In the eOMP, the remineralization of NO 3 and PO 4 is numerically related to an oxygen consumption rate, which, if multiplied with a pseudo-age, is similar to apparent oxygen utilization (AOU) (Poole and Tomczak, 1999). OMP and eOMP analyses have previously been used to describe in  Table 2, including their standard deviations. These values were determined from the 10 % of data in the relevant density class that were closest to the property maximum or minimum used to delineate the SWT. For example, if a SWT was defined as a salinity minimum, all data points within a specific potential density (σ) range were sorted by salinity and the mean and standard deviation over the first 10 % of the data points, gave the salinity properties for that SWT. The approximate locations of all SWTs are shown DSOW is the densest water mass in the Irminger Sea and defined as an O 2 maximum at σ 2 levels denser than 37.10 kg m −3 . ISOW is defined as a salinity maximum between 36.89 and 37.10 kg m −3 (σ 2 ) and θ between 2. σ 2 range, effectively separating uLSW and cLSW. IW is a saline water mass, depleted in O 2 , and of southern origin (Sarafanov et al., 2008). IW was identified by O 2 values below 250 µmol kg −1 at σ 0 between 27.45 and 27.65 kg m −3 . ISW is fresh, elevated in O 2 , and found between 4 and 5 • C. 15 Finally the Subpolar Mode Water (SPMW), the Subarctic Intermediate Water (SAIW) and the North Atlantic Central Water (NACW) are all typically found in the upper Irminger Sea. SPMW is oxygenated and of subpolar origin. It was defined as a salinity maximum in the 7-8 • C θ range. SAIW is a salinity minimum in the 6.5-7.5 • C θ range. NACW was defined as the salinity maximum for θ above 9 • C. 20 The seven hydrographic parameters, describing the SWTs, limit the number of SWTs included in one eOMP analysis to a maximum of seven. In addition, the required mass conservation over-determines the system of linear equations. However, NO 3 and PO 4 might be locally correlated, therefore one degree of freedom for the eOMP analysis is potentially lost. Therefore, the Irminger Sea was divided into four regions, defined such that each contained a maximum number of 5 SWTs to be determined  While the density boundaries were used to identify the set of SWTs potentially present, the eOMP analyses were performed at each sampling point. The equations were normalized and weighted, accounting for differences in measurement accuracies and potential environmental variability. Weights were assigned according to the variability and accuracy of the parameters following García-Ibáñez et al. (2015). The highest weight was assigned to mass to ensure its conservation. The second highest weights were assigned to θ and salinity, because they are the most accurate. PV was weighted high as well due to its good 35 accuracy and to enable resolution of both LSW classes. The eOMP results in a ratio r ij , this describes the contribution of each SWT, i, to each data point in space and time, j.
In order to determine SWT concentrations of time-varying species such as DIC, DIC nat and C ant , the mixing-weighted concentration or archetypal concentration, C i , of these species was calculated for each SWT, i (Álvarez-Salgado et al., 2013): Here, the concentration in each sampling point j, C j , is multiplied with the ratio of the SWT in that point r ij . When summed over all points and divided by the total fraction the SWT occupies, this estimates C i . Further, the layer thickness, Th, of each SWT at each station is estimated according to: Here, is the fraction of the total water column that each SWT occupies at each station, with n k describing the number of sampling points per station. The fraction is unitless and needs to be scaled to the total water column height, i.e. multiplied by the maximum depth of station, d k,max . The average layer thickness over all stations of the Irminger Sea transect is the mean layer thickness. As each water parcel is a mixture of different water masses represented by r ij , Eq. (5) allows to convert each composite to a measure of height.

Uncertainty analysis
Uncertainties for the distribution of water masses result from measurement uncertainties and errors of the eOMP analysis.
Here, the largest source of errors is the definition of the SWTs. The SWT matrix needs to represent the known features of the circulation (Tanhua et al., 2005), but temporal shifts in SWT characteristics cannot be accounted for with the eOMP analysis.
A measure of uncertainty is given by the difference between measured and eOMP calculated values, the residual R in Eq.
(3). 20 The total residual, calculated by taking the square of the largest parameter residual at each sampling point (García-Ibáñez et al., 2015), and the individual parameter residuals are shown in Fig. 4. Below 1200 m, the total residual is close to zero, as are the residuals of θ, salinity and O 2 . In the intermediate and surface ocean these residuals increase, in particular for O 2 , which might be a consequence of gas exchange. The residuals of PO 4 , NO 3 and SiO 2 are larger, as expected due to their lower weights in the eOMP, and do not show any trend with depth. The mean error for each parameter is listed in Table 2. These are of similar 25 magnitude as the errors determined for the Irminger Sea eOMP analysis by Tanhua et al. (2005).
In order to test how robust the results of the eOMP analysis are, a Monte-Carlo simulation was performed (Tanhua et al., 2005). The properties of the SWT matrix were randomly perturbed, within the standard deviation of each parameter. 100 of such perturbed SWT matrices were created and the eOMP was solved for each perturbed system. This allows for quantification 30 of the sensitivity of the eOMP to potential temporal variations of the SWT properties. The standard deviation of the mean SWT contribution over all 100 perturbations is shown in the last column of Table 2. The uncertainties are generally low, hence the robustness of the eOMP analysis is high.
The layer thickness uncertainties were estimated by scaling the averaged standard deviation of each SWT, which were quantified with the Monte-Carlo simulation, to the width and depth of the Irminger Sea. The uncertainty of C ant concentra-5 tions is 5 µmol kg −1 , for DIC and for DIC nat it is 4 µmol kg −1 . Errors for the inventories were estimated by propagating the uncertainties of the layer thicknesses and the concentrations through the water column.

Results
Over the 24-year period considered here, the frequency and the amplitude of the mean-winter NAO changed significantly (Hurrell and Deser, 2009). The convection in the subpolar gyre during winter, which is driven by the large-scale atmospheric Therefore, the annual change in the DIC inventory is mainly driven by the large C ant storage rate, but partially offset by the loss in DIC nat inventory. The variability of the DIC and C ant inventories over the 24-year period is of similar magnitude, as indicated by the error of the slope in Fig. 6, whereas the DIC nat inventory varies slightly more. DIC nat is increasing as water masses age and DIC nat from remineralisation of organic matter accumulate, while it decreases during water mass renewal / ventilation, which brings water with preformed, relatively low, DIC nat concentrations into the ocean interior. It is notable that the C ant inventory increased sharply from 2012 to 2015, while there was a comparably large decline in the DIC nat inventory.

5
This is not an artefact of the method, but can be explained by the fact that the 2015 data were obtained during active convection in the Irminger Sea (Fröb et al., 2016). During the strong convection, older water masses enriched in DIC nat were replaced by water masses high in Cant due to their most recent contact to the atmosphere, and relatively low in DIC nat as DIC nat from remineralisation has not yet accumulated. In contrast to that, the peak in 2005 in the C ant inventory could be related to the advection of C ant enriched water masses formed in the Labrador Sea or in the region south of Cape Farewell into the Irminger 10 Sea (Straneo et al., 2003;Palter et al., 2016), but this strong signal may also reflect the true error or reveal measurement bias.

SWT distribution
The layer thickness of the Irminger Sea SWTs from 1991 to 2015 is presented in Fig. 7. Because their individual contributions are small, the upper ocean SWTs, i.e. NACW, ISW, and SPMW, are combined and titled upper waters (UW). For reasons of simplicity, the number of SWTs were reduced from 11 to 9, by performing composite analyses for uNEADW and IcSW. 15 The uNEADW was determined to be a composite of 26 % ISOW, 14 % LSW, 58 % lNEADW and 2 % Mediterranean Water analysis. Therefore, not all 11 SWTs used for the eOMP analysis are shown, but only UW, IW, uLSW, cLSW, DSOW and ISOW. Since the FOUREX section, occupied in 1997, was located further south than the AR07E section, covered that year by 06MT19970707, the SWT distribution differs slightly between the two cruises. At the FOUREX section, the ISOW layer is on average 50 m and the SAIW layer about 15 m thicker than further north, while the ISW layer is about 47 m and the SPMW layer 19 m thicker for 06MT19970707 than at the FOUREX section. For the other SWTs, differences are smaller than 8 m. 25 Relatively speaking, 50 m is less than 2 % of the entire water column, so that the discrepancy between the two cruises is small compared to the mean depth of the Irminger Sea.
The distribution of SWTs in the Irminger Sea, as shown in Fig. 7, changes substantially from 1991 to 2015. The overall trend is indicated, but the rates of change are often larger, if sub-periods are considered. For example, the IW layer thickens the layer thickness of uLSW increases substantially. The cLSW layer shows the largest changes in thickness at a loss rate of -54.0±4.3 m yr −1 from 1991 to 2015. The convective activity in the subpolar gyre in 1991 to 1997 led to extensive production of cLSW. After that, the cLSW layer was not renewed and strongly diminished until 2015. In contrast to that, the DSOW layer decreases only little in thickness over the 24-year period covered by the data. The ISOW layer thins by -9.9±3.1 m yr −1 . Over the entire period from 1991 to 2015 the UW layer thickness increases at a rate of 26.2±4.8 m yr −1 , but it has essentially a constant thickness from 1991 to 2000, then thickens until 2007 and decreases in thickness after that. Overall, the change of the 5 distribution of the main SWTs is well captured by the eOMP analysis. Especially the transformation of the two LSW classes seems to match the observations well .
The mean layer thickness for each SWT in the time periods considered here are summarized in Table 3. The main variability in the distribution of water masses is created by layer thickness changes of cLSW, uLSW and UW. With only small changes 10 over the 24-year period, DSOW, ISOW and IW occupy a little more than a third of the Irminger Sea sections. In the mid-1990s, the cLSW layer occupied close to 50 % of the entire water column, which left thin uLSW and UW layers, corresponding to less than 7 % and around 8 %, respectively, of the entire water column. As cLSW was advected out of the Irminger Sea, while not being re-formed by convection, that layer was replaced mainly by UW in the early 2000s. In 2004, cLSW occupied 37 % of the water column and UW 24 %, while the uLSW layer only accounted for 4 %. With the recurring convection events between 15 2008 and 2015, uLSW was formed more frequently and displacing UW as well as the remainder of the cLSW layer. The measurements from winter 2015 reveal that by then the fraction of cLSW was as low as 17 %, and that of UW only 5 %, but that uLSW occupied 45 % of the entire water column.

Archetypal concentration changes
Within the SWTs, the archetypal concentrations of carbon species change over time due to e.g. remineralization of organic 20 matter, which add DIC nat , or air-sea gas exchange, which increases C ant in water masses that are in contact with the atmosphere. The archetypal concentrations of DIC, DIC nat and C ant from 1991 to 2015 are shown in Fig. 8, for the same SWTs as in Fig. 7, and are also summarized in Table 3. In all SWTs, DIC is similar, except for UW, which have a more variable DIC that is generally lower compared to all other SWTs. The older SWTs in the deep ocean have higher DIC nat concentrations and lower C ant concentrations compared to SWTs that have been ventilated more recently. Therefore, the deep SWTs, namely ISOW and  For all SWTs, the rate of change in the archetypal DIC concentration is similar, except for ISOW, where it is slightly smaller.
The archetypal concentration change over time for DIC nat is not statistically different from zero for most SWTs, apart from cLSW, where the DIC nat concentration increases by 0.30±0.06 µmol kg −1 yr −1 . In contrast, the archetypal C ant concentration increases significantly in all SWTs over time. Further, the increase rate of the archetypal DIC concentration is not statistically different from the rate of increase of the archetypal C ant concentration. This indicates that the increase in DIC can be explained by the input of C ant to the entire water column. This is true for all SWTs, except cLSW. For this water mass the increase in the archetypal C ant concentration contributes by 0.31±0.07 µmol kg −1 yr −1 to the increase in DIC, which is only half of the DIC concentration increase. This is because DIC nat accumulates as cLSW ages, while at the same time, a smaller fraction of C ant is added to this water mass due to less frequent ventilation.

DIC storage rate decomposition
The decomposition of the inventory changes reveals the contribution that changes in the SWT distribution and changes in the concentration within these SWTs have on the total storage rate of DIC and its natural and anthropogenic components. Figure 9 summarizes the storage rates of DIC, DIC nat and C ant over the entire 24-year time period and for the periods 1991-1997, 2000-2007 and 2008-2015. The first bar shows the total storage rate summed over all SWTs. The second bar shows the concentration 10 driven storage rate and the third bar the layer thickness-driven storage rate. In theory, the first bar should be the sum of the last two bars. However, since the storage rates have been calculated using a linear regression over only a small number of data points, the residuals can become quite large, this is especially the case for the shorter time periods. Nevertheless, some conclusions can be drawn. The increase of the DIC inventory from 1991 to 2015 is driven by the increase in the C ant inventory, partially offset by a decrease in DIC nat inventory. The rise in the C ant inventory is primarily due to a rise in C ant concentration, 15 but the contribution from layer thickness-driven changes is also significantly positive (Fig. 9a). While the rise in C ant concentration occurs in all SWTs (Fig. 8), the contribution from the latter factor appears mostly driven by the increase in the thickness of uLSW (Fig. 7), which is rich in C ant (Fig. 8). The decrease in the DIC nat inventory is the result of the layer thickness-driven reduction, which is larger than the concentration driven increase (Fig. 9a). This can be attributed to the replacement of cLSW with uLSW and UW, which both have lower concentrations of DIC nat . 20 The variations at subdecadal timescales can be understood in terms of the convective activity in the Irminger Sea, although with larger uncertainty. Figure 9b shows the decomposed trends from 1991-1997, a period when convective activity was high in the subpolar gyre. The total DIC storage rate is driven by the C ant storage rate, while the changes in the DIC nat inventory are not significantly different from zero. The C ant storage rate is mainly driven by increasing concentrations, which occur in all 25 SWTs (Fig. 8). indicative of the ventilation that occurred. As evaluated from Fig. 8, the C ant increase is greatest in IW, uLSW and UW, while the loss of DIC nat occurred primarily in the IW and uLSW. However, also cLSW appears to be affected by the most recent 10 event in 2015.

Discussion
The data that have been collected in the Irminger Sea over the past decades provide unequivocal evidence for climate forcing of the carbon cycle in this oceanic region. Over long time scales, the steady trend due to uptake of anthropogenic CO 2 clearly dominates, but at shorter time scales, it varies and can also be significantly masked by variability of DIC nat . In particular, this 15 was the case in the period from 2000 to 2007, when the negative DIC nat storage rate partially offset the increasing C ant storage, resulting in a DIC storage rate that was only about a third of the storage rates in the preceding and later time periods considered here. In that period, when convection was shallow, the replacement of relatively DIC nat rich cLSW with relatively DIC nat -poor UW and uLSW led to the loss of DIC nat . Ageing might have increased the DIC nat concentration, but the water masses in which these processes occur were flushed out the study area. Climate feedback mechanisms that involve natural carbon cycling in the The results presented here do not indicate a consistent response in the storage rates for anthropogenic CO 2 to the NAO. The storage rates were similar for the first two periods, with predominantly high and low values of the NAO index, while it was clearly larger for the last period, with predominantly high NAO index. In that time period, the C ant storage rate increased by 30 a factor of 1.6 compared to the period from 2000 to 2007 (Fig. 9). The lack of change in C ant storage rates between the two first periods, contrasts with the results of Pérez et al. (2008), who found that the storage rates were low from 1997-2006. This is a result of the differences in data: when using the same cruises and the same periods of time (i.e 1997-2006 for the middle period) as in Pérez et al. (2008) to estimate C ant storage rates, we estimated a significant decline in C ant storage rates for the middle period compared to the first (C ant storage rates: 2.78±0.28 mol m −2 yr −1 for 1991-1997: 1.06±0.47 mol m −2 yr −1 for 1997-2006).
The DIC storage, on the other hand, do show a consistent response to the NAO index, it is larger in the periods with pre-5 dominantly high NAO index winters (the first and the last) than in the middle period, with a predominance of low NAO index winters. This response was also found when using the same cruises and time periods as Pérez et al. (2008) to calculate the storage rates. Altogether, this shows that while calculation of C ant storage rates over small time periods is sensitive to data selection, calculation of DIC storage rates is not. This is not unreasonable, as estimates of C ant inventories involves more variables (for example AOU and alkalinity), which increases the risk of introducing sampling or measurement biases. Regardless, convective who also used an independent approach for estimating C ant . Most likely because the thick layer of C ant rich cLSW was mostly established in 1991 when the first data were collected, no large storage rate is observed for the early NAO positive period included here . Keeling (2012) is not significantly different from zero. The entire increase in DIC inventory is explained by an increase of the remineralised fraction, as determined from AOU and their assumed C:O ratio. This implies either that there is no storage of C ant in the region or that any increase in C ant is completely offset by reduced CO 2 solubility. 25 In order to compare to the Tanhua and Keeling (2012) results here, the DIC abio storage rates as well as O 2 and AOU storage rates are shown in Fig. 10 for the same periods as in Fig. 9. For the Irminger Sea cruise data from 1991-2015, the DIC abio storage rate is 1.40±0.17 mol m −2 yr −1 , and explains the DIC storage rate almost entirely. The negligible trend in AOU explains the lack of difference between the DIC and DIC abio storage rates. The estimated DIC abio storage rate is lower than the C ant 30 storage rate (1.92±0.17 mol m −2 yr −1 , Fig. 6), probably due to a decline in preformed DIC values, i.e. a loss of solubility as a consequence of the positive surface temperature trends in the Irminger Sea from 1991-2015 (Stendardo and Gruber, 2012;Maze et al., 2012). The long-term warming trends in the surface ocean and thus the decreasing O 2 solubility leads to a deoxygenation of -0.8±0.3 mol m −2 yr −1 over the 24-year period, which is consistent with the suggested loss in preformed DIC.
Similar to the C ant storage rates, both periods before and after 1997 show a loss in O 2 , despite stronger convection in the early 1990s. At the same time, the small increase in AOU does neither reflect the constant inventory rate in DIC nat before 1997, nor does the constant AOU storage rate reflect the loss in DIC nat after 1997. Again, this is attributed to the fact that 10 cLSW was mostly formed before 1991  AOU matches the loss in DIC nat as the water column is reventilated so that remineralized organic matter is replaced by a larger fraction of C ant . In fact, the convection in 2015 was strong enough to restore the O 2 inventory levels so that the mean inventory from 2008 to 2015 was 749±4 mol m −2 , the same level as the well-ventilated early 1990s.

Conclusions
The repeat observations within the Irminger Sea show significant changes in total, natural and anthropogenic CO 2 inventories 20 from 1991 to 2015 with large interannual variability in the natural component. The eOMP method results and the decomposition of the inventory changes give valuable insight into the driving mechanisms to interpret the observed variability. Overall, changes in layer thickness of the main water masses appear most important for the DIC nat inventory, while concentration change within these water masses is the key factor for C ant . C ant is typically more important for changes in total DIC inventories than DIC nat .
While the DIC inventory changes show a clear signal associated to the NAO, for C ant the signal is less robust, especially before 25 and after 1997, likely because the data used here does not cover the period before 1991, when the thick layer of cLSW was formed. From 1991 to 2015 the mean C ant saturation over the entire water column increased from 52 % to 67 %, increasing the C ant inventory from 53±3 mol m −2 to 117±3 mol m −2 , mainly driven by the most recent convection in 2015. Despite the negative trend in O 2 inventory from 1991 to 2015, the convection in 2015 was strong enough to replenish O 2 levels at depth, leading to a mean saturation of 89 % and an O 2 inventory of 749±4 mol m −2 , which was as high as in 1991. C ant is sensitive 30 to time period considered, while DIC appears more robust to sampling and measurement bias. Therefore, for a comprehensive view on carbon cycle feedback mechanisms, not only C ant , but also natural and total DIC should be taken into account.

Appendix A: Salinity -Alkalinity relationship
The application of the surface relationship between salinity, temperature and A T by Lee et al. (2006) is tested for Irminger Sea cruise data in the entire water column. Further, the linear relationship between salinity and A T by Nondal et al. (2009) is applied as well and compared to the Lee et al. (2006) relationship. For all cruise data between 1991 and 2015, the difference between measured and calculated A T is presented in Fig. A1. Both relationships perform reasonably well. For the Lee et al. 5 (2006) relationship there is no bias with depth, but measured A T is slightly overestimated. The Nondal et al. (2009)  Further, the impact of the overestimated calculated A T on C ant , DIC nat and DIC storage rates was tested. For that, 4.5 µmol kg −1 , which was the mean difference between measured and calculated A T based on the Lee et al. (2006) relationship, were added to the calculated A T for the 3 cruises, where only DIC was measured. The archetypal C ant concentrations for all SWTs were less than 1 µmol kg −1 smaller after the correction of A T . No significant difference between the storage rates was evident.
Appendix B: Location SWT 15 As a result of the eOMP, the fraction of SWTs in each sampling point is estimated. Figure B1 shows this fraction for the 1991 (a) and the 2015 (b) data for IW, uLSW, cLSW, DSOW, ISOW and UW, which is the sum over ISW, NACW and SPMW.
Overall, the position of the water masses is well represented through time.