Reconstruction of the history of anthropogenic CO2 concentrations in the ocean

The release of fossil fuel CO2 to the atmosphere by human activity has been implicated as the predominant cause of recent global climate change. The ocean plays a crucial role in mitigating the effects of this perturbation to the climate system, sequestering 20 to 35 per cent of anthropogenic CO2 emissions. Although much progress has been made in recent years in understanding and quantifying this sink, considerable uncertainties remain as to the distribution of anthropogenic CO2 in the ocean, its rate of uptake over the industrial era, and the relative roles of the ocean and terrestrial biosphere in anthropogenic CO2 sequestration. Here we address these questions by presenting an observationally based reconstruction of the spatially resolved, time-dependent history of anthropogenic carbon in the ocean over the industrial era. Our approach is based on the recognition that the transport of tracers in the ocean can be described by a Green’s function, which we estimate from tracer data using a maximum entropy deconvolution technique. Our results indicate that ocean uptake of anthropogenic CO2 has increased sharply since the 1950s, with a small decline in the rate of increase in the last few decades. We estimate the inventory and uptake rate of anthropogenic CO2 in 2008 at 140 ± 25 Pg C and 2.3 ± 0.6 Pg C yr-1, respectively. We find that the Southern Ocean is the primary conduit by which this CO2 enters the ocean (contributing over 40 per cent of the anthropogenic CO2 inventory in the ocean in 2008). Our results also suggest that the terrestrial biosphere was a source of CO2 until the 1940s, subsequently turning into a sink. Taken over the entire industrial period, and accounting for uncertainties, we estimate that the terrestrial biosphere has been anywhere from neutral to a net source of CO2, contributing up to half as much CO2 as has been taken up by the ocean over the same period.

The release of fossil fuel CO 2 to the atmosphere by human activity has been implicated as the predominant cause of recent global climate change 1 . The ocean plays a crucial role in mitigating the effects of this perturbation to the climate system, sequestering 20 to 35 per cent of anthropogenic CO 2 emissions [2][3][4] . Although much progress has been made in recent years in understanding and quantifying this sink, considerable uncertainties remain as to the distribution of anthropogenic CO 2 in the ocean, its rate of uptake over the industrial era, and the relative roles of the ocean and terrestrial biosphere in anthropogenic CO 2 sequestration. Here we address these questions by presenting an observationally based reconstruction of the spatially resolved, time-dependent history of anthropogenic carbon in the ocean over the industrial era. Our approach is based on the recognition that the transport of tracers in the ocean can be described by a Green's function, which we estimate from tracer data using a maximum entropy deconvolution technique. Our results indicate that ocean uptake of anthropogenic CO 2 has increased sharply since the 1950s, with a small decline in the rate of increase in the last few decades. We estimate the inventory and uptake rate of anthropogenic CO 2 in 2008 at 140 6 25 Pg C and 2.3 6 0.6 Pg C yr 21 , respectively. We find that the Southern Ocean is the primary conduit by which this CO 2 enters the ocean (contributing over 40 per cent of the anthropogenic CO 2 inventory in the ocean in 2008). Our results also suggest that the terrestrial biosphere was a source of CO 2 until the 1940s, subsequently turning into a sink. Taken over the entire industrial period, and accounting for uncertainties, we estimate that the terrestrial biosphere has been anywhere from neutral to a net source of CO 2 , contributing up to half as much CO 2 as has been taken up by the ocean over the same period.
A key challenge for estimating anthropogenic CO 2 (C ant ) in the ocean is that C ant is not a directly measurable quantity. Existing estimates of C ant are thus based on indirect techniques, such as socalled 'back calculation' methods that attempt to separate the small anthropogenic perturbation (of the order of a few per cent) from the large background distribution of carbon by correcting the measured total dissolved inorganic carbon (DIC) concentration for changes due to biological activity and air-sea disequilibrium 5,6 . The recent availability of a high quality, global tracer data set 7 and significant improvements in methodology, notably the development of the DC* approach 8 , made it possible to apply these ideas, and led to one of the first observationally based global estimates of the distribution of C ant in the ocean 2 . Although a major advance in our understanding of C ant in the ocean, it has been suggested that this estimate suffers from a number of biases and limitations [9][10][11] , one of them being that it provides only a snapshot for the mid-1990s.
Our method for estimating anthropogenic CO 2 builds on previous work 11, 12 , but extends it in several significant ways. Following previous studies, we exploit the fact that the anthropogenic perturbation can be treated as a conservative tracer 13 transported by ocean circulation from the surface mixed layer into the interior. The transport of water that carries C ant into the interior ocean involves considerable dispersion and mixing of water masses of different ages and of different surface origin. Anthropogenic CO 2 at an interior location x and at time t is therefore related to its history in the surface mixed layer, C s ant , through a convolution equation involving a kernel, G, which partitions each water parcel according to where and when it was last in contact with the sea surface: The time integral in equation (1) is over the industrial era and the space integral is over the whole ocean surface. The kernel, G x, t; x', t' ð Þ , is an intrinsic property of the ocean circulation and not of any particular tracer. As such, it can be used to propagate the surface boundary condition of any conservative tracer into the full three-dimensional concentration field. We exploit this fact by using a suite of well sampled oceanic tracers such as chlorofluorocarbons, natural 14 C, temperature, and salinity from the GLODAP and World Ocean Atlas databases (see Methods), to provide constraints analogous to equation (1) from which we deconvolve G. We recognize that there is no single tracer that perfectly emulates the C ant transient, necessitating the use of multiple tracers with distinct time histories to constrain G. To regularize the under-determined deconvolution problem, we use a maximum entropy method 14 which is well suited for problems with positive kernels. We note that the C ant estimated via equation (1) is relatively insensitive to errors in G as the latter only appears via a convolution 12 .
Apart from G, we need an estimate of C s ant in order to apply equation (1) to compute C ant . We obtain this boundary condition from the known atmospheric CO 2 history by requiring that the rate of change of the inventory of C ant (obtained by integrating equation (1) over the volume of the ocean) must, by mass conservation, be equal to its net flux into the ocean. The air-sea flux of C ant is proportional to the change in surface disequilibrium of CO 2 (see Methods). To make further progress, we exploit the empirical result from ocean carbon cycle models that the change in disequilibrium is, to a very good approximation, proportional to the (known) anthropogenic perturbation in the atmospheric partial pressure of CO 2 , p CO2 . To estimate the unknown, spatially variable proportionality constant, we combine the above constraint with the requirement that our solution match observed p CO2 values averaged over a discrete set of surface patches, subject to the CO 2 -system equilibrium chemistry. Once the proportionality constants are known, the history of C ant on each surface patch is readily obtained.
We note that our inversion method provides improvements to reduce the three main biases of most previous techniques 10 , namely: (1) the air-sea disequilibrium is allowed to evolve and is estimated by solving a nonlinear optimization problem, (2) the mixing of waters of different ages is accounted for by using multiple transient and steady tracers to constrain the age distribution, and (3) the mixing of different end-member water types is accounted for by constraining G using multiple steady and transient tracers. Tests using tracers simulated with a global ocean circulation model show that our method is able to successfully recover the complex spatial distribution of C ant in the model with a maximum error in the column inventory of ,2 mol m 22 and an error in the global inventory of less than 2% (Supplementary Information).
Using our new inverse method, we have reconstructed the first three-dimensional, time-varying, history of anthropogenic carbon in the ocean from AD 1765 to AD 2008. Figure 1a shows the column inventory of C ant in 2008. The total inventory in that year was ,140 6 25 Pg C. This estimate excludes the Arctic Ocean and marginal seas not covered by the GLODAP database. Using a recent estimate based on CFC-11 for the former 15 , and an area scaling approach for the latter 2 , would increase our estimate of the global inventory by ,11 Pg C. We have also partitioned this inventory according to where at the surface the anthropogenic CO 2 penetrated the ocean. As is evident from Fig. 1b, the high latitude oceans, driven by intermediate and deep water formation, constitute the most intense sinks of C ant . In particular, the Southern Ocean is by far the largest conduit by which anthropogenic CO 2 enters the ocean: roughly 40% of the C ant residing in the ocean in 2008 entered the ocean south of 40u S.
It is useful to compare our result for 1994 to previous estimates that are available for that year. The inventory estimated using the socalled DC* method 2 for 1994 is 106 6 21 Pg C. The 1994 inventory estimated using the transit-time distribution (TTD) method 11 is 107 Pg C, with a range of 94-121 Pg C. Both estimates are consistent with our 1994 estimate of 114 6 22 Pg C. However, the previous TTD-based estimate includes a 20% downward correction for the fact that air-sea disequilibrium was incorrectly treated as being constant. Such a correction, based on model simulations, is not necessary for our estimate because our inverse method explicitly accounts for changing air-sea disequilibrium. The spatial distribution we obtain is quite different from that obtained with the DC* method, particularly in the Southern Ocean. Relative to the DC*-based estimate, our estimate of C ant , like the TTD-based estimate, is generally lower in the upper ocean and higher in the deep ocean 11 . Two key aspects of the DC* method are the use of a single tracer age to characterize transport and the assumption of constant disequilibrium. These give rise to competing biases 10 that largely cancel out, leading to the close, but fortuitous, agreement between our estimate of the total inventory and that derived using the DC* method. Figure 2 shows the uptake history over the industrial era (AD 1765 to AD 2008) computed from the time-varying inventory. (The corresponding space-and time-varying change in surface disequilibrium of CO 2 driving this uptake is also estimated by our inversion method.) There has been a sharp increase in ocean uptake since the 1950s in response to a higher growth rate of atmospheric CO 2 , although the rate of increase has decreased somewhat in the last few decades. Our estimated uptake rate for the 1990s, 2.0 6 0.6 Pg C yr 21 , agrees well with the IPCC consensus estimate based on independent methods 4 ( Fig. 2 and Table 1).
Our work has important implications for the terrestrial carbon budget, here computed as a residual between the main fossil fuel source and the ocean and atmosphere sinks. (Including the highly  The shaded area represents the error envelope (see Fig. 1 legend). Also shown are the decadal average uptake rates adopted by the IPCC fourthassessment report (AR4) 4 (blue circles; vertical error bars are 61 s.d. and horizontal error bars span the averaging period of years) and the atmospheric CO 2 mixing ratio 29 used for the inversion (red dashed line).

NATURE | Vol 462 | 19 November 2009
LETTERS uncertain source due to changes in land use 16 provides a different perspective; Supplementary Fig. 3.) There remains considerable uncertainty regarding the relative partitioning of anthropogenic emissions between the ocean and land biosphere 3 . Our time-evolving estimate of the ocean uptake provides a more precise and detailed view of the land sink. Table 1 compares our estimates for the uptake rate of C ant by the land biosphere with other estimates 4,17 . There is generally good agreement between the two, although the latter only go back to the 1980s, whereas our approach covers the entire industrial period. Figure 3 shows the evolution of the various sources and sinks of anthropogenic CO 2 between AD 1765 and AD 2005. Our results indicate that the terrestrial biosphere was a source of C ant until the 1940s, roughly in line with previous model-based estimates 3,18 , after which it turned into a sink of anthropogenic CO 2 . Taken over the entire industrial period, and accounting for uncertainties, we estimate that the terrestrial biosphere has been anywhere from neutral to a net source of CO 2 , contributing up to half as much C ant as has been taken up by the ocean over the same period.
One potential source of error we have neglected is variability in ocean circulation, especially long-term trends. Such variability, however, remains poorly constrained, and its impact on anthropogenic CO 2 uptake is still debated. For example, ocean models show a slight decline in the rate-of-growth of C ant uptake by the Southern Ocean over the past few decades due to an increase in the strength of the meridional overturning circulation (MOC) in response to strengthened westerly winds 19,20 . However, a recent study 21 finds no observational evidence for such a change in the MOC, suggesting that the simulated changes in the MOC and uptake may be an artefact of the eddy-parametrization used in coarse-resolution ocean models 21,22 . Nevertheless, it is useful to place these modelled changes within the context of our inverse calculations, which assume a cyclo-stationary circulation. The simulated decrease 19 , of ,0.08 Pg C yr 21 per decade between 1981 and 2004, would imply a reduction in C ant uptake of roughly 2.1 Pg C over that period (equivalent to a 1 p.p.m. increase in atmospheric CO 2 concentration). This value should be compared with a total ocean uptake of 46 6 6 Pg C over the same period, as estimated by our inverse method (17.5 6 2.5 Pg C for the Southern Ocean)-that is, a 5% reduction in the global ocean sink of anthropogenic CO 2 , but still well within the uncertainty of our purely data-based estimate.
To conclude, we have presented an observationally based estimate of the time-evolving distribution of anthropogenic CO 2 in the ocean over the industrial era. Unlike other recent inverse calculations 23 , we do not rely on an ocean model for transport information. Errors in transport simulated by such models [23][24][25] are a large source of uncertainty in estimates of C ant based on those inverse calculations. Instead, we constrain the transport by using tracer observations. Our results can thus be used to assess and improve global ocean biogeochemical models, initialize high-resolution simulations, and provide boundary conditions for atmospheric inversions of anthropogenic sources and sinks.

METHODS SUMMARY
The combined GLODAP/WOA05 databases (see Methods) provide six constraints, one for each tracer, of the form: where C is the observed tracer concentration, and l its radioactive decay rate (non-zero only for 14 C). To reduce the indeterminacy, we assume that the ocean circulation is stationary except for a cyclo-stationary seasonal cycle and we discretize the sea surface position variable, x9, into a discrete set of 26 surface patches (Fig. 1b). With these assumptions, the kernel for each interior position, x, can be written as a discrete function, G i, k, m ð Þin which i is the surface-patch index, and k and m are respectively the number of years since, and the month when, the water was last in the surface mixed layer. To further regularize the deconvolution, we used a maximum entropy approach in which we maximize an entropy functional, J G ½ ~{ Þ , subject to the discretized tracer constraints. M is a prior estimate of G, which we take to be an analytical solution to the one-dimensional advection-diffusion equation known as the inverse Gaussian 26 . Simulations in ocean general circulation models 27 show that the inverse Gaussian form captures well the general characteristics of G. The inverse Gaussian is characterized by two parameters, a mean age C and width D. Consistent with tracer observations 26,28 , we set D 5 C, which leaves us with one free parameter. To specify C, we make it a function of depth, increasing linearly from 10 yr in the surface layer to 2,000 yr in the deepest layer. The resulting variational problem is solved using the method of Lagrange multipliers, and yields a solution of the form (shown for clarity using the continuous time variables): where a j (x) is the Lagrange multiplier required to enforce the jth observational constraint. The above yields a system of nonlinear algebraic equations for the Lagrange multipliers, which we solve using standard methods.
Full Methods and any associated references are available in the online version of the paper at www.nature.com/nature.    30 (including a small contribution from cement production) is the only source considered here, and is shown as positive values. Sinks, shown as negative values, include the atmosphere, ocean, and land biosphere. Error envelope (as in Fig. 1 legend), indicated by broken lines and the shaded area, includes a 5% uncertainty in fossil fuel emissions 17 .