Interannual sea-air CO2 flux variability from an observation-driven ocean mixed-layer scheme

Interannual anomalies in the sea-air carbon diox- ide (CO2) exchange have been estimated from surface-ocean CO2 partial pressure measurements. Available data are suffi- cient to constrain these anomalies in large parts of the trop- ical and North Pacific and in the North Atlantic, in some ar- eas covering the period from the mid 1980s to 2011. Global interannual variability is estimated as about 0.31 Pg C yr 1 (temporal standard deviation 1993-2008). The tropical Pa- cific accounts for a large fraction of this global variabil- ity, closely tied to El Nino-Southern Oscillation (ENSO). Anomalies occur more than 6 months later in the east than in the west. The estimated amplitude and ENSO response are roughly consistent with independent information from atmo- spheric oxygen data. This both supports the variability esti- mated from surface-ocean carbon data and demonstrates the potential of the atmospheric oxygen signal to constrain ocean biogeochemical processes. The ocean variability estimated from surface-ocean carbon data can be used to improve land CO2 flux estimates from atmospheric inversions.


Introduction
The ocean currently accounts for about half the sink of excess atmospheric CO 2 (Sarmiento et al., 2010). Long-term changes in this sink capacity therefore affect the climate change trajectory. In order to be able to relate changes in sea-air CO 2 exchanges to driving influences and to test concepts or models, contemporary variations need to be quantified from suitable data sets. pressure (pCO 2 ) and sea-air CO 2 flux presented in the companion paper Rödenbeck et al. (2013) (run SFC). As documented in detail there, these fields were estimated by fitting a simple diagnostic model of mixed-layer biogeochemistry to SOCAT pCO 2 data in the following way: (1) A mixed-layer carbon budget equation (including simple parameterizations of sea-air CO 2 exchange (Wanninkhof, 1992), solubility (Weiss, 10 1974), and carbonate chemistry (Sarmiento and Gruber, 2006), driven by observationbased environmental fields listed in Table 1) was used to express surface-ocean pCO 2 and sea-air CO 2 flux as a function of ocean-internal carbon sources and sinks ( Fig. 1) for each pixel. (2) A cost function was formed to measure the mismatch between the individual pCO 2 data points in SOCAT and the model's pCO 2 field at the corresponding 15 location and time. (3) The ocean-internal carbon sources and sinks were adjusted as to minimize this model-data mismatch. In order to interpolate areas/periods without data, additional smoothness constraints were applied (similar to those in the atmospheric inversion of Rödenbeck, 2005). The only difference of the present estimates with respect to Rödenbeck et al. (2013) is that we updated the data source to the ship-board ob-Introduction The calculation has been done on a global grid of ≈ 4 • × 5 • pixels and daily time steps over 1985-2012 (inclusive). The first and last year will be discarded from any plots to avoid edge effects. Statistical analysis is done over 1993-2008 (inclusive) only, as this is a period with largest data coverage both in SOCAT v2 and in atmospheric oxygen records. To obtain interannual variations, we de-seasonalize the time series by 5 subtracting the mean seasonal cycle (1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008) and remove variations faster than about 1 year by a spectral filter corresponding to a Gaussian smoothing kernel (Rödenbeck, 2005).

Process contributions to variability
For any given location, the (interannual) variability of the sea-air CO 2 flux is the sum 10 of two contributions: Contribution OIS (Ocean-Interior Sources/sinks). Variations in biological conversion and/or ocean transport into the mixed layer lead to anomalous ocean-internal sources/sinks f DIC int of dissolved inorganic carbon (DIC) (Fig. 1 bottom), which then 15 lead to enhanced/reduced sea-air CO 2 fluxes f CO 2 ma ( Fig. 1 top). The responses in f CO 2 ma are however of smaller magnitude than and delayed with respect to the causes in f DIC int , because the buffering capacity of carbonate chemistry (Revelle factor) strongly enhances the limiting effect of the finite gas exchange velocity. 20 Contribution TE (Thermally induced Exchange). Even in the absence of ocean-internal sources or sinks, the sea-air flux varies in response to temperature-induced changes in solubility and chemical equilibrium 2 (plus minor contributions from changes in atmo-Introduction Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | spheric pressure, freshwater effects, and alkalinity). Contribution TE also comprises a secular CO 2 uptake (and consequent rise of pCO 2 ) induced by the prescribed rising atmospheric CO 2 content. Both contributions are somewhat modulated by variations in gas exchange (Supplement, Sect. S5) and mixed-layer depth. 5 In the diagnostic scheme, contribution TE is considered known from the employed process parametrizations driven by observed variables (predominantly Sea Surface Temperature (SST)). In contrast, the ocean-internal DIC sources/sinks f DIC int causing contribution OIS are taken as unknowns, to be adjusted in such a way that the total variations in the CO 2 partial pressure (sum of contributions OIS and TE) are compatible 10 with the SOCAT data (for illustration see also Fig. 3 of Rödenbeck et al., 2013). The Bayesian prior is chosen to be f DIC int = 0, such that any positive or negative anomalies in the ocean-internal sources/sinks are a-priori equally likely, i.e., we do not assume any prior knowledge on the ocean-internal processes. Consequently, the prior of the sea-air flux coincides with contribution TE. 15

Reduction of uncertainty (RoU)
As the main diagnostic to identify which regions are well constrained by the data, we use the "reduction of uncertainty" (RoU), a standard diagnostic of Bayesian estimation defined as which counteract the direct temperature effect. As a consequence of this balance, peaks in the TE contribution tend to be related to temperature changes and thus occur earlier than the actual temperature peaks. with σ pri and σ post denoting the formal a-priori and a-posteriori uncertainties of the flux. High values (RoU towards 1) indicate strong data constraints, while low values (close to 0) indicate that the data are not able to move the estimates away from the prior. The uncertainty intervals σ pri and σ post have been calculated statistically from an ensemble of n = 50 inversion runs with pseudo-random realizations of a-priori errors and 5 model-data mismatch errors (drawn according to their specified uncertainties). Each member i of the pseudo-random ensemble gives a realization ∆ i f post of the a-posteriori error of the flux field. The uncertainties σ pri and σ post are then given by the root mean square of the error fields across the ensemble, σ = i (∆ i f ) 2 /n (as the expectation value of the error is zero by construction). As we want to obtain RoU specifically for 10 interannual variations (and possibly for regional averages), the a-priori and a-posteriori error fields of each ensemble member are interannually filtered (and regionally averaged) in the same way as the actual flux estimates, before the root mean square is calculated. That way, any covariances between the a-posteriori errors at different locations/times are automatically taken into account. 15 Of course, the root mean square across the finite ensemble only gives an approximation to the uncertainty intervals. For our sample size of 50, standard deviations could be a factor 0.79 smaller or a factor 1.34 larger (confidence interval for a confidence level of 99 %). Calculated RoU values can therefore only give the 1st order pattern of the strength of the data constraint. When calculating the overall performance over a given 20 time period, however, we increased the sample size (and thus improved the statistics) by calculating the root mean square error not only across the ensemble but also across all monthly values of the flux error realizations within that period (method proposed by Chevallier et al., 2007). 25 Most driving variables of the mixed-layer scheme (in particular SST) also contain interannual variations. If a region is well-constrained with respect to IAV, however, then the Introduction IAV of the estimated pCO 2 field should only depend on the signals in the pCO 2 data, "overwriting" any influence of IAV in the drivers. We thus performed a test run in which any IAV has been removed from the driver variables (except for a linear rising trend in atmospheric CO 2 ). In well-constrained regions the result will be similar to the standard result. We use this as a qualitative confirmation to the RoU diagnostic.

5
Remark: for quantities other than the pCO 2 field, IAV in the drivers can of course be important also in well-constrained regions, for example IAV in wind speed influencing sea-air CO 2 fluxes.

Synthetic data test
Another traditional method to assess the strength of the data constraint (also used in 10 Rödenbeck et al., 2013, for seasonality) is "synthetic data tests" where (1) a modelled pCO 2 field is sampled at the locations/times of the SOCAT data to create a synthetic data set, and (2) a pCO 2 field is retrieved by the diagnostic mixed-layer scheme from the synthetic data and compared to the original modelled field as the known "true" answer. However, this method is dependent on the particular time course of the cho-15 sen synthetic "truth", confounding any information about the temporal changes in the strength of the constraint. Despite the higher computational demand, we therefore prefer the RoU diagnostic (Sect. 2.3.1) here. Nevertheless, synthetic data tests have been performed to check consistency with the RoU diagnostic (not shown). 20 In addition to gaps in data coverage, results may be affected by uncertainties in parameters or input data sets used in the calculation. To assess these errors, we performed a series of sensitivity runs varying those parameters that are considered most uncertain: (i) increase and decrease of the a-priori uncertainty by a factor 2, leaving more/less freedom to inverse adjustments, (ii) decrease of the a-priori uncertainty of Introduction Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | lengths in latitude direction by a factor 2, (iv) increase and decrease of the global mean piston velocity by 3.2 cm h −1 (range given by Naegler, 2009), or using a cubic dependence of piston velocity on wind speed (still keeping the global mean piston velocity at the value from Naegler, 2009), (v) increase and decrease of the mixed-layer depth h by a factor 2. As all these changes are still considered reasonable values, the envelope of 5 these sensitivity results (range between smallest and largest value at each time) gives a lower bound of the error (to be plotted as grey bands).

Residuals and comparison to independent time series
As a consistency check of the fit, residuals between the estimated pCO 2 field and the SOCAT data are assessed. According to the Supplement (Sect. S3.1), the data 10 are fit within 10 to 15 µatm (root-mean-squared residuals), consistent with assumed mismatch uncertainty. Any remaining interannual signals in the residuals are negligible, confirming that the data information is used completely. Short-term variations are further validated against high-frequency time series from various moorings not used in the fit (Supplement, Sect. S3.2). For testing interannual 15 variations, independent continuous long-term surface-ocean carbon time series are used (Supplement, Sect. S3.3).

Lagged correlation analysis
To assess the relation of estimated interannual flux anomalies to ENSO, we correlate them to the Multivariate El Niño Index (MEI) by Wolter and Timlin (1993). MEI has 20 been interannually filtered in the same way as the flux (Sect. 2.1). The correlation is calculated for different time lags between flux and MEI; we always report the maximum correlation coefficient. If local maxima exist for several lags, we chose the smallest absolute lag. Correlations have been calculated over the 1993-2008 analysis period. Assuming each of these 16 years to be statistically independent, correlations are moderately significant (at the 90 % level) if their correlation coefficient exceeds 0.50.

Calculating sea-air oxygen fluxes
In order to relate the SOCAT-based variability estimates to inverse estimates by Rö-5 denbeck et al. (2008) based on atmospheric oxygen data, we extend the diagnostic mixed-layer scheme in three steps:

Ocean-internal oxygen sources/sinks
The ocean-internal processes (biology, transport) that add or remove dissolved inorganic carbon to/from the mixed-layer also add or remove dissolved oxygen (symbolized 10 as box "TT" in Fig. 2). Biological respiration/photosynthesis lead to opposite, mutually proportional changes in carbon and oxygen that we assume to follow a fixed Redfield stoichiometry. Mixing of water from the ocean interior into the mixed layer can change its carbon and oxygen concentrations as well. These transport-induced changes do not follow a simple stoichiometry, but analysis of vertical tracer gradients suggests that 15 their mutual relation is also roughly in Redfield proportions. We therefore assume with the Redfield ratio r O 2 :C = −150/106 ≈ −1.4 (Anderson, 1995) as an approximation.

Sea-air oxygen fluxes
Sea-air oxygen fluxes f the mixed-layer oxygen budget and sea-air oxygen exchange (Fig. 2, middle column). These oxygen parameterizations are analogous to those of carbon (left column) detailed in Rödenbeck et al. (2013) Rödenbeck et al., 2013), (ii) there is no equivalent to carbonate chemistry involved, (iii) sea-air oxygen exchange contains an additional contribution from air injection through bubbles 3 .

Atmospheric Potential Oxygen (APO) fluxes
The use of atmospheric oxygen data involves the complication that the atmospheric 5 oxygen abundance is also influenced by oxygen exchanges from the land biosphere. To account for that, the inversion (Rödenbeck et al., 2008) has not been applied to oxygen 3 Gas injection by bubbles is parametrized by with the oxygen mixing ratio X O 2 = 209 392 ppm. This assumes complete dissolution of bubbles (shown to be the dominating process, Stanley et al., 2009). The air injection flux (positive for injection into the ocean) is taken to depend on wind speed u according to (Monahan and Torgersen, 1990) f Air bubb = with the threshold wind speed u 0 = 2.27 m s −1 , barometric pressure p baro , gas constant R gas = 8.3144 J mol −1 K −1 , and absolute temperature T +T 0 . The global constant B has been calculated such that the global diffusive gas exchange f O 2 ge calculated from the World Ocean Atlas (WOA) climatology of surface-ocean dissolved oxygen (Garcia et al., 2006) equals the global bubble injection flux, such that the long-term global sea-air oxygen flux is zero.
It should be noted that the bubble flux has only very little influence on the sea-air oxygen fluxes fluxes calculated here from internal O 2 sources and sinks: due to the fast equilibration of dissolved oxygen with the atmosphere (within less than a month), the oxygen input by bubbles is essentially compensated (on interannual time scales) by a corresponding enhancement in the diffusive sea-air oxygen exchange. Introduction

Tables Figures
Back Close

Full Screen / Esc
Printer-friendly Version

Interactive Discussion
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | itself but to the conceptual tracer "Atmospheric Potential Oxygen" (APO) introduced by Stephens et al. (1998). APO is defined as a combination of O 2 and CO 2 abundances in such a way that its surface-to-atmosphere fluxes are As O 2 exchanges from the land biosphere are, to good approximation, −1.1 times 5 the CO 2 exchanges, they cancel each other in the APO flux. There is a remaining contribution from fossil fuel burning as it has a slightly different stoichiometry (about −1.4, Keeling, 1988), but this contribution has been accounted for in the APO inversion based on fuel-use statistics (see Rödenbeck et al., 2008). Thus, the APO inversion yields an estimate of the oceanic APO flux f APO ma .

10
The comparison is therefore done in terms of f APO ma . The SOCAT-based APO flux has been calculated according to Eq. (5) from the (dominant) oxygen contribution f O 2 ma obtained as described above, and the (small) carbon contribution f CO 2 ma directly available from the SFC calculation. The atmospheric inversion of Rödenbeck et al. (2008) has been updated by adding more recent observations to extend the time period to end of 15

Overview
While the sea-air CO 2 flux varies on interannual, seasonal, and day-to-day time scales (Fig. 3, top), this paper focuses on interannual anomalies of the sea-air CO 2 flux 20 around its mean (filtered flux, bottom). The largest contributions to the global interannual variability (IAV) are estimated to come from the subpolar North Pacific and Atlantic, the tropical Pacific, and parts of the Southern Ocean (Fig. 4, top). Due to its large size and its spatially coherent variations, the tropical Pacific ( , Feely et al., 1999;Inoue et al., 2001;Feely et al., 2002), the sea-air CO 2 flux anomalies are strongly tied to El Niño/Southern Oscillation (ENSO), with a reduced CO 2 source during El Niño phases (see Sect. 3.4 below).

5
To identify which parts of the ocean have been sampled frequently enough to reliably estimate interannual flux variations from the SOCAT data, Fig. 4 (bottom) shows the "reduction of uncertainty" (RoU) defined in Sect. 2.3.1. It indicates the strength of the data to detect deviations of the flux from its climatological mean. Good constraints (larger RoU) are found in the tropical and Northern Pacific and in the North Atlantic, 10 highlighting the large long-running observational programs by US, Japanese, and European research groups. Specifically, the locations of the repeated cruises crossing the Equator are seen, as well as the ship routes traversing the North Atlantic. In contrast, poor constraints (low RoU, even below the arbitrary threshold of 0.2 left uncolored in the map) prevail in the Indian Ocean and most parts of the Southern-hemisphere ex- 15 tratropical Ocean. Figure 4 (bottom) gives the overall performance within 1993-2008. Depending on the region considered, this period may comprise both well constrained and badly constrained years (Fig. S1). For example, the RoU in the extratropical North Atlantic clearly increased after the implementation of the pCO 2 observing network in the first half of 20 the 2000s (Watson et al., 2009) and has remained steady at between approximately 0.6 and 0.8 since then.
In well-constrained regions/periods (high RoU), estimated interannual anomalies are essentially independent of the prior. For example, the data are able to almost reverse the a-priori anomalies in the tropical Pacific (  Pacific are further confirmed by a test run without interannual variations in all driving variables (Sect. 2.3.2): despite the corresponding substantial change in the prior, the estimated pCO 2 anomalies hardly change (not shown). In contrast, poorly constrained regions (low RoU) stay close to the prior (Fig. S1). In areas where the temperaturerelated variability (contribution TE contained in the prior, Sect. 2.2) dominates IAV, the 5 results of the diagnostic mixed-layer scheme may still capture part of the real IAV, otherwise the variability from the ocean-internal fluxes (contribution OIS) is missing there.
The pattern of RoU shown in Fig. 4 is similar to the map of number of data points ( Fig. S7.4 of Rödenbeck et al., 2013), but additionally takes into account that any data 10 point in an already well-covered area has less individual effect. Though RoU is to some extent also influenced by our a-priori uncertainty settings (e.g., a smaller a-priori uncertainty would lower the achievable RoU, and the applied smoothness constraints carry information into unconstrained areas close to constrained ones), the patterns of Fig. 4 (bottom) should be broadly representative for SOCAT's information content on 15 interannual anomalies, largely independently of the method used here.

Robustness
To exclude that the results are dominated by uncertain parameters or input data sets used in the calculation, Fig. 5 shows the envelope of a set of sensitivity results (Sect. 2.3.4), using the tropical Pacific as example region. Reflecting the fact that pCO 2 20 (middle panel) is the quantity directly constrained by the data, alterations in parameters only have a very small effect (sensitivity band much more narrow than the interannual variations). There are short periods of enhanced sensitivity in 2009 and before about 1992, when the tropical Pacific is less well constrained by the SOCAT data such that it stays closer to the prior and thus the effects of our sensitivity cases on the prior cannot 25 be completely overwritten.
The sensitivity of the sea-air CO 2 fluxes (Fig. 5 despite some change in amplitude, the time course of the interannual flux variations is very robust. Larger sensitivity is found in the ocean-internal DIC sources/sinks (Fig. 5,bottom). This is expected as any alteration in the parameterized relation between f DIC int and pCO 2 enforces compensating changes in f DIC int . The biggest effect comes from the alterations 5 in mixed-layer depth (MLD). Note however that our alterations of MLD by factors 2 and 0.5 are ad-hoc; this large range has been chosen to account for the missing interannual variations in MLD values used, but may overestimate the actual MLD errors.
In regions less well constrained than the tropical Pacific, estimates are less robust. Besides the above-mentioned effect that the alterations of the prior have more impact 10 on the result, estimates in unconstrained regions critically depend on the choice of spatial correlation length: the longer the correlations, the farther variability from dataconstrained locations is spread into unconstrained locations (see Fig. 5 of Rödenbeck et al., 2013). In particular, longer correlation carry the strong variability of the tropical Pacific farther into the poorly constrained areas south of it (Fig. 4, bottom), thus in- 15 creasing the areal extent of this variability and leading to higher global IAV 5 . This leads to a sensitivity band for the global sea-air CO 2 fluxes (Fig. 3, bottom) that is somewhat wider than in the tropical Pacific. As the band is asymmetric, it suggests a somewhat lower amplitude of global IAV than in the standard case.

Interannual variability in the tropical Pacific
As the presented material establishes that the tropical Pacific is both a well-constrained and a globally important contributor to interannual sea-air CO 2 flux variability 6 , we consider the contributions to the variability at different longitudes along the Equator (Fig. 6).
During the 1993-2008 analysis period, fluxes at most longitudes are continuously constrained (left panel); the Western part even throughout the earlier years. The reduction of the flux during El Niño phases (middle, red colors) extends throughout considered part of the equatorial Pacific (200 • W to 85 • W), but with the tendency to be of smaller amplitude and to occur several months later in the East than in the West. This striking "propagation" can be confirmed statistically: the flux is significantly anti-correlated  Fig. S4 for illustration). Such a slow West-to-East "propagation" of sea-air CO 2 flux anomalies as con- 15 strained by SOCAT is not present in physical surface ocean variables: for example, SST is positively correlated to MEI, with almost no lag and propagating much faster (within about 1 month consistent with an equatorial Kelvin wave, not shown). 7 No slow "propagation" is thus found in the SST-dominated TE contribution either: consistent with Footnote 2 on page 6, the TE contribution is positively correlated to MEI as well (mid-20 6 Ocean regions other than the tropical Pacific either have much less IAV, or are not well constrained by the available pCO 2 observations (Fig. S1). Comparison with independent data at station HOT in the North Pacific broadly confirms the SOCAT-based anomalies, including the low IAV there (Sect. S3.3). Comparison at BATS is not successful in terms of interannual features, but still confirms the low amplitude in the North Atlantic as well (Sect. S3.3). Despite their lower IAV, various extratropical regions can be expected to play an important role in the long-term flux trend. A more detailed consideration warrants a separate study. 7 In weak ENSO events, SST anomalies can even propagate westward, see Santoso et al. (2013). shows a West-to-East "propagation" of only 1-2 months. This suggests that the slow "propagation" of the total anomalies is only an apparent one and does not arise from any actual propagation mechanism. Rather, as the amplitudes of the opposite TE and 5 OIS contributions are markedly changing from West to East (bottom right panel), the temporal phase of their superposition is shifting (in the same way as sine and cosine Fourier terms of different relative weight add up to modes with different phase). Can the remaining "propagation" in the OIS contribution arise from memory effects in the mixed-layer carbon budget? Anomalies in contribution OIS occur several months 10 later than the causing anomalies in the ocean-interior DIC sources/sinks (not shown, but see the Fig. S4), due to the finite sea-air exchange rate enhanced by the buffer effect of carbonate chemistry. Anomalies in f DIC int thus even lead MEI by 5-6 months, i.e., are associated with build-up or decline of the ENSO states. However, this delay between f DIC int and f CO 2 ma is essentially the same throughout the equatorial Pacific, i.e., it 15 does not contribute to the propagation.

Consistency to atmospheric oxygen observations
In order to relate the SOCAT-based IAV estimates to the information from atmospheric oxygen, the oceanic APO flux has been calculated from run SFC (Sect. 2.5). The information flow is illustrated in Fig. 7. Sea-air CO 2 flux (f The APO flux inferred from SOCAT is compared to that inferred from atmospheric oxygen data by a transport inversion (Rödenbeck et al., 2008, updated) (Fig. 7 top). Both of these independent estimates show enhanced APO outgassing during El Niño phases (grey stripes). While the SOCAT-based anomalies consistently lead ENSO events (due to the lead of f DIC int , compare Sect. 3.4), the APO-based anomalies co-5 incide with the 1995 and 1997/98 ENSO events but also lead the 2009 event. Both estimates are of roughly similar amplitudes, though the relative heights of the anomalies in individual ENSO events are different.
Various conceivable error influences have the potential to distort the comparison: (1) Both the SOCAT and the APO constraints suffer from incomplete spatial coverage, 10 where the regions of good coverage do not necessarily coincide. An APO inversion of "synthetic data" suggests that this may explain part of the differences in timing (Supplement, Sect. S4).
(2) The information flow from the SOCAT data ( Fig. 2) involves various uncertain model steps: even if run SFC would lead to a perfectly constrained p (3) The APO inversion involves various uncertainties both from necessary assumptions in the estimation and from observational error (Rödenbeck et al., 2008). Also, interannual variations in the sea-air N 2 (influencing the APO inversion because the atmospheric O 2 abundance is measured relative to the N 2 abundance) are neglected. 25 In the light of all these error influences, but also considering that the SOCAT and APO constraints are fully independent from each other in terms of both data and model, the partial agreement in the interannual APO flux variations is remarkable. Even though BGD Introduction

Tables Figures
Back Close

Full Screen / Esc
Printer-friendly Version

Interactive Discussion
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | the comparison is not conclusive enough to provide a quantitative validation, we take it as a confirmation of plausibility of the mixed-layer scheme. Conversely, the comparison also confirms that atmospheric oxygen observations contain information about ocean-internal carbon sources/sinks and thus sea-air CO 2 exchange. To use this information as an additional constraint, we need to tackle the er-5 ror influences listed above, in particular (1) increase the information from atmospheric oxygen by adding existing or new observation sites, (2) use a more realistic representation of carbon-oxygen coupling (TT) which may involve further data sources. In any case, a crucial role is played by the availability of atmospheric oxygen records that cover the entire time period under consideration without interruption.

Comparison to simulations by an ocean process model -IAV and trend
The data-based estimates presented here can be used to challenge the results of process simulations for well-constrained features, but conversely, comprehensive process models can be used to challenge our results where constraints are weak such that we rely on our simple parameterizations. Therefore, sea-air CO 2 fluxes have been 15 calculated from pCO 2 simulations by an ocean biogeochemical process model run (NEMOv2.3 with PlankTOM5, Buitenhuis et al., 2010) (but using the same gas exchange parametrizations as in our diagnostic scheme). As in the SOCAT-based run SFC, the total ocean flux variability (Fig. 8) is strongly tied to ENSO. Though the process model gives a smaller amplitude than SFC (about half), this model is the one 20 with the largest interannual variability among the suite of process models considered in Wanninkhof et al. (2013) (their Fig. 6), and agrees by far the best with our SOCATbased estimate.
In both the model simulation and our data-based run SFC, total ocean flux variability is dominated by the tropical Pacific, where both estimates are in closest agreement 25 (Supplement, Fig. S2). The estimates further agree that all other regions have smaller IAV, though there is hardly any correspondence in the detailed features. The regional comparison reveals that the differences in global amplitude are to large extent related 3186 Introduction to areas south of the tropical Pacific, reinforcing the possibility that our amplitude is overestimated due to the tropical variability being spread over too great an area (end of Sect. 3.3). Rising atmospheric CO 2 content leads to rising CO 2 -undersaturation of the global ocean and thus an increasing oceanic sink. The magnitude of this trend, however, 5 cannot be predicted from simple considerations as it critically depends on the rate by which the sequestered carbon is passed on from the surface ocean into the interior. Therefore process model simulations are needed. From a suite of models, Le Quéré et al. (2013) (Global Carbon Project, GCP) give a sink time series increasing by −0.032 (Pg C yr −1 ) yr −1 over 1960-2012 (linear fit, converted to atmospheric sign convention). While the period of good data coverage is too short to compare this long-term trend, both our data-based estimates, NEMO-PlankTOM5, and the GCP model mean agree that the 1990-1999 period saw a negligible or even reversed trend, followed by the 2000-2009 period with a trend almost 50 % stronger than long-term (trend lines in Fig. 8, also see Fig. 5c of Le Quéré et al., 2013). Unfortunately, a more quantitative 15 analysis of these short-term trends is difficult due to the large interannual variability and region-to-region differences (e.g., Fay and McKinley, 2013).

Combination with an atmospheric CO 2 inversion -implication for land flux estimates
Interannual variations of regional sea-air CO 2 fluxes can also be estimated from atmo-20 spheric CO 2 mixing ratio measurements by atmospheric transport inversion (Newsam and Enting, 1988;Rayner et al., 1999). However, the atmospheric signal is dominated by the much larger variability of land-atmosphere CO 2 fluxes. Though the total (land+ocean) flux within latitude bands is relatively well-constrained from atmospheric CO 2 data (due to atmospheric tracer mass conservation and the predominantly longitu- estimates and an atmospheric inversion (Fig. 9, top): the ocean total from the atmospheric inversion is almost anti-correlated to run SFC, but rather more in phase with land biosphere variability (see Fig. 9, bottom), suggesting that these variations are spuriously spilling over from there. This is in line with a relatively large spread in the ocean fluxes from different atmospheric inversions (see the RECCAP ensemble, Peylin 5 et al., 2013), confirming a limited constraint of the atmospheric data on land/ocean flux partitioning. Similar discrepancies between SOCAT-based and atmospheric inversion estimates are also found in regional fluxes (Supplement, Fig. S3), though at least the relative amplitudes of flux IAV between the regions broadly agree. Surprisingly similar variations are found in the tropical Indian.

10
As the SOCAT-based estimates show more plausible ocean flux IAV than the atmospheric inversion at least in the region contributing the largest variability (Equatorial Pacific), it would be beneficial to add pCO 2 as further constraint to the atmospheric inversion: due to the atmospheric mass conservation, land flux estimates profit from any improvement in the poorly constrained ocean fluxes. Even though interannual vari- 15 ations outside the Equatorial Pacific are not well constrained from the pCO 2 data in many places and thus stay close to the prior of the diagnostic scheme (Sect. 3.2), they still contain the temperature-related part as represented in the parametrizations. In practical terms, even in areas where this is a bad approximation, using the SOCATbased estimates will not worsen the land fluxes from the atmospheric inversion much 20 (compared to an inversion with a freely adjustable ocean flux) due to their small amplitude. Only in the Southern Ocean where influence from adjacent land regions is smallest, the constraint from atmospheric CO 2 data may be powerful, and able to compensate the weakness of the pCO 2 constraint there.
The inclusion of the pCO 2 constraint into the atmospheric inversion can be imple- 25 mented by adding together the cost function contributions of pCO 2 data (from SFC) and of atmospheric CO 2 (see the companion paper Rödenbeck et al., 2013, Appendix A2.3). However, as the ocean fluxes are much more strongly constrained by the pCO 2 data rather than the atmospheric data, such a joint inversion gives ocean fluxes almost identical to run SFC (test run not shown) 8 . Given this, we may simply use the results of run SFC as a fixed ocean prior in a subsequent "classical" atmospheric inversion in which only the land fluxes are adjusted, which is much more practical 9 .
Figure 9 (bottom) compares total land CO 2 fluxes from atmospheric inversions with adjustable or fixed ocean fluxes. Using the SOCAT-based run SFC as prior slightly increases the land flux IAV. This is largely due to differences in South America (Fig. S3) where the density of atmospheric measurement sites is low, such that the distinction of land and ocean fluxes is not well constrained in the "classical" atmospheric inversion. Differences also occur in the Asian regions. However, the impact is still roughly within the range of many other uncertainties in global atmospheric CO 2 inversions (e.g., see the spread of results by various inversion studies participating in the RECCAP ensemble, Peylin et al., 2013). Future developments of carbon cycle quantification should use the good constraint of pCO 2 data on the tropical Pacific and much of the Northern extratropics, but also the constraint of atmospheric CO 2 data on the Southern Ocean. Further, a joint CO 2 and 15 APO inversion, linking sea-air CO 2 and O 2 fluxes through the diagnostic mixed-layer scheme, would exploit the potential of atmospheric oxygen observations to constrain the ocean-internal processes driving the sea-air CO 2 flux variability as demonstrated in Sect. 3.5. The APO constraint applies on the same large spatial scales that drive the air-sea CO 2 flux variability, and thus would be a valuable complement to the pCO 2 20 and atmospheric CO 2 data in the previously underconstrained areas.

Conclusions
Based on the SOCAT v2 data set of pCO 2 observations, we estimated interannual variations of the sea-air CO 2 exchange, using a data-driven diagnostic scheme of mixed-layer carbon biogeochemistry as a space-time interpolator. The scheme links sea-air CO 2 exchange to ocean-internal DIC sources and sinks, thereby also allowing 5 to relate carbon anomalies to signals in oxygen (or nutrient) observations.
-SOCAT pCO 2 data constrain interannual variations in the sea-air CO 2 flux in parts of the ocean, notably in the tropical Pacific where the largest interannual variations are found.
-The tropical Pacific shows a reduced CO 2 outgassing during El Niño phases. In 10 the East, this anomaly occurs more than 6 months later than in the West, likely due to different relative contributions of temperature-related and biologically/physically caused responses.
-The SOCAT-based estimates of the interannual variations in tropical sea-air CO 2 flux and ocean-internal DIC sources/sinks are broadly consistent with the inde-15 pendent constraint from atmospheric oxygen measurements.
-In qualitative agreement with ocean process model simulations, the ocean sink in the SOCAT-based estimates increases less than expected from atmospheric CO 2 rise during 1990-1999 and more than expected during 2000-2009.
-Surface-ocean pCO 2 data constrain interannual variations in sea-air CO 2 fluxes 20 better than atmospheric CO 2 data, at least for the dominating variations in the tropical Pacific. The pCO 2 -based estimates can be used as ocean prior in atmospheric CO 2 inversions to improve their land flux estimates.
The presented gridded sea-air CO 2 flux estimates can be downloaded in digital form from the Jena inversion web-site, http://www.bgc-jena.mpg.de/~christian.roedenbeck/ Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Atlantic sink for atmospheric CO 2 , Science, 326, 1391Science, 326, -1393Science, 326, , doi:10.1126Science, 326, /science.1177394, 2009 Weiss, R.: Carbon dioxide in water and seawater: the solubility of a non-ideal gas, Mar. Chem., 2, 203-205, 1974. 3171 Wolter, K. and Timlin, M.: Monitoring ENSO in COADS with a seasonally adjusted principal 5 component index, in: Proceedings of the 17th Climate Diagnostics Workshop, Norman, OK, NOAA/NMC/CAC, NSSL, Oklahoma Clim. Survey, CIMMS and the School of Meteor., Univ. of Oklahoma, 52-57, 1993. 3176, 3183, 3204 Yu, L. and Weller, R. A.: Objectively Analyzed air-sea heat Fluxes (OAFlux) for the global oceans, B. Am. Meteorol. Soc., 88, 527-539, 2007. 3198 BGD 11,2014 Interannual sea-air CO 2 flux variations C. Rödenbeck et al.  Fig. 1. Illustration of the inverse procedure in our main case (run SFC): Boxes denote the process representations causally linking quantities from bottom to top: Mixed-layer carbon budget (CB), Carbonate chemistry (CC), Solubility and sea-air gas exchange (GE). The double arrow symbolizes the matching between observed and modelled CO 2 partial pressure, as gauged by the cost function J p . The thin arrow indicates the adjustments of the unknown oceaninternal carbon sources/sinks done to minimize the model-data mismatch (see Rödenbeck et al. (2013) for full details).  . The double arrow symbolizes the matching between observed and modelled CO 2 partial pressure, as gauged by the cost function J p . The thin arrow indicates the adjustments of the unknown ocean-internal carbon sources/sinks done to minimize the model-data mismatch (see Rödenbeck et al., 2013, for Fig. 3): Top: Sea-air CO 2 flux anomalies (f CO 2 ma ), middle: surface CO 2 partial pressure anomalies (p CO 2 ), bottom: anomalies in the ocean-internal DIC sources and sinks (f DIC int , note vertical scale different from f CO 2 ma ). The dark grey lines give the a-priori state of the diagnostic scheme only responding to changes in the driving variables (predominantly temperature).   Top: Interannual anomalies of the sea-air flux of atmospheric potential oxygen (APO) in the tropical Pacific, inferred from SOCAT (run SFC, blue) or independently estimated from atmospheric O2/N2 ratios and CO2 mixing ratios (atmospheric APO inversion, Rödenbeck et al., 2008, updated, orange). The background shading indicates El Niño (MEI index). Middle: Over selected years for illustration, sea-air fluxes of CO2 and O2 composing the APO flux (Eq. (5)). Bottom: Ocean-internal sources and sinks of DIC and O2. Vertical axes of all panels span the same range on a molar basis.  Fig. 9. Top: Interannual anomalies of the total ocean CO2 flux estimated by run SFC (blue) and an atmospheric CO2 inversion (magenta). Bottom: Total land CO2 flux estimated by an atmospheric CO2 inversion (magenta -standard Jena inversion (s90 v3.5) with constant priors and ocean flux IAV being adjusted from the atmospheric data, dashed magenta/blue -atmospheric inversion using run SFC as fixed ocean prior). Note threefold vertical range in the land flux panel relative to ocean panel. Fig. 7. Top: interannual anomalies of the sea-air flux of atmospheric potential oxygen (APO) in the tropical Pacific, inferred from SOCAT (run SFC, blue) or independently estimated from atmospheric O 2 /N 2 ratios and CO 2 mixing ratios (atmospheric APO inversion, Rödenbeck et al., 2008, updated, orange). The background shading indicates El Niño (MEI index). Middle: over selected years for illustration, sea-air fluxes of CO 2 and O 2 composing the APO flux (Eq. 5). Bottom: ocean-internal sources and sinks of DIC and O 2 . Vertical axes of all panels span the same range on a molar basis.  Fig. 9. Top: Interannual anomalies of the total ocean CO 2 flux estimated by run SFC (blue) and an atmospheric CO 2 inversion (magenta). Bottom: Total land CO 2 flux estimated by an atmospheric CO 2 inversion (magenta -standard Jena inversion (s90 v3.5) with constant priors and ocean flux IAV being adjusted from the atmospheric data, dashed magenta/blue -atmospheric inversion using run SFC as fixed ocean prior). Note threefold vertical range in the land flux panel relative to ocean panel. Fig. 9. Top: interannual anomalies of the total ocean CO 2 flux estimated by run SFC (blue) and an atmospheric CO 2 inversion (magenta). Bottom: total land CO 2 flux estimated by an atmospheric CO 2 inversion (magenta -standard Jena inversion (s90_v3.5) with constant priors and ocean flux IAV being adjusted from the atmospheric data, dashed magenta/blue -atmospheric inversion using run SFC as fixed ocean prior). Note threefold vertical range in the land flux panel relative to ocean panel.