# Current estimates of K1* and K2* appear inconsistent with measured CO2 system parameters in cold oceanic regions

Olivier Sulpis, Siv Lauvset, Mathilde Hagens

► To cite this version:

Olivier Sulpis, Siv Lauvset, Mathilde Hagens. Current estimates of K1* and K2* appear inconsistent with measured CO2 system parameters in cold oceanic regions. Ocean Science, 2020, 16 (4), pp.847-862. 10.5194/os-16-847-2020. hal-04700515

HAL Id: hal-04700515

https://hal.science/hal-04700515v1

Submitted on 19 Sep 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

# Current estimates of K₁* and K₂* appear inconsistent with measured CO₂ system parameters in cold oceanic regions

Olivier Sulpis¹,², Siv K. Lauvset³, and Mathilde Hagens⁴

¹Department of Earth Sciences, Utrecht University, Utrecht, the Netherlands

²Department of Earth and Planetary Sciences, McGill University, Montreal, Canada

³NORCE Norwegian Research Centre, Bjerknes Centre for Climate Research, Bergen, Norway

⁴Soil Chemistry and Chemical Soil Quality, Wageningen University and Research, Wageningen, the Netherlands

**Correspondence:** Olivier Sulpis (o.j.t.sulpis@uu.nl)

Received: 5 March 2020 – Discussion started: 25 March 2020

Revised: 11 June 2020 – Accepted: 16 June 2020 – Published: 21 July 2020

**Abstract.** Seawater absorption of anthropogenic atmospheric carbon dioxide (CO₂) has led to a range of changes in carbonate chemistry, collectively referred to as ocean acidification. Stoichiometric dissociation constants used to convert measured carbonate system variables (pH, pCO₂, dissolved inorganic carbon, total alkalinity) into globally comparable parameters are crucial for accurately quantifying these changes. The temperature and salinity coefficients of these constants have generally been experimentally derived under controlled laboratory conditions. Here, we use field measurements of carbonate system variables taken from the Global Ocean Data Analysis Project version 2 and the Surface Ocean CO₂ Atlas data products to evaluate the temperature dependence of the carbonic acid stoichiometric dissociation constants. By applying a novel iterative procedure to a large dataset of 948 surface-water, quality-controlled samples where four carbonate system variables were independently measured, we show that the set of equations published by Lueker et al. (2000), currently preferred by the ocean acidification community, overestimates the stoichiometric dissociation constants at temperatures below about 8 °C. We apply these newly derived temperature coefficients to high-latitude Argo float and cruise data to quantify the effects on surface-water pCO₂ and calcite saturation states. These findings highlight the critical implications of uncertainty in stoichiometric dissociation constants for future projections of ocean acidification in polar regions and the need to improve knowledge of what causes the CO₂ system inconsistencies in cold waters.

## 1 Introduction

In the last decades, oceans have absorbed over a quarter of the anthropogenic carbon dioxide (CO₂) emitted to the atmosphere (Le Quéré et al., 2018; Gruber et al., 2019). Upon dissolution in seawater, this CO₂ triggers a suite of reactions that lead to a range of chemical changes jointly termed ocean acidification (Zeebe and Wolf-Gladrow, 2001; Gattuso and Hansson, 2011). To accurately calculate the magnitude of these changes, it is crucial to understand the chemical behaviour of CO₂ in seawater.

Upon dissolution, CO₂ takes the form of solvated CO₂ (CO₂(aq), CO₂·H₂O) or carbonic acid (H₂CO₃), which are here both represented by H₂CO₃*, since they can only be readily distinguished by infra-red spectrometry (Zeebe and Wolf-Gladrow, 2001), and the following series of reactions occurs:

$$CO_{2(aq)} + H_2O \leftrightarrow H_2CO_3^*, \qquad (R1)$$

$$H_2CO_3^* \leftrightarrow H^+ + HCO_3^-, \qquad (R2)$$

$$HCO_3^- \leftrightarrow H^+ + CO_3^{2-}. \qquad (R3)$$

Together, these three reactions and their species constitute the marine CO₂–H₂O system, which is responsible for about 95 % of the acid–base buffering capacity of seawater and maintains the pH of the ocean within a narrow range (Bates, 2019; Zeebe and Wolf-Gladrow, 2001).

Each of these reversible reactions is associated with a thermodynamic equilibrium constant, a number that expresses the relationship between the activities of products and reactants present at equilibrium at a given temperature and pressure. Reactions (R2) and (R3) describe the first and second
step in the dissociation of carbonic acid. Their equilibrium constants are therefore termed the first and second dissociation constants, K<sub>1</sub> and K<sub>2</sub>, respectively. To avoid the use of activity coefficients, which are not straightforward to derive in seawater, marine scientists have developed a set of stoichiometric (or apparent) equilibrium constants to represent the state of the system at a given pressure (P), temperature (T), and salinity (S). To describe the carbonate system, two stoichiometric constants (K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup>, conventionally denoted by a star) are defined in terms of the concentrations of the different species:

$$K_1^* = \frac{[HCO_3^-][H^+]}{[H_2CO_3]^*}, \qquad (1)$$

$$K_2^* = \frac{[H^+][CO_3^{2-}]}{[HCO_3^-]}. \qquad (2)$$

Using these stoichiometric equilibrium constants, we can calculate the relative quantities of the dissolved inorganic carbon (DIC = [H<sub>2</sub>CO<sub>3</sub><sup>*</sup>] + [HCO<sub>3</sub><sup>-</sup>] + [CO<sub>3</sub><sup>2-</sup>]) species. With improved analytical techniques, measurement accuracy of carbonate system variables has substantially increased over the past decades. As a result, uncertainty in carbonate system calculations is currently dominated by the uncertainty of K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> values (Orr et al., 2018), justifying the need to investigate whether these uncertainties can be reduced.

Using popular software for carbonate system calculations, e.g. CO2SYS (Lewis and Wallace, 1998; Pierrot et al., 2006; van Heuven et al., 2011) or seacarb (Gattuso et al., 2019), and recently published literature as references, roughly 15 different expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> are currently in use, some of which are (partly) based on refitting data from earlier experiments (e.g. Dickson and Millero, 1987; Lueker et al., 2000); see Table S1 in the Supplement for an overview. Some expressions are based on measurements in artificial seawater of various compositions, while others were carried out in natural seawater. The vast majority of expressions were obtained in the laboratory under controlled conditions, using electrochemical cells either with (e.g. Millero et al., 2006; Millero, 2010) or without (e.g. Roy et al., 1993; Tishchenko et al., 2013) liquid junction. Within these cells, electromotive force readings of equilibrated seawater are used to compute equilibrium constants. Each expression is valid over its own range of T and S.

The various expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> obtained this way generally agree well, but discrepancies at low salinities have been highlighted (Cai and Wang, 1998; Millero, 2010; Dinauer and Mucci, 2017; Orr et al., 2018). In addition, the temperature range covered by various K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> expressions, although generally broad, only extends below 0°C in a few studies (Millero et al., 2002; Goyet and Poisson, 1989; Padimitriou et al., 2018). In fact, Mehrbach et al. (1973), who provided experimental data used by several authors to derive expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> (e.g. Dickson and Millero, 1987; Lueker et al., 2000), used data obtained at only four

different temperatures (2, 13, 25, and 35°C), which brings into question the accuracy of the temperature dependency of these constants. Bailey et al. (2018) recently suggested that the same bias exists for the dissolution of CO<sub>2</sub> in seawater and showed that previous expressions of Henry's Law constant for CO<sub>2</sub> underestimate the CO<sub>2</sub> solubility below 0°C due to a lack of samples in cold waters. As explained by Raimondi et al. (2019), because the only carbonate system variables currently measured by in situ sensor technologies are pH and the partial pressure of CO<sub>2</sub> in seawater (pCO<sub>2</sub>), relating laboratory or on-board measurements that are usually performed at temperatures ~25°C to these in situ measurements requires an accurate knowledge of the K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> temperature dependency. About 40% of the ocean volume is at an average temperature lower than 2°C, outside of the temperature range for which the Mehrbach et al. (1973) and derived constants are valid (from the data of Lauvset et al., 2016). An example of this are high-latitude cold waters, which are a critical component of the current global oceanic carbon cycle, as the Southern Ocean surface waters account for ~40% of the annual anthropogenic CO<sub>2</sub> uptake by the ocean (Landschützer et al., 2015). Given past difficulties to obtain direct pCO<sub>2</sub> measurements from ships in the Southern Ocean (Bakker et al., 2016), a number of autonomous floats have been deployed in recent years (see, e.g., Williams et al., 2017; Takeshita et al., 2018). Since these floats estimate pCO<sub>2</sub> from a pH measurement and a calculated total alkalinity (TA), our knowledge of surface pCO<sub>2</sub> in the Southern Ocean strongly relies on the accuracy of dissociation constants in these cold waters.

Best practices for oceanic carbonate system measurements generally recommend the Lueker et al. (2000) constants (Dickson et al., 2007), but the choice for a set of constants may depend on the environment and/or measured carbonate system variables. Only two of the measurable variables are required to characterize the whole carbonate system, except under conditions where substantial impact of dissolved organic carbon on TA is expected (i.e. significant organic alkalinity). Overdetermination of the carbonate system, i.e. the concomitant measurement of at least three of the carbonate system variables (1) pCO<sub>2</sub>, (2) DIC, (3) TA, and (4) pH, is often used as a tool to identify the best pair of input variables for carbonate system calculations under specific environmental conditions, e.g. in sea-ice brines (Brown et al., 2014) or in systems with substantial organic alkalinity (Koeve and Oschlies, 2012). We refer the reader to Raimondi et al. (2019) for an overview of internal consistency studies, i.e. the agreement between measured and calculated variables. Disagreement between measured and computed values may arise from uncertainties in measurements and, more importantly, equilibrium constants (Orr et al., 2018) but can also result from the choice of relationship between total boron and salinity, as well as organic alkalinity (Fong and Dickson, 2019).

Field measurements are rarely used to derive stoichiometric equilibrium constants because of their interdependence.
For example, ship-based measurements of pH are normally conducted at fixed temperature (commonly 25 °C) and converted to in situ temperature using a second input parameter as well as a set of stoichiometric equilibrium constants (Hunter, 1998). Similarly, using measured TA to calculate the contribution of the carbonate system to total alkalinity (carbonate alkalinity, CA) requires that the proton concentration and thus pH be known (Dickson et al., 2007). To the best of our knowledge, only two studies have so far used overdeterminations of the carbonate system to derive expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> (Millero et al., 2002; Papadimitriou et al., 2018). Both studies used concurrent measurements of pCO<sub>2</sub>, TA, DIC, and pH over a range of temperatures and salinities to calculate K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup>. Millero et al. (2002) used over 6000 sets of pressure-corrected field measurements. They argued that determinations of stoichiometric dissociation constants measured in natural seawater are preferable over those determined in artificial seawater and concluded that the value of K<sub>2</sub><sup>*</sup> depends on pCO<sub>2</sub>, possibly linked to organic alkalinity, which is not accounted for in carbonate system calculations. Papadimitriou et al. (2018), who focussed especially on highly saline brines down to their freezing points, used the same methods as Millero et al. (2002) for their calculations. However, instead of using field measurements, they overdetennined their system under controlled laboratory temperatures and salinities, thus avoiding temperature corrections of the pH measurements. Their work, like Orr et al. (2018), confirmed the high uncertainties associated with extrapolating expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> beyond the investigated salinity and temperature ranges.

In the present study, we use the Global Ocean Data Analysis Project version 2 (GLODAPv2; Key et al., 2015; Olsen et al., 2016) and the Surface Ocean CO₂ Atlas (SOCAT; Bakker et al., 2016) global data products to constrain stoichiometric equilibrium constants based on surface-water field measurements. Using an iterative procedure that takes into account the lack of independence of CA and pH, we quantify the temperature dependence of the stoichiometric equilibrium constants. We then use these constants to recommend input pairs for pCO₂ and CaCO₃ saturation state determinations over various temperature ranges and apply them onto a high-latitude dataset.

## 2 Materials and methods

### 2.1 Expressions for K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> as a function of carbonate system variables

Aside from recent advances that allow spectrophotometric determinations of CO₂⁻⁻⁻ concentrations (Byrne and Yao, 2008; Easley et al., 2013; Sharp and Byrne, 2019), the concentrations of H₂CO₃⁻, HCO₃⁻, and CO₂⁻⁻⁻ are normally not directly measured in seawater. Instead, at least two of the four parameters pCO₂, pH, DIC, and TA are measured and the

concentrations of the individual species inferred from them. In practical terms, TA is the sum of all bases that are titratable with a strong acid to an equivalence point corresponding to the conversion of HCO₃⁻ to H₂CO₃⁻. Here it is defined as

$$TA = CA + BA + PA + SiA + [OH⁻] - [H⁺], \quad (3)$$

with CA = [HCO₃⁻] + 2[CO₂⁻⁻⁻] and where BA is the carbonate alkalinity ([B(OH)₄⁻], PA is the phosphate alkalinity ([HPO₄²⁻] + 2[PO₄³⁻] - [H₃PO₄]), SiA is the silicate alkalinity ([SiO(OH)₃⁻], [OH⁻] is the hydroxide ion concentration, and [H⁺] is the hydrogen ion concentration. Equation (3) approximates the definition of TA provided by Dickson (1981) but does not take into account the hydrogen sulfide and ammonia acid-base systems. The terms BA, PA, and SiA can all be expressed in terms of stoichiometric equilibrium constants, total concentrations, and [H⁺]. Hence, knowing TA, the total concentrations of dissolved silicate ([DSi]), soluble reactive phosphate ([SRP]), and boron, as well as [OH⁻] and [H⁺], CA can be calculated.

To estimate K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> as a function of salinity and temperature based solely on independent measurements, we first need to define expressions that define both constants as functions of CA, DIC, pCO₂, and [H⁺]. Both K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> are normally defined in terms of proton concentration, [H⁺], and the acid-base species they describe; see Eqs. (1) and (2). In this work, we replace [H₂CO₃⁻], [HCO₃⁻], and [CO₂⁻⁻⁻] by expressions that only contain the four variables present in the dataset and Henry's constant, K₀, taken from Weiss (1974). This leads to the following set of equations, which are equivalent to those presented in Millero et al. (2002) and Papadimitriou et al. (2018):

$$K_{1}^{*} = \frac{[H^{+}](2\text{DIC} - \text{CA} - 2K_{0}p\text{CO}_{2})}{K_{0}p\text{CO}_{2}}, \quad (4)$$

$$K_{2}^{*} = \frac{[H^{+}](\text{CA} - \text{DIC} + K_{0}p\text{CO}_{2})}{2\text{DIC} - \text{CA} - 2K_{0}p\text{CO}_{2}}. \quad (5)$$

Note that similar expressions can also be derived when only three independently measured variables are available in the dataset. In this case, either K<sub>1</sub><sup>*</sup> or K<sub>2</sub><sup>*</sup> remains in the expression, in addition to any three variables of the set CA, DIC, pCO₂, and [H⁺]. Derivations of all these expressions can be found in the Supplement.

### 2.2 Data

Data for T, P, practical salinity (S<sub>p</sub>), DIC, TA, pH, [SRP], and [DSi] were taken from GLODAPv2 (Key et al., 2015; Olsen et al., 2016). Only data associated with a WOCE (World Ocean Circulation Experiment) flag of 2 were retained for this analysis. WOCE flags are associated with each GLODAPv2 variable during quality control. Data associated with a flag of 2 were assessed as “acceptable” by quality controllers of the original dataset. Recalculated or estimated
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Spatial coverage of the dataset containing GLODAPv2 and SOCAT samples for which DIC, TA, pH, and $pCO_2$ are available from independent, high-quality measurements.</figcaption></figure>

variables and samples with missing *T* or *S*<sub>P</sub> were always discarded. In total, we obtained 98326 samples for which TA, DIC, and pH are available from independent, high-quality measurements.

$pCO_2$ values were obtained from SOCAT version 3 (Bakker et al., 2016). Only data associated with a WOCE flag of 2 are used. When available, a $pCO_2$ value was selected and added to corresponding surface GLODAPv2 samples. To select the most accurate $pCO_2$ value, we only merged GLODAPv2 and SOCAT samples from the same cruise and taken within the same hour; in most cases within the same 20 min. As a result, we assembled 1024 samples for which TA, DIC, pH, and $pCO_2$ are all available from independent high-quality measurements. As underway $pCO_2$ measurements available in the SOCAT database are all from the surface ocean, it was not possible to assign measured $pCO_2$ values to samples at depth. Note that we discarded data from two cruises (EXPOCODES 33AT20120419 and 49NZ20010828) for reasons explained in the Supplement, ultimately using data from 948 samples for this analysis. Samples within this dataset were taken between 1993 and 2012 over 26 different research cruises. These samples were taken at the ocean surface, always in the top 5 m. They cover a range of practical salinities from 30.73 to 37.57 and temperatures from -1.67 to 31.80 °C at locations shown in Fig. 1.

## 2.3 Iterative methods and underlying assumptions

To the 948 samples for which independent, high-quality measurements of pH, DIC, TA, [DSi], [SRP], and $pCO_2$ are available, we would preferably directly apply Eqs. (4) and (5). This was, however, not possible given the inter-dependence of pH and CA, both of which are necessary to compute $K_1^*$ and $K_2^*$ and, in turn, other carbonate system parameters. Rather than estimating the temperature dependence of pH from $\Delta\text{pH}/\Delta T$ as done by Millero et al. (2002), we used a novel iterative fitting procedure. This procedure is based on an initial estimate of both pH and CA using the Lueker et al. (2000) constants, followed by a recomputation at each iteration using the values of $K_1^*$ and $K_2^*$ from the previous iteration. The calculations were executed in R (R Core Team, 2019) and detailed below.

The code and data files can be downloaded from an online repository (https://doi.org/10.5281/zenodo.3725889, Sulpis et al., 2020).

Since the objective of the present study is to obtain independent measurements of $K_1^*$ and $K_2^*$, we could not directly use GLODAPv2 in situ pH data because the majority of these data were obtained on board by potentiometric or spectrophotometric methods at an equilibrium temperature often higher than the surface seawater (usually 25 °C, occasionally 20 or 13 °C). In addition, pH data were not always delivered to GLODAPv2 on the total pH scale. Consequently, the in situ pH measurements available in GLODAPv2 are all recalculated using measured TA and the Lueker et al. (2000) stoichiometric dissociation constants and converted to the total pH scale if necessary (Olsen et al., 2016). To obtain the pH values delivered to GLODAPv2, which should be independent of $K_1^*$ and $K_2^*$, we converted the GLODAPv2 pH values back to their measurement temperatures (labT) and pH scales, as recorded in the cruise reports. For this reconciliation, we used the bias-corrected TA value from GLODAPv2 rather than the measured TA values that were used in the GLODAPv2 conversion process. Bias correction of TA was done through crossover and inversion analysis of the data (Olsen et al., 2016); for the 26 research cruises we selected, bias correction resulted in TA adjustments of -1 to 10 µmol kg⁻¹. These adjustments however affected the recalculated pH values by less than 0.0001. We then converted all recalculated pH values to the free pH scale ($pH_F^{\text{labT}}$) using the default settings of the *pHconv* function in the seacarb R package (Gattuso et al., 2019). The free pH scale was used during the fitting procedure to avoid further complications with the sulfate and fluoride acid-base systems. Nevertheless, final results are presented on the total pH scale.

Carbonate alkalinity was not directly measured but is instead a pH-dependent quantity computed from TA; see Eq. (3). As a first approximation, we calculated CA from the measured TA by subtracting the contributions of the carbonate, silicate, and phosphate acid-base systems, as well as the auto-dissociation of water, using [SRP], [DSi], and the in situ pH from GLODAPv2. We estimated the total boron concentration from salinity using the Uppström (1974) relationship and calculated its acid-base speciation using the equilibrium constants of Dickson (1990). For the silicate and phosphate acid-base speciation, the equilibrium constants of Yao and Millero (1995) were used. All of these expressions are only valid for temperatures above 0 °C; thus, extrapolation to lower temperatures yields an additional uncertainty to the method. All equilibrium constants were corrected for pressure following Millero (1995), but given that all samples were taken at depth shallower than 5 m depth, this correction is negligible. Using $H_F^{\text{labT}}$, the proton concentration computed from $pH_F^{\text{labT}}$, we also calculated carbonate alkalinity at the temperature of pH measurements ($CA^{\text{labT}}$). This variable was used during the iteration procedure.
Since there are two pH-independent parameters (DIC, pCO₂), we can use these two parameters and one pH-dependent parameter (either pH itself, or CA) to initialize the iterative procedure. This implies that either K₁* or K₂* must be assigned an initial value before starting the iterations. Here, we initially set the in situ K₂* to the value calculated from the Lueker et al. (2000) expressions. The alternative case in which in situ K₁* was initially set to the Lueker et al. (2000) value is described in the Supplement and would have no appreciable impact on the results presented here; i.e. whether the first or the second dissociation constant is assigned an initial value does not affect the results. Each iteration consisted of four different steps:

1. First, K₁* was computed from in situ DIC, CA, pCO₂, and K₂* from Lueker et al. (2000) using Eq. (4) and a Newton-Raphson technique (function *uniroot.all* from package rootSolve; Soetaert and Herman, 2009). These calculated K₁* values were subsequently fitted to a general expression as a function of temperature and salinity of the form:

$$pK_{1\text{ or }2}^* = a_1 + a_2 S_p + a_3 S_p^2 + \frac{a_4}{T} + a_5 \ln(T), \quad (6)$$

where pK₁* corresponds to -log₁₀(Kᵢ) and aᵢ are fitting coefficients determined using the Levenberg-Marquardt algorithm for non-linear least-squares estimates (function *nlsLM* from the minpack.lm package; Elzhov et al., 2016). This expression is of a similar form as Lueker et al. (2000), to facilitate the comparison. Because the salinity range in the sub-dataset where four carbonate system variables are available is narrow (30.73 to 37.57), it was not possible to obtain converging iterations where all the coefficients in Eq. (6) were resolved. Thus, we kept a₂ and a₃ fixed to the Lueker et al. (2000) values, assuming that the salinity dependence of K₁* and K₂* is correct for the salinity range of our dataset, and only solved for a₁, a₄, and a₅.

2. Second, this new expression for K₁*, as well as CA and the expression for K₂* used in step 1, were used to compute pH at in situ temperature. For this, both K₁* and K₂* were calculated at the temperature of pH measurement (K₁*, labT and K₂*, labT). These were used together with the free proton concentration at lab temperature (Hₐ<sup>labT</sup>) and the calculated carbonate alkalinity (CA<sup>labT</sup>), both of which do not change during the iterative procedure, to calculate [H⁺], the free proton concentration at in situ temperature. We expressed DIC as a function of CA, [H⁺], K₁* and K₂*, and assumed that the value of DIC is independent of temperature. Thus,

$$CA \frac{[H^+]^2 + [H^+] K_1^* + K_1^* K_2^*}{[H^+] K_1^* + 2K_1^* K_2^*} = CA^{\text{labT}} \frac{(H_F^{\text{labT}})^2 + H_F^{\text{labT}} K_1^{*, \text{labT}} + K_1^{*, \text{labT}} K_2^{*, \text{labT}}}{H_F^{\text{labT}} K_1^{*, \text{labT}} + 2K_1^{*, \text{labT}} K_2^{*, \text{labT}}}. \quad (7)$$

This equation was rewritten into a quadratic equation, solved analytically for [H⁺], and converted to pH.

3. Third, CA – which is dependent on pH – was updated based on the new [H⁺], as per Eq. (3) and the method outlined for the initial calculation of CA.

4. Fourth, we used Eq. (5) to calculate K₂* as a function of pCO₂, DIC, the new pH, and CA and fit these in situ computed constants to an equation of the form of Eq. (6).

These four steps were repeated, and at each iteration, K₂*, CA, and pH from the previous iteration, were used as initial values. Note that this method assumes that the uncertainty in K₀ is minor compared to that in K₁* and K₂*. We also assumed that no acid-base systems other than the carbonate, borate, silicate, and phosphate acid-base systems contributed to TA – this point will be elucidated later – and that uncertainties in the calculated contributions of the latter three acid-base systems to TA were also minor compared to the uncertainties in K₁* and K₂*.

## 2.4 Uncertainty propagation

The overall uncertainty on the final K₁* and K₂* values is a combination of the uncertainties associated with measurement errors (hereafter termed “analytical uncertainty”) and the uncertainties resulting from the fitting procedures (hereafter termed “fitting uncertainty”) that are propagated throughout the iterations. The analytical uncertainty (σK<sup>ana.</sup>) was computed using the predefined accuracy limits (here, for simplicity, denoted σ) used for the GLODAPv2 secondary quality control procedures. This accuracy limit reflects the minimum bias that can be detected with reasonable certainty (Tanhua et al., 2010) and is based on an objective analysis of systematic biases in ship-based data. Within the GLODAP context the accuracy limit should be interpreted as “the range within which we can realistically expect measurements from the deep ocean to be reproducible”. For each variable the corresponding value is taken from Table 2 in Olsen et al. (2016); i.e. σSₚ = 0.005, σ[DSi] = 2%, σ[SRP] = 2%, σDIC = 4 µmol kg⁻¹, σTA = 6 µmol kg⁻¹. σpH is set to 0.01 following Table 3 in Olsen et al. (2019) and σpCO₂ is set to 2 µatm, corresponding to the minimum accuracy of SOCAT quality control flags A or B. While referred to as accuracy, this number is actually a measure of overall measurement uncertainty and includes uncertainties due to environmental factors (Pierrot et al., 2009). σCA was computed as the square root of the sum of the squares of σCO₂⁻ and σHCO₃⁻. In turn, these were computed using TA, pH, [SRP], [DSi], P, T, and Sₚ as input variables, as well as their respective aforementioned uncertainties, using the error propagation code of Orr et al. (2018). The analytical uncertainty on both K₁* and K₂* was then estimated following the standard
rules of error propagation, as per the following equations:

$$ \sigma K_1^{\text{ana.}} = K_1^* \sqrt{\left(\frac{\sigma[H^+]}{[H^+]}\right)^2 + \left(\frac{\sigma f CO_2}{f CO_2}\right)^2 + \left(\frac{\sqrt{(2\sigma\text{DIC})^2 + (\sigma\text{CA})^2 + (2K_0\sigma f CO_2)^2}}{2\text{DIC}-\text{CA}-2K_0 f CO_2}\right)^2} \quad (8) $$

$$ \sigma K_2^{\text{ana.}} = K_2^* \sqrt{\left(\frac{\sigma[H^+]}{[H^+]}\right)^2 + \left(\frac{\sqrt{(\sigma\text{CA})^2 + (\sigma\text{DIC})^2 + (K_0\sigma f CO_2)^2}}{\text{CA}-\text{DIC}+K_0 f CO_2}\right)^2 + \left(\frac{\sqrt{(2\sigma\text{DIC})^2 + (\sigma\text{CA})^2 + (2K_0\sigma f CO_2)^2}}{2\text{DIC}-\text{CA}-2K_0 f CO_2}\right)^2} \quad (9) $$

The fitting uncertainty ($\sigma K^{\text{fit.}}$) was obtained using a Monte Carlo simulation technique that propagates errors in the fitting coefficients to the predicted K values. At the end of the iterations, the non-linear least-square model fits obtained with the nlsLM function were used as an input in the predictNLS function, from the propagate R package (Spiess, 2018), to calculate $\sigma K^{\text{fit.}}$, neglecting any error in the temperature measurements. The overall uncertainty on $K_1^*$ and $K_2^*$ was then assumed to be the square root of the sum of the squares of the analytical and fitting uncertainties. The 95% confidence intervals for each of the fitting coefficients, i.e. $a_i$ in Eq. (6), shown in Table 1, were extracted from the result of the non-linear least-squares model fits in R using the summary function. Note that, because we did not solve for the salinity coefficients in Eq. (6) due to the limited salinity range of the four carbonate system variable datasets, $a_2$ and $a_3$ are set to the Lueker et al. (2000) values and no confidence interval is computed for these coefficients.

### 3 Results

Using all the data from which T, S<sub>p</sub>, DIC, pH, CA, and pCO<sub>2</sub> are available as high-quality, independent measurements, we were able to derive expressions for $K_1^*$ and $K_2^*$ with a new temperature dependence. The coefficients $a_i$ for $pK_1^*$ and $pK_2^*$, in an equation of the form of Eq. (6), and after 30 iterations, are reported in Table 1, along with their respective 95% confidence intervals. In both expressions, all coefficients are significantly different from zero ($p$ values < 0.001). As CA and pH are being updated throughout the iterations, CA, pH, $K_1^*$, and $K_2^*$ all evolved until the 10th iteration, before converging to a final value (Fig. 2). Thus, we stopped the iterative process after 30 iterations. Between the initial GLODAPv2 value and the 30th iteration, [H<sup>+</sup>] values vary by up to 6.6%. CA values are only weakly affected by these in situ pH updates as they change by a maximum of 0.2% throughout the first 30 iterations. Largest pH and CA changes occur at the colder end of the temperature range. $pK_1^*$ and $pK_2^*$ values both shift upward throughout the iterations, the $pK_1^*$ increase being higher than that of $pK_2^*$, especially in cold waters (Fig. 2).

As shown in Fig. 3, the $pK^*$ values obtained with the iterative procedure are statistically indistinguishable from the $pK^*$ values of Lueker et

al. (2000) (i.e. $pK^{\text{this study}} - \sigma pK^{\text{this study}} < pK^{\text{Lueker}} < pK^{\text{this study}} + \sigma pK^{\text{this study}}$) over most of the temperature range. Nevertheless, in cold waters, below temperatures 8.1 and 9.2 °C respectively, $pK_1^*$ and $pK_2^*$ are significantly higher than the values reported by Lueker et al. (2000). In Fig. 3, we also provide a comparison of the $pK^*$ values of this study with those of Papadimitriou et al. (2018), as the latter study focuses on low-temperature waters. As shown in Fig. 3c and d, both $pK_1^*$ and $pK_2^*$ from Papadimitriou et al. (2018) are slightly lower than the fitted values from this study, for a similar salinity, and are also lower than the Lueker et al. (2000) $pK^*$ values except for temperatures below 4.0 and 1.2 °C, respectively. Thus, the fact that the Lueker et al. (2000) study underestimates both $pK_1^*$ and $pK_2^*$ in waters near freezing point seems to be a consistent feature across studies.

### 4 Discussion

Using underestimated $pK_1^*$ and $pK_2^*$ values implies that, for a given state, computed [H<sub>2</sub>CO<sub>3</sub><sup>*</sup>] or pCO<sub>2</sub> would be underestimated and [CO<sub>3</sub><sup>2-</sup>] overestimated. This potentially has strong implications for our representation of seawater carbonate chemistry in low-temperature marine environments, such as polar regions. Hence, we highlight the implications of this work for the estimation of two carbonate system variables in polar regions, i.e. pCO<sub>2</sub> and the saturation state of seawater with respect to calcite ($\Omega_{\text{Ca}}$). But, first, we examine error propagation and the dependence of $pK_1^*$ and $pK_2^*$ to salinity and discuss the influence of organic alkalinity and the quality of pH measurements on the results presented here.

#### 4.1 Influence of pH and TA measurement quality

The pH of samples that are used to derive the $K^*$ fits presented here was measured using electrodes or spectrophotometrically, between 1993 and 2012. During this period, it was shown that impurities present in commercially available dyes could generate a systematic bias in the measured pH (Yao et al., 2007). Recently, Carter et al. (2018) pointed to systematic discrepancies resulting from differing approaches to pH measurements. Yao et al. (2007) noted that these impurities can contribute to pH offsets as large as 0.01 pH units, which corresponds to the analytical uncertainty in pH ($\sigma$ pH) that we use here, taken from Olsen et al. (2019). We must therefore investigate whether the fact that most pH measurements in the pH dataset are not from spectrophotometric measurements with purified dyes could alter the conclusion of underestimated $pK^*$ values in cold waters. To answer this question, we gathered a sub-dataset of the more recent GLODAPv2.2019 data product (Olsen et al., 2019), composed of samples from nine different cruises (EXPOCODES 320620140320, 06AQ20150817, 33AT20120324, 33AT20120419, 33HQ20150809,
<table><caption>**Table 1.** Comparison of the coefficients for $pK_1^*$ and $pK_2^*$ between this study and Lueker et al. (2000), using an equation of the form of Eq. (6). Coefficients are given as value ±95 % confidence interval. In bold are the coefficients ($a_1$, $a_4$ and $a_5$) that were solved for, while the two salinity coefficients ($a_2$ and $a_3$) were fixed to the Lueker et al. (2000) values.</caption><thead><tr><th rowspan="2"></th><th colspan="2">This study</th><th colspan="2">Lueker et al. (2000)</th></tr><tr><th>pK<sub>1</sub><sup>*</sup></th><th>pK<sub>2</sub><sup>*</sup></th><th>pK<sub>1</sub><sup>*</sup></th><th>pK<sub>2</sub><sup>*</sup></th></tr></thead><tbody><tr><td>a<sub>1</sub></td><td>-172.4493±26.131</td><td>-59.4636±24.016</td><td>-61.2172</td><td>25.9290</td></tr><tr><td>a<sub>2</sub></td><td>-0.011555</td><td>-0.01781</td><td>-0.011555</td><td>-0.01781</td></tr><tr><td>a<sub>3</sub></td><td>0.0001152</td><td>0.0001122</td><td>0.0001152</td><td>0.0001122</td></tr><tr><td>a<sub>4</sub></td><td>8510.63±1139.8</td><td>4226.23±1050.8</td><td>3633.86</td><td>471.78</td></tr><tr><td>a<sub>5</sub></td><td>26.32996±3.9161</td><td>9.60817±3.5966</td><td>9.67770</td><td>-3.16967</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Evolution of pH<sub>T</sub>, CA, $pK_1^*$, and $pK_2^*$ as a function of in situ temperature and iterations. (a) Differences between recomputed pH<sub>T</sub> and pH<sub>T</sub> from GLODAPv2. (b) Differences between recomputed CA and CA estimated from GLODAPv2 data and the Lueker et al. (2000) constants. Differences between recomputed (c) $pK_1^*$ and (d) $pK_2^*$ and the Lueker et al. (2000) values at in situ temperature and at a practical salinity of 35.</figcaption></figure>

33RO20150410, 33RO20150525, 33RO20161119, and 33RR20160208) for which T, P, S<sub>p</sub>, DIC, TA, [SRP], and [DSi] are available and associated with a WOCE flag of 2, for which pH was measured spectrophotometrically using purified dyes only, and for which an associated SOCAT pCO<sub>2</sub> value is available. Although this independent dataset is too small to apply the iterative procedure and obtain acceptable pK<sup>*</sup> fits, we can use it to compare the $K_1^*$ and $K_2^*$ values obtained with Eqs. (4) and (5) and this purified-dye independent dataset to the $K_1^*$ and $K_2^*$ values obtained with Eqs. (4) and (5) and the regular dataset.

For both the purified-dye independent dataset and the regular dataset, samples were sorted according to in situ temperature and grouped into bins of 0.5 °C. Temperature bins containing a single sample were not used. For each tempera-

ture bin, the mean pK<sup>*</sup> values, obtained with Eqs. (4) and (5), along with the associated uncertainties, were computed. Plotting the difference between the pK<sup>*</sup> values computed using the regular dataset and those computed using the purified-dye independent dataset as a function of seawater in situ temperature (Fig. 4), we do not see any clear systematic bias caused by the use of purified dyes. This means that pK<sup>*</sup> values computed from a dataset with purified-dye pH measurements only are not higher or lower than pK<sup>*</sup> values computed from a dataset with pH measured using primarily impure dye. More importantly, in colder waters (T < ~ 2 °C), the differences between the pK<sup>*</sup> values from this study and those from Lueker et al. (2000) (black line in Fig. 4) are larger than what can be explained by the choice of dye for spectrophotometric pH measurements. Thus, the use of impure vs. purified
<figure><img src="image_6.png" /><figcaption>**Figure 3.** (a) pK<sub>1</sub><sup>*</sup> and (b) pK<sub>2</sub><sup>*</sup> as a function of temperature where the colour represents practical salinity (S<sub>p</sub>) and the fits are fixed for a S<sub>p</sub> of 35. Comparison of (c) pK<sub>1</sub><sup>*</sup> and (d) pK<sub>2</sub><sup>*</sup> as a function of temperature from this study (red lines), Lueker et al. (2000, blue line) and Papadimitriou et al. (2018, black line). The solid blue line represents the pK<sup>*</sup> fits from Lueker et al. (2000), the solid red line the pK<sup>*</sup> from this study computed with the coefficients presented in Table 1. Dashed red lines are overall uncertainties as defined in Sect. 2.4.</figcaption></figure>

dye in pH measurements should not affect the conclusions presented here.

Another issue in GLODAPv2 carbonate system measure-
ments may be the fact that some seawater samples contain
measurable amounts of organic bases (Fong and Dickson,
2019; Patsavas et al., 2015; Yang et al., 2015). This organic
alkalinity is unaccounted for in the definition of total alka-
linity of Eq. (3), thus causing biased, overestimated, com-
puted carbonate alkalinity values. This does not only con-
cern coastal waters but also open-ocean waters, where the
total concentration of these organic bases could be in the or-
der of a few micromols per kilogram (Fong and Dickson,
2019). This, however, should not substantially alter the re-
sults presented here, due to the small amount of these organic

bases and consequently small impact on computed pK* val-
ues. If anything, subtracting the contributions of these un-
accounted bases to the total alkalinity measurements would
have a unidirectional effect on the dissociation constant esti-
mates, shifting the pK* values upwards – see Eqs. (4) and (5)
– further away from the Lueker et al. (2000) values in cold
waters.

## 4.2 Uncertainties in carbonate system calculations

The relative overall uncertainties (σK*/K*) were ~2.5% for both K₁* and K₂*. In both cases, the analytical uncertainty (σKana.) was more than twice as high as the fitting uncertainty (σKft). The overall pH measurement uncertainty of
<figure><img src="image_7.png" /><figcaption>**Figure 4.** (a) pK\*₁ and (b) pK\*₂ differences between those obtained with the regular dataset and those obtained with the purified-dye independent dataset, as a function of in situ seawater temperature. Orange circles represent the mean computed pK\* values within a temperature bin, and the vertical black bars stand for the associated uncertainties deviations. The solid black line is the difference between the pK\* fit from this study and that from Lueker et al. (2000).</figcaption></figure>

the GLODAP dataset ($\sigma$ pH = 0.01) is relatively high, causing the pH term in Eqs. (8) and (9) to be the dominant factor in $\sigma$ K<sup>ana.</sup>. Thus, what dominates the overall uncertainties for K\* estimates from this study is the uncertainty in pH, which explains the fact that the relative uncertainties of K\*₁ and K\*₂ are similar. Converting from $\sigma$ K\* to $\sigma$ pK\*, we find $\sigma$ pK\*₁=$\sigma$ pK\*₂ = 0.011. The overall uncertainty for pK\*₁ is higher than that reported by Orr et al. (2018), but the overall uncertainty for pK\*₂ is smaller. These values are quite high relative to the uncertainties of previous expressions for K\*₁ and K\*₂ as reported in Table 2 of Millero (2007), which we attribute to the fact that they reflect both the uncertainty from the fits and from the measurements that we use.

Using the pK\* values from Table 1, setting the analytical uncertainties for each variable to the values reported in Sect. 2.4, and using $\sigma pK\*₁=\sigma pK\*₂ = 0.011$, we use the ExCEL version of CO2SYS from Orr et al. (2018) and analyse the propagation of uncertainties on two computed variables, pCO₂ and the saturation state of seawater with respect to calcite, ΩCa (Mucci, 1983; see discussion in Sect. 4.4). For this purpose, we use all data points from GLODAPv2 that contain T, P, S<sub>p</sub>, DIC, TA, pH, [SRP], and [DSi], in the top 10 m of the water column. The quality criterion remains unchanged; i.e. we use only data associated with a WOCE flag of 2. The obtained dataset contains 3392 samples, including the 948 samples of the regular dataset, and covers a salinity range from 3.46 to 37.57 and a temperature range from -1.91 to 31.80 °C. Depending on which carbonate system pair of variables is used, both the magnitude and the uncertainties of computed variables can differ (Orr et al., 2018; Ribas-Ribas et al., 2014). Here, we use three different pairs of variables, i.e. TA–DIC, TA–pH, and DIC–pH, to compute pCO₂ and ΩCa and their associated propagated uncertainties, $\sigma$ pCO₂ and $\sigma$ ΩCa, respectively.

Relative uncertainties generated with the TA–DIC pair appear to be particularly sensitive to salinity (Fig. 5), increasing from ~5% to ~15% for both pCO₂ and ΩCa as the salin-

ity decreases from 35 to 5. The overall relative uncertainty for both pCO₂ and ΩCa is less dependent on salinity when pH is used as an input variable with either DIC or TA. For both the TA–pH and the DIC–pH pairs, the overall relative uncertainty on pCO₂ increases with an increasing temperature, while the overall relative uncertainty on ΩCa decreases with an increasing temperature (Fig. 5).

As depicted in Fig. 6, in agreement with Raimondi et al. (2019), we conclude that the DIC–pH pair offers the lowest overall relative uncertainty for computed pCO₂ over the range of salinities and temperatures investigated. Using the DIC–pH pair also has the important benefit of not having to make any assumption regarding organic alkalinity or the boron-to-salinity ratio (Fong and Dickson, 2019). Conversely, the TA–pH pair is the one generating the highest overall relative uncertainties on computed pCO₂. As for ΩCa, the TA–pH pair provides the lowest overall relative uncertainty below a temperature of ~20°C, whereas the TA–DIC pair should be preferred in warmer waters (Fig. 6).

## 4.3 Implications for surface ocean pCO₂

To evaluate the implications of the revised temperature dependence of the carbonic acid dissociation constants, we compare ocean carbonate chemistry as calculated with the constants from this study, those of Lueker et al. (2000), and those of Millero et al. (2002). Whereas the constants from Lueker et al. (2000) are the most commonly used by the oceanographic community, as recommended by Dickson et al. (2007), the constants from Millero et al. (2002) were derived in an approach similar to that presented here, using a large dataset of in situ measurements. Thus, it appears relevant to include a comparison with Millero et al. (2002) in the present discussion. The major differences between our approach and the approach of Millero et al. (2002) are the calculation of CA from measured TA and pH (iteratively versus direct) and the conversion of pH measurements (iteratively versus estimating $\Delta$pH/$\Delta$T from the constants of Mehrbach et al., 1973).

For this comparison, we use data from the Southern Ocean Carbon and Climate Observations and Modelling (SOCOM) project (https://soccom.princeton.edu/, last access: 5 March 2020). The SOCCOM project has deployed more than 100 Argo floats equipped with biogeochemical sensors in the Southern Ocean. These sensors include pH, and SOCCOM routinely calculates the full carbon chemistry (including pCO₂ and ΩCa) using a combination of measured T, S<sub>p</sub>, pH, O₂, and empirical algorithms for TA (Carter et al., 2018). The SOCCOM data used here, both measured and calculated, were downloaded as a Matlab file from https://library.ucsd.edu/dc/object/bb0515927k (last access: 5 March 2020).

The method used to calculate pCO₂ is detailed in the data file (within the *FloatViz* structure). Briefly, the SOCCOM project uses the Lueker et al. (2000) K\*₁ and K\*₂, Perez and Fraga (1987) for K<sub>F</sub>, Dickson (1990) for KSO₄,
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Overall relative uncertainty for (a, b, c) pCO₂ and (d, e, f) ΩCa as a function of in situ temperature (horizontal axis) and salinity (colour bar). Panels (a, d) are for the TA–DIC pair, panels (b, e) are for the TA–pH pair, and panels (c, f) are for the DIC–pH pair. This is computed using all data points from GLODAPv2 that contain T, P, S<sub>p</sub>, DIC, TA, pH, [SRP], and [DSi], in the top 10 m of the water column (3392 samples).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Pair of carbonate system variables (TA–DIC, TA–pH, or DIC–pH) providing the lowest overall relative uncertainty, as a function of in situ temperature and practical salinity, for pCO₂ and ΩCa. This is computed using all data points from GLODAPv2 that contain T, P, S<sub>p</sub>, DIC, TA, pH, [SRP], and [DSi] in the top 10 m of the water column (3392 samples).</figcaption></figure>

and Lee et al. (2010) for total boron estimates. Both [DSi] and [SRP] are estimated from the measured nitrate concentration using stoichiometric ratios of 2.5 and 1/16, respectively. We applied the same method but substituted the Lueker et al. (2000) K<sub>1</sub><sup>*</sup> and K<sub>2</sub><sup>*</sup> constants with either the constants from this study or the Millero et al. (2002) constants. We were then able to compare surface ocean (defined as the upper 10 m of the water column) pCO₂ obtained using Lueker et al. (2000) or Millero et al. (2002) with pCO₂ obtained using the constants from this study. The analytical uncertainties were set to σS<sub>p</sub> = 0.005 and σT = 0.005 °C (Olsen et al., 2016), σ[DSi] = 0.9 µmol kg<sup>-1</sup> and σ[SRP] = 0.5 µmol kg<sup>-1</sup> (combination of uncertainty in nitrate concentration from Argo data, i.e. 0.5 µmol kg<sup>-1</sup> as given in Johnson et al. (2017) and a 30% uncertainty in stoichiometric ratios), σTA = 5.6 µmol kg<sup>-1</sup> (Carter et al., 2018), and σpH = 0.005 (Johnson et al., 2017). σpK<sub>1</sub><sup>*</sup> and σpK<sub>2</sub><sup>*</sup> were set to 0.011, respectively, when the constants from this study were used. For both Lueker et al. (2000) or Millero et al. (2002), they were set to the default values given by Orr et al. (2018), i.e. σpK<sub>1</sub><sup>*</sup> = 0.0075 and σpK<sub>2</sub><sup>*</sup> = 0.015. Uncertainties on the computed pCO₂ were propagated using the Matlab version
of the Orr et al. (2018) CO2SYS software with error propa-
gation.

pCO₂ values obtained with the constants derived from this study are clearly higher than the Lueker et al. (2000)-based values in the southernmost regions, where temperatures are lowest (Fig. 7a, b), with a maximum difference (ΔpCO₂=pCO₂Lueker−pCO₂this study) of −55 ± 17 µatm when the surface ocean is near the freezing point. The uncertainty on ΔpCO₂ (σΔpCO₂, *grey* lines in Fig. 7) is computed as the square root of the sum of the squares of σpCO₂Lueker and σpCO₂this study. Given the large uncertainties, the pCO₂ difference between values based on constants derived from this study and values based on Lueker et al. (2000) is only statistically significant (i.e. ΔpCO₂ + σΔpCO₂this study < 0) for temperatures below ~8 °C (Fig. 7a). pCO₂ values obtained using Millero et al. (2002) constants appear to be midway between pCO₂ values obtained using Lueker et al. (2000) constants and pCO₂ values based on the constants from this study (Fig. 7c, d).

Recently the SOCCOM Argo array was used to re-evaluate
the Southern Ocean carbon sink (Gray et al., 2018). Tra-
ditional ship-based observations indicate a strong CO₂ up-
take in the entire Southern Ocean, but these observations are
known to have a strong seasonal bias (Bakker et al., 2016),
as well as a smaller spatial bias due to many areas being
severely undersampled (Takahashi et al., 2012). Using pCO₂
calculated by the above method, Gray et al. (2018) showed
that the Southern Ocean CO₂ uptake is considerably smaller
than previously estimated. In parallel, Bailey et al. (2018)
showed that the CO₂ solubility constant from Weiss (1974)
used in the majority of studies, including this one, was un-
derestimated in waters below 0°C, which implies that sur-
face pCO₂ is underestimated. In this study, using the new
constants in Table 1, the computed Southern Ocean pCO₂
is also higher than when computed using the constants of
Lueker et al. (2000) or the constants of Millero et al. (2002),
as shown in Fig. 7. The Southern Ocean is a net CO₂ sink
because the pCO₂ in surface waters is on average lower than
the atmospheric pCO₂. If the surface-water pCO₂ is revised
upward, the resulting flux of CO₂ from the atmosphere to the
surface waters becomes smaller. Thus, results from Gray et
al. (2018), Bailey et al. (2018), and the present study all advo-
cate for a weaker CO₂ sink in the Southern Ocean. The ocean
CO₂ sink is immensely important and currently estimated to
remove ~25% of anthropogenic CO₂ emissions (Le Quéré
et al., 2018). If the CO₂ uptake by the Southern Ocean is
much smaller than previously estimated, there must be miss-
ing sinks elsewhere in the Earth system, be it in the oceanic
or terrestrial realm. This highlights the need for a better un-
derstanding of the dynamics of the ocean carbon sink, in-
cluding its regional and temporal variability. To validate our
results, the high uncertainties associated with stoichiometric
constants (Orr et al., 2018), coupled to the low spatial and
temporal resolution of measurements in high latitudes, need

to be addressed. Whether in the laboratory or in the field, fu-
ture work should focus on a better understanding of seawater
carbonate chemistry in cold waters.

4.4 Implications for calcium carbonate chemistry

In seawater undersaturated with respect to calcite or arago-
nite, (ΩCa or ΩAr < 1), the CaCO3 phase of interest should
dissolve if present. ΩCa depends on the Ca²⁺ concentration
in seawater (a function of salinity and therefore nearly invari-
ant at depth), the solubility product of calcite (Mucci, 1983),
and the CO₃²⁻ concentration in seawater. Because of the lat-
ter, computed ΩCa values are impacted by the choice of car-
bonic acid dissociation constants. Note that, as reported in
Orr et al. (2018), the relative uncertainty on the solubility
product of calcite is about 5%. Below, we test the implica-
tions of the K⁺₁ and K⁺₂ values from this study on predictions
of calcite saturation state in seawater.

Naviaux et al. (2019) recently observed discrepancies be-
tween ΩCa computed with the TA–DIC pair and ΩCa com-
puted with the TA–pH pair, which they attributed to the in-
ternal inconsistency of the carbonate system, i.e. the fact that
measured pH does not correspond to calculated pH. Instead,
or in addition, the calcite saturation depth calculated by Navi-
aux et al. (2019) could be erroneously too shallow due to an
overestimated K⁺₂ and, consequently, overestimated seawater
[CO₃²⁻] and ΩCa.

Here, we used data from a cruise (33RO20071215, GLO-
DAPv2 cruise no. 345) along the CLIVAR repeat section
P18 that took place in 2007, following a latitudinal transect
in the south-eastern Pacific Ocean, in which the carbonate
chemistry variables DIC, TA, and pH were measured (see
Olsen et al., 2016, for details about the data). All calcula-
tions were carried out using the discrete data, but for pur-
poses of visualization (Fig. 8) we used a nearest-neighbour
interpolation (function *griddata* in Matlab). In Fig. 8, we
compare ΩCa as computed using K<sup>*</sup><sub>1</sub> and K<sup>*</sup><sub>2</sub> from this study,
with ΩCa based on Lueker et al. (2000) constants. We also
compare ΩCa computed from TA–DIC with ΩCa computed
from TA–pH. From Fig. 8, it can be seen that the ΩCa dif-
ference between two different carbonate system pairs (TA–
DIC, TA–pH) is ~5 times smaller than the ΩCa difference
that is due to the set of dissociation constants. Thus, the ap-
parent dissolution observed by Naviaux et al. (2019) may
be explained by overestimated dissociation constants atop
inconsistencies arising from the choice of carbonate vari-
ables used in the calculations. We also note, based on Fig. 8,
that ΩCa overestimation is largest in the southernmost part
of the Pacific surface waters, where the temperature is the
lowest. Nevertheless, the maximum calculated ΩCa differ-
ences, i.e. ΔΩCa = 0.06 with TA–DIC and ΔΩCa = 0.07
with TA–pH, are 2–3 times lower than the average combined
uncertainty, i.e. ((σΩ<sup>Lueker</sup><sub>Ca</sub>)² + (σΩ<sup>this study</sup><sub>Ca</sub>)²)⁰.⁵ = 0.20 and
0.17, respectively. These high uncertainties are attributed to
the high measurement uncertainties that we use (those from
<figure><img src="image_10.png" /></figure>

**Figure 7.** Difference in surface *pCO*₂ obtained from pH, temperature, practical salinity, and dissolved oxygen measured by the SOCCOM Argo array using (**a**, **b**) the Lueker et al. (2000) constants and the constants in Table 1 or (**c**, **d**) the Millero et al. (2002) constants and the constants in Table 1. Plots (**a**, **c**) represents the *pCO*₂ difference (Δ*pCO*₂) as a function of in situ temperature. Solid black lines are the mean Δ*pCO*₂, while solid grey lines are plus or minus σ Δ*pCO*₂ (square root of the sum of the squares of σ*pCO*₂<sub>Lueker</sub> and σ*pCO*₂<sub>this study</sub>). Maps (**b**, **d**) depict the spatial distribution of Δ*pCO*₂ in the Southern Ocean, where each point corresponds to an Argo float measurement.

<figure><img src="image_11.png" /></figure>

**Figure 8.** Difference between (**a**) the calcite saturation state (Ω<sub>Ca</sub>) computed from TA and DIC using the dissociation constants of Lueker et al. (2000) and the constants derived from this study, (**b**) Ω<sub>Ca</sub> computed from pH and TA using the dissociation constants of Lueker et al. (2000) and the constants derived from this study, and (**c**) Ω<sub>Ca</sub> computed from the DIC and TA pair and the pH and TA pair and the dissociation constants derived from this study.
Olsen et al., 2016, see Sect. 2.4), the overall uncertainty on the dissociation constants from this study, the uncertainty on Lueker et al. (2000) constants, and especially to the uncertainty in the calcite solubility product. As noted by Orr et al. (2018), the uncertainty on the calcite solubility product (~5%) causes the total uncertainty on ΩCa to be considerably larger than the uncertainty in seawater [CO₃²⁻]. In this study, only the K\*₁ and K\*₂ constants have been re-evaluated, but the large overall uncertainties on calculated saturation states clearly indicate that more work is necessary to define the solubility products of calcium carbonate minerals in the ocean. While beyond the scope of the present study, the results presented here show that proper assessments of present and future ocean acidification are highly sensitive to the present knowledge gaps regarding the thermodynamics of ocean carbon chemistry.

## 5 Conclusion

An iterative procedure allowed us to estimate the temperature dependence of the first and second carbonic acid stoichio-metric dissociation constants (K\*₁ and K\*₂, respectively) from a large dataset of high-quality oceanographic measurements. Both K\*₁ and K\*₂ were similar to the constants of Lueker et al. (2000) that are currently used by most of the oceanographic community, as recommended by Dickson (2007), but the K\*₁ and K\*₂ values were lower in cold seawater, below a temperature of ~8–9 °C. Consequently, at these temperatures, pCO₂ computed using the constants of Lueker et al. (2000) may be underestimated and [CO₃²⁻] overestimated, meaning that the cold oceans are more undersaturated with respect to CaCO₃ minerals than expected. We also used a GLODAP sub-dataset to study the internal consistency of the carbonate system and found that the DIC–pH carbonate system pair provides the smallest overall uncertainty when computing seawater pCO₂. When calculating the saturation state of seawater with respect to calcite, the TA–DIC pair should be used to minimize the overall uncertainty when seawater is warmer than ~20 °C, whereas the TA–pH pair should be preferred below ~20 °C. These results are of critical importance for scientists contemplating studies of high-latitude marine carbonate chemistry and underline that improved knowledge of what causes the CO₂ system inconsistencies in cold waters is key to improve our understanding of the marine carbon budget.

*Data availability.* The R code and data files are available on Zenodo (https://doi.org/10.5281/zenodo.3725889, Sulpis et al., 2020).

*Supplement.* The supplement related to this article is available online at: https://doi.org/10.5194/os-16-847-2020-supplement.

*Author contributions.* OS conceived the original idea, OS and MH designed the research, MH wrote the R code, and SKL advised on the use of the GLODAPv2, SOCAT, and Argo datasets. All authors contributed to paper writing, with OS taking the lead.

*Competing interests.* The authors declare that they have no conflict of interest.

*Acknowledgements.* We thank Dorothee C. E. Bakker and Mariana Ribas-Ribas for fruitful discussions, Alfonso Mucci and Andrew Dickson for helpful comments and suggestions on an earlier version of the paper, the journal editor Mario Hoppema for handling the paper, and two anonymous reviewers for their comments. We also thank Anna de Kluijver for her help in setting up the Monte Carlo simulation. We thank all who contributed to the creation of GLODAPv2 and SOCAT. Olivier Sulpis acknowledges the Department of Earth and Planetary Sciences at McGill University for financial support during his residency in the graduate programme and acknowledges funding from the Dutch Ministry of Education via the Netherlands Earth System Science Centre (NESSC). Siv K. Lauvset acknowledges funding from the Research Council of Norway.

*Financial support.* This research has been supported by Dutch Ministry of Education via the Netherlands Earth System Science Centre (NESSC (grant no. 024.002.001)) and the Research Council of Norway (ICOS-Norway (grant no. 245927) and NorArgo2 (grant no. 269753)).

*Review statement.* This paper was edited by Mario Hoppema and reviewed by two anonymous referees.

## References

Bailey, N., Papakyriakou, T. N., Bartels, C., and Feiyue, W.: Henry’s Law constant for CO₂ in aqueous sodium chloride solutions at 1 atm and sub-zero (Celsius) temperatures, Mar. Chem., 207, 26–32, https://doi.org/10.1016/j.marchem.2018.10.003, 2018.

Bakker, D. C. E., Pfeil, B., Landa, C. S., Metzl, N., O’Brien, K. M., Olsen, A., Smith, K., Cosca, C., Harasawa, S., Jones, S. D., Nakaoka, S., Nojiri, Y., Schuster, U., Steinhoff, T., Sweeney, C., Takahashi, T., Tilbrook, B., Wada, C., Wanninkhof, R., Alin, S. R., Balestrini, C. F., Barbero, L., Bates, N. R., Bianchi, A. A., Bonou, F., Boutin, J., Bozec, Y., Burger, E. F., Cai, W.-J., Castle, R. D., Chen, L., Chierici, M., Currie, K., Evans, W., Featherstone, C., Feely, R. A., Fransson, A., Goyet, C., Greenwood, N., Gregor, L., Hankin, S., Hardman-Mountford, N. J., Harlay, J., Hauck, J., Hoppema, M., Humphreys, M. P., Hunt, C. W., Huss, B., Ibánhez, J. S. P., Johannessen, T., Keeling, R., Kitidis, V., Körtzinger, A., Kozyr, A., Krasakopoulou, E., Kuwata, A., Landschützer, P., Lauvset, S. K., Lefèvre, N., Lo Monaco, C., Manke, A., Mathis, J. T., Merlivat, L., Millero, F. J., Monteiro, P. M. S., Munro, D. R., Murata, A., Newberger, T., Omar, A. M., Ono, T., Paterson, K., Pearce, D., Pierrot, D., Robbins, L. L., Saito, S.,
Salisbury, J., Schlitzer, R., Schneider, B., Schweitzer, R., Sieger, R., Skjelvan, I., Sullivan, K. F., Sutherland, S. C., Sutton, A. J., Tadokoro, K., Telszewski, M., Tuma, M., van Heuven, S. M. A. C., Vandemark, D., Ward, B., Watson, A. J., and Xu, S.: A multi-decade record of high-quality fCO₂ data in version 3 of the Surface Ocean CO₂ Atlas (SOCAT), Earth Syst. Sci. Data, 8, 383–413, https://doi.org/10.5194/essd-8-383-2016, 2016.

Bates, N. R.: Ocean Carbon Cycle, in: Encyclopedia of Ocean Sciences, 3rd edn., edited by: Cochran, J. K., Bokuniewicz, H. J., and Yager, P. L., Elsevier, 1, 418–428, 2019.

Brown, K. A., Miller, L. A., Davelaar, M., Francois, R., and Tortell, P. D.: Over-determination of the carbonate system in natural seawater brine and assessment of carbonic acid dissociation constants under low temperature, high salinity conditions, Mar. Chem., 165, 36–45, https://doi.org/10.1016/j.marchem.2014.07.005, 2014.

Byrne, R. H. and Yao, W.: Procedures for measurement of carbonate ion concentrations in seawater by direct spectrophotometric observations of Pb(II) complexation, Mar. Chem., 112, 128–135, 2008.

Cai, W.-J. and Wang, Y.: The chemistry, fluxes, and sources of carbon dioxide in the estuarine waters of the Satilla and Altamaha Rivers, Georgia, Limnol. Oceanogr., 43, 657–668, 1998.

Carter, B. R., Feely, R. A., Williams, N. L., Dickson, A. G., Fong, M. B., and Takeshita, Y.: Updated methods for global locally interpolated estimation of alkalinity, pH, and nitrate, Limnol. Oceanogr.-Meth, 16, 119–131, 2018.

Dickson, A. G.: An exact definition of total alkalinity and a procedure for the estimation of alkalinity and total inorganic carbon from titration data, Deep-Sea Res., 28A, 609–623, 1981.

Dickson, A. G.: Thermodynamics of the dissociation of boric acid in synthetic seawater from 273.15 to 318.15 K, Deep-Sea Res., 37, 755–766, 1990.

Dickson, A. G. and Millero, F. J.: A comparison of the equilibrium constants for the dissociation of carbonic acid in seawater media, Deep-Sea Res., 34, 1733–1743, 1987.

Dickson, A. G., Sabine, C. L., and Christian, J. R.: Guide to best practices for ocean CO₂ measurement, Sidney, British Columbia, North Pacific Marine Science Organization, PICES Special Publication 3, 191 pp., http://hdl.handle.net/11329/249, 2007.

Dinauer, A. and Mucci, A.: Spatial variability in surface-water pCO₂ and gas exchange in the world's largest semi-enclosed estuarine system: St. Lawrence Estuary (Canada), Biogeosciences, 14, 3221–3237, https://doi.org/10.5194/bg-14-3221-2017, 2017.

Easley, R. A., Patsavas, M. C., Byrne, R. H., Liu, X., Feely, R. A., and Mathis, J. T.: Spectrophotometric measurement of calcium carbonate saturation states in seawater, Environ. Sci. Technol., 47, 1468–1477, 2013.

Elzhov, T. V., Mullen, K. M., Spiess, A.-N., and Bolker, B.: minpack.lm: R Interface to the Levenberg-Marquardt Nonlinear Least-Squares Algorithm Found in MINPACK, Plus Support for Bounds, R package version 1.2-1, available at: https://CRAN.R-project.org/package=minpack.lm (last access: 5 March 2020), 2016.

Fong, M. B. and Dickson, A. G.: Insights from GO-SHIP hydrography data into the thermodynamic consistency of CO₂ system measurements in seawater, Mar. Chem., 211, 52–63, https://doi.org/10.1016/j.marchem.2019.03.006, 2019.

Gattuso, J.-P. and Hansson, L.: Ocean Acidification, Oxford University Press, Oxford, 2011.

Gattuso, J.-P., Epitalon, J.-M., Lavigne, H., and Orr, J.: seacarb: Seawater Carbonate Chemistry. R package version 3.2.12, available at: https://CRAN.R-project.org/package/seacarb (last access: 5 March 2020), 2019.

Goyet, C. and Poisson, A.: New determination of carbonic acid dissociation constants in seawater as a function of temperature and salinity, Deep-Sea Res., 36, 1635–1654, https://doi.org/10.1016/0198-0149(89)90064-2, 1989.

Gray, A. R., Johnson, K. S., Bushinsky, S. M., Riser, S. C., Russell, J. L., Talley, L. D., Wanninkhof, R., Williams, N. L., and Sarmiento, J. L.: Autonomous Biogeochemical Floats Detect Significant Carbon Dioxide Outgassing in the High-Latitude Southern Ocean, Geophys. Res. Lett., 45, 9049–9057, 2018.

Gruber, N., Clement, D., Carter, B. R., Feely, R. A., van Heuven, S., Hoppema, M., Ishii, M., Key, R. M., Kozyr, A., Lauvset, S. K., Lo Monaco, C., Mathis, J. T., Murat, A., Olsen, A., Perez, F., Sabine, C. L., Tanhua, T., and Wanninkhof, R.: The oceanic sink for anthropogenic CO₂ from 1994 to 2007, Science, 363, 1193–1199, 2019.

Hunter, K. A.: The temperature dependence of pH in surface seawater, Deep-Sea Res. Pt. I, 45, 1919–1930, 1998.

Johnson, K. S., Plant, J. N., Coletti, L. J., Jannasch, H. W., Sakamoto, C. M., Riser, S. C., Swift, D. D., Williams, N. L., Boss, E., Haentjens, N., Talley, L. D., and Sarmiento, J. L.: Biogeochemical sensor performance in the SOCCOM profiling float array, J. Geophys. Res.-Oceans, 122, 6416–6436, 2017.

Key, R. M., Olsen, A., van Heuven, S., Lauvset, S. K., Velo, A., Lin, X., Schirnick, C., Kozyr, A., Tanhua, T., Hoppema, M., Jutterström, S., Steinfeldt, R., Jeansson, E., Ishii, M., Perez, F. F., and Suzuki, T.: Global Ocean Data Analysis Project, Version 2 (GLODAPv2), ORNL/CDIAC-162, NDP-P093, Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, US Department of Energy, Oak Ridge, Tennessee, 2015.

Koeve, W. and Oschlies, A.: Potential impact of DOM accumulation on fCO₂ and carbonate ion computations in ocean acidification experiments, Biogeosciences, 9, 3787–3798, https://doi.org/10.5194/bg-9-3787-2012, 2012.

Landschützer, P., Gruber, N., Haumann, F. A., Rödenbeck, C., Bakker, D. C. E., van Heuven, S., Hoppema, M., Metzl, N., Sweeney, C., Takahashi, T., Tilbrook, B., and Wanninkhof, R.: The reinvigoration of the Southern Ocean carbon sink, Science, 349, 1221–1224, 2015.

Lauvset, S. K., Key, R. M., Olsen, A., van Heuven, S., Velo, A., Lin, X., Schirnick, C., Kozyr, A., Tanhua, T., Hoppema, M., Jutterström, S., Steinfeldt, R., Jeansson, E., Ishii, M., Perez, F. F., Suzuki, T., and Watelet, S.: A new global interior ocean mapped climatology: the 1° × 1° GLODAP version 2, Earth Syst. Sci. Data, 8, 325–340, https://doi.org/10.5194/essd-8-325-2016, 2016.

Lee, K., Kim, T. W., Byrne, R. H., Millero, F. J., Feely, R. A., and Liu, Y. M.: The universal ratio of boron to chlorinity for the North Pacific and North Atlantic oceans, Geochim. Cosmochim. Ac., 74, 1801–1811, 2010.

Le Quéré, C., Andrew, R. M., Friedlingstein, P., Sitch, S., Pongratz, J., Manning, A. C., Korsbakken, J. I., Peters, G. P., Canadell, J. G., Jackson, R. B., Boden, T. A., Tans, P. P., Andrews, O. D., Arora, V. K., Bakker, D. C. E., Barbero, L., Becker, M.,
Betts, R. A., Bopp, L., Chevallier, F., Chini, L. P., Ciais, P., Cosca, C. E., Cross, J., Currie, K., Gasser, T., Harris, I., Hauck, J., Haverd, V., Houghton, R. A., Hunt, C. W., Hurtt, G., Ilyina, T., Jain, A. K., Kato, E., Kautz, M., Keeling, R. F., Klein Goldewijk, K., Körtzinger, A., Landschützer, P., Lefèvre, N., Lenton, A., Lienert, S., Lima, I., Lombardozzi, D., Metzl, N., Millero, F., Monteiro, P. M. S., Munro, D. R., Nabel, J. E. M. S., Nakaoka, S., Nojiri, Y., Padin, X. A., Peregon, A., Pfeil, B., Pierrot, D., Poulter, B., Rehder, G., Reimer, J., Rödenbeck, C., Schwinger, J., Séférian, R., Skjelvan, I., Stocker, B. D., Tian, H., Tilbrook, B., Tubiello, F. N., van der Laan-Luijkx, I. T., van der Werf, G. R., van Heuven, S., Viovy, N., Vuichard, N., Walker, A. P., Watson, A. J., Wiltshire, A. J., Zaehle, S., and Zhu, D.: Global Carbon Budget 2017, Earth Syst. Sci. Data, 10, 405–448, https://doi.org/10.5194/essd-10-405-2018, 2018.

Lewis, E. and Wallace, D. W. R.: Program Developed for CO₂ Sys-
tem Calculations, ORNL/CDIAC-105, Carbon Dioxide Informa-
tion Analysis Center, Oak Ridge National Laboratory, U.S. De-
partment of Energy, Oak Ridge, Tennessee, 38 pp., 1998.

Lueker, T. J., Dickson, A. G., and Keeling, C. D.: Ocean pCO₂ calculated from dissolved inorganic carbon, alkalinity, and equa- tions for K₁ and K₂: validation based on laboratory measure- ments of CO₂ in gas and seawater at equilibrium, Mar. Chem., 70, 105–119, 2000.

Mehrbach, C., Culberson, C. H., Hawley, J. E., and Pytkowicz, R. M.: Measurement of the apparent dissociation constants of carbonic acid in seawater at atmospheric pressure, Limnol. Oceanogr., 18, 897–907, 1973.

Millero, F. J.: Thermodynamics of the carbon dioxide system in the oceans, Geochim. Cosmochim. Ac., 59, 661–677, 1995.

Millero, F. J.: The Marine Inorganic Carbon Cycle, Chem. Rev., 107, 308–341, 2007.

Millero, F. J.: Carbonate constants for estuarine waters, Mar. Fresh-
water Res., 61, 139–142, 2010.

Millero, F. J., Pierrot, D., Lee, K., Wanninkhof, R., Feely, R., Sabine, C. L., Key, R. M., and Takahashi, T.: Dissociation constants for carbonic acid determined from field measurements, Deep-Sea Res. Pt. I, 49, 1705–1723, 2002.

Millero, F. J., Graham, T. B., Huang, F., Bustos-Serrano, H., and Pierrot, D.: Dissociation constants of carbonic acid in seawater as a function of salinity and temperature, Mar. Chem., 100, 80–94, https://doi.org/10.1016/j.marchem.2005.12.001, 2006.

Mucci, A.: The solubility of calcite and aragonite in seawater at various salinities, temperatures and one atmosphere total pressure, Am. J. Sci., 283, 780–799, https://doi.org/10.2475/ajs.283.7.780, 1983.

Naviaux, J. D., Subhas, A. V., Dong, S., Rollins, N. E., Liu, X., Byrne, R. H., Berelson, W. M., and Adkins, J. F.: Calcite dissolution rates in seawater: lab vs. in situ measurements and inhibition by organic matter, Mar. Chem., 215, 103684, https://doi.org/10.1016/j.marchem.2019.103684, 2019.

Olsen, A., Key, R. M., van Heuven, S., Lauvset, S. K., Velo, A., Lin, X., Schirnick, C., Kozyr, A., Tanhua, T., Hoppema, M., Jutterström, S., Steinfeldt, R., Jeansson, E., Ishii, M., Pérez, F. F., and Suzuki, T.: The Global Ocean Data Analysis Project version 2 (GLODAPv2) – an internally consistent data product for the world ocean, Earth Syst. Sci. Data, 8, 297–323, https://doi.org/10.5194/essd-8-297-2016, 2016.

Olsen, A., Lange, N., Key, R. M., Tanhua, T., Álvarez, M., Becker,
S., Bittig, H. C., Carter, B. R., Cotrim da Cunha, L., Feely,
R. A., van Heuven, S., Hoppema, M., Ishii, M., Jeansson, E.,
Jones, S. D., Jutterström, S., Karlsen, M. K., Kozyr, A., Lau-
vset, S. K., Lo Monaco, C., Murata, A., Pérez, F. F., Pfeil,
B., Schirnick, C., Steinfeldt, R., Suzuki, T., Telszewski, M.,
Tilbrook, B., Velo, A., and Wanninkhof, R.: GLODAPv2.2019 –
an update of GLODAPv2, Earth Syst. Sci. Data, 11, 1437–1461,
https://doi.org/10.5194/essd-11-1437-2019, 2019.

Orr, J. C., Epitalon, J.-M., Dickson, A. G., and Gatuso, J.-P.: Routine uncertainty propagation for the marine carbon dioxide system, Mar. Chem., 207, 84–107, https://doi.org/10.1016/j.marchem.2018.10.006, 2018.

Papadimitriou, S., Loucaides, S., Rérolle, V. M. C., Kennedy, P., Achterberg, E. P., Dickson, A. G., Mowlem, M., and Kennedy, H.: The stoichiometric dissociation constants of carbonic acid in seawater brines from 298 to 267 K, Geochim. Cosmochim. Ac., 220, 55–70, https://doi.org/10.1016/j.gca.2017.09.037, 2018.

Patsavas, M. C., Byrne, R. H., Wanninkhof, R., Feely, R. A., and Cai, W.-J.: Internal consistency of marine carbonate system measurements and assessments of aragonite saturation state: Insights from two U.S. coastal cruises, Mar. Chem., 176, 9–20, 2015.

Perez, F. F. and Fraga, F.: Association constant of fluoride and hy-
drogen ions in seawater, Mar. Chem., 21, 161–168, 1987.

Pierrot, D., Lewis, E., and Wallace, D. W. R.: MS Excel Pro-
gram Developed for CO₂ System Calculations, Tech. rep.,
ORNL/CDIAC-105a, Carbon Dioxide Information Analysis
Center, Oak Ridge National Laboratory, U.S. Department of En-
ergy, Oak Ridge, Tennessee, 2006.

Pierrot, D., Neill, C., Sullivan, K., Castle, R., Wanninkhof, R., Lüger, H., Johannessen, T., Olsen, A., Feely, R. A., and Cosca, C. E.: Recommendations for autonomous underway pCO₂ measuring systems and data-reduction routines, Deep-Sea Res. Pt. II, 56, 512–522, https://doi.org/10.1016/j.dsr2.2008.12.005, 2009.

R Core Team: R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria, available at: http://www.R-project.org/ (last access: 5 March 2020), 2019.

Raimondi, L., Matthews, J. B. R., Atamanchuck, D., Azetsu-Scott, K., and Wallace, D.: The internal consistency of the marine carbon dioxide system for high latitude shipboard and in situ monitoring, Mar. Chem., 213, 49–70, 2019.

Ribas-Ribas, M., Rérolle, V. M. C., Bakker, D. C. E., Kitidis, V., Lee, G. A., Brown, I., Achterberg, E. P., Hardman-Mountford, N. J., and Tyrrell, T.: Intercomparison of carbonate chemistry measurements on a cruise in northwestern European shelf seas, Biogeosciences, 11, 4339–4355, https://doi.org/10.5194/bg-11-4339-2014, 2014.

Roy, R. N., Roy, L. N., Vogel, K. M., Porter-Moore, C., Pearson, T., Good, C. E., Millero, F. J., and Campbell, D. M.: The dissociation constants of carbonic acid in seawater at salinities 5 to 45 and temperatures 0 to 45 °C, Mar. Chem., 44, 249–267, 1993.

Sharp, J. D. and Byrne, R. H.: Carbonate ion concentrations in seawater: Spectrophotometric determination at ambient temperatures and evaluation of propagated calculation uncertainties, Mar. Chem., 209, 70–80, 2019.

Spiess, A.-N.: propagate: Propagation of uncertainty using higher-order Taylor expansion and Monte Carlo simulation. R-
package version 1.0-6, available at: https://CRAN.R-project.org/
package=propagate (last access: 5 March 2020), 2018.

Soetaert, K. and Herman, P. M. J.: A Practical Guide to Ecological Modelling. Using R as a Simulation Platform, Springer, 372 pp., 2009.

Sulpis, O., Lauvset, S. K., and Hagens, M.: R script belonging to “Current estimates of K1* and K2* appear inconsistent with measured CO₂ system parameters in cold oceanic regions”, Zenodo, https://doi.org/10.5281/zenodo.3725889, 2020.

Takahashi, T., Sweeney, C., Hales, B., Chipman, D., Newberger, T., Goddard, J., Iannuzzi R., and Sutherland, S.: The Changing Carbon Cycle in the Southern Ocean, Oceanography, 25, 26–37, 2012.

Takeshita, Y., Johnson, K. S., Martz, T. R., Plant, J. N., and Sarmiento, J. L.: Assessment of Autonomous pH Measurements for Determining Surface Seawater Partial Pressure of CO₂, J. Geophys. Res.-Oceans, 123, 4003–4013, 2018.

Tanhua, T., van Heuven, S., Key, R. M., Velo, A., Olsen, A., and Schirnick, C.: Quality control procedures and methods of the CARINA database, Earth Syst. Sci. Data, 2, 35–49, https://doi.org/10.5194/essd-2-35-2010, 2010.

Tishchenko, P. Y., Wong, C. S., and Johnson, W. K.: Measurements of Dissociation Constants of Carbonic Acid in Synthetic Seawater by Means of a Cell Without Liquid Junction, J. Solution Chem., 42, 2168–2186, https://doi.org/10.1007/s10953-013-0094-7, 2013.

Uppström, L. R.: The boron/chlorinity ratio of deep-sea water from the Pacific Ocean, Deep-Sea Research and Oceanographic Abstracts, 21, 161–162, 1974.

van Heuven, S., Pierrot, D., Rae, J. W. B., Lewis, E., and Wallace, D. W. R.: MATLAB program developed for CO₂ system calculations, ORNL/CDIAC-105b, Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, U.S. Department of Energy, Oak Ridge, Tennessee, 2011.

Weiss, R.: Carbon dioxide in water and seawater: the solubility of a non-ideal gas, Mar. Chem., 2, 203–215, 1974.

Williams, N. L., Juranek, L. W., Feely, R. A., Johnson, K. S., Sarmiento, J. L., Talley, L. D., Dickson, A. G., Gray, A. R., Wanninkhof, R., Russell, J. L., Riser, S. C., and Takeshita, Y.: Calculating surface ocean pCO₂ from biogeochemical Argo floats equipped with pH: An uncertainty analysis, Global Biogeochem. Cy., 31, 591-604, 2017.

Yang, B., Byrne, R. H., and Lindemuth, M.: Contributions of organic alkalinity to total alkalinity in coastal waters: A spectrophotometric approach, Mar. Chem., 176, 199–207, 2015.

Yao, W. and Millero, F. J.: The chemistry of the anoxic waters in the Framvaren Fjord, Norway, Aquat. Geochem., 1, 53–88, 1995.

Yao, W., Liu, X., and Byrne, R. H.: Impurities in indicators used for spectrophotometric seawater pH measurements: Assessment and remedies, Mar. Chem., 107, 167–172, 2007.

Zeebe, R. E. and Wolf-Gladrow, D.: CO₂ in seawater: equilibrium, kinetics, isotopes, Elsevier Science, Amsterdam, 2001.