2.1 CMS-Flux inversion system

The CMS-Flux framework is summarized in Fig. 1. The center of the system is the CMS-Flux inversion system, which optimizes NBE and air–sea net carbon exchanges with a 4D-Var inversion system (Liu et al., 2014). In the current system, we assume no uncertainty in fossil fuel emissions, which is a widely adopted assumption in global flux inversion systems (e.g., Crowell et al., 2019), since the uncertainty in fossil fuel emissions at regional scales is substantially less than the NBE uncertainties. The 4D-Var minimizes a cost function that includes two terms:

The first term measures the differences between the optimized fluxes and the prior fluxes normalized by the prior flux error covariance B. The second term measures the differences between observations (y) and the corresponding model simulations (h(x)) normalized by the observation error covariance R. The term h(⋅) is the observation operator that calculates observation-equivalent model-simulated X${}_{{\mathrm{CO}}_{\mathrm{2}}}$. The 4D-Var uses the adjoint (i.e., the backward integration of the transport model) (Henze et al., 2007) of the GEOS-Chem transport model to calculate the sensitivity of the observations to surface fluxes. The configurations of the inversion system are summarized in Table 1. We run both the forward model and its adjoint at 4^∘ × 5^∘ spatial resolution and optimize monthly NBE and air–sea carbon fluxes at each grid point from January 2010 to December 2018. Inputs for the system include prior carbon fluxes, meteorological drivers, and the satellite X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ (Fig. 1). Sect. 2.2 (Table 2) describes the prior flux and its uncertainties, and Sect. 2.3 (Table 3) describes the observations and the corresponding uncertainties.

Table 2Description of the prior fluxes and assumed uncertainties in the inversion system.

Download Print Version | Download XLSX

Table 3Description of observation and evaluation dataset. Data sources are listed in Table 7.

Download Print Version | Download XLSX

2.2 The prior CO₂ fluxes and uncertainties

The prior CO₂ fluxes include NBE, air–sea carbon exchange, and fossil fuel emissions (see Table 2). The data sources for the prior fluxes are listed in Table 7 and provided in the gridded fluxes. Methods to generate prior ocean carbon fluxes and fossil fuel emissions are documented in Brix et al. (2015), Caroll et al. (2020), and Oda et al. (2018). The focus of this dataset is optimized terrestrial biosphere fluxes, so we briefly describe the prior terrestrial biosphere fluxes and their uncertainties.

We construct the NBE prior using the CARDAMOM framework (Bloom et al., 2016). The CARDAMOM data assimilation system explicitly represents the time-resolved uncertainties in the NBE. The prior estimates are already constrained with multiple data streams accounting for measurement uncertainties following a Bayesian approach similar to that used in the 4D-variational approach. We use the CARDAMOM setup as described by Bloom et al. (2016, 2020) resolved at monthly timescales; data constraints include Global Ozone Monitoring Experiment 2 (GOME-2) solar-induced fluorescence (Joiner et al., 2013), MODIS leaf area index (LAI), and biomass and soil carbon (details on the data assimilation are provided in Bloom et al., 2020). In addition, mean GPP and fire carbon emissions from 2010–2017 are constrained by FLUXCOM RS + METEO version 1 GPP (Tramontana et al., 2016; Jung et al., 2017) and GFEDv4.1s (Randerson et al., 2018), respectively, both assimilated with an uncertainty of 20 %. We use the Olsen and Randerson (2004) approach to downscale monthly GPP and respiration fluxes to 3-hourly timescales, based on ERA-Interim re-analysis of global radiation and surface temperature. Fire fluxes are downscaled using the GFEDv4.1 daily- and diurnal-scale factors on monthly emissions (Giglio et al., 2013). Posterior CARDAMOM NBE estimates are then summarized as NBE mean and standard deviation values.

The NBE from CARDAMOM shows net carbon uptake of 2.3 Gt C/year over the tropics and close to neutral in the extra-tropics (Fig. B1). The year-to-year variability (i.e., interannual variability, IAV) estimated from CARDAMOM is generally less than 0.1 g C/m²/d outside of the tropics (Fig. B1). Because of the weak interannual variability estimated by CARDAMOM, we use the same 2017 NBE prior for 2018.

CARDAMOM generates uncertainty along with the mean state. The relative uncertainty over the tropics is generally larger than 100 %, and the magnitude is between 50 % and 100 % over the extra-tropics (Fig. B2). We assume no correlation in the prior flux errors in either space or time. The temporal and spatial error correlation estimates can in principle be computed by CARDAMOM. We anticipate incorporating these error correlations in subsequent versions of this dataset.

2.3 Column CO₂ observations from GOSAT and OCO-2

We use the satellite-column CO₂ retrievals from the Atmospheric Carbon Observations from Space (ACOS) team for both GOSAT (version 7.3) and OCO-2 (version 9) (Table 3). The use of the same retrieval algorithm and validation strategy adopted by the ACOS team to process both GOSAT and OCO-2 spectra maximizes the consistency between these two datasets. Both GOSAT and OCO-2 satellites carry high-resolution spectrometers optimized to return high-precision measurements of reflected sunlight within CO₂ and O₂ absorption bands in the shortwave infrared (Crisp et al., 2012). Both satellites fly in a sun-synchronous orbit. GOSAT has a 13:00±0.15 h local passing time and a 3 d ground track repeat cycle. The footprint of GOSAT is ∼10.5 km in diameter in sun nadir view (Crisp et al., 2012). The daily number of soundings processed by the ACOS-GOSAT retrieval algorithm is between a few hundred to ∼2000. Further quality control and filtering reduce the ACOS-GOSAT X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ retrievals to ∼100–300 daily (Fig. B5 in Liu et al., 2017). We only assimilate ACOS-GOSAT land nadir observations flagged as being good quality, which are the retrievals with quality flag equal to zero.

OCO-2 has a 13:30 local passing time and 16 d ground track repeat cycle. The nominal footprints of the OCO-2 are 1.25 km wide and ∼2.4 km along the orbit. Because of its small footprints and sampling strategy, OCO-2 has many more X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ retrievals than ACOS-GOSAT. To reduce the sampling error due to the resolution differences between the transport model and OCO-2 observations, we generate super observations by aggregating the observations within ∼100 km (along the same orbit) (Liu et al., 2017). The super-obing strategy was first proposed in numerical weather prediction (NWP) to assimilate dense observations (Lorenc, 1981) and is still broadly used in NWP (e.g., Liu and Rabier, 2003). More detailed information about OCO-2 super observations can be found in Liu et al. (2017). OCO-2 has four observing modes: land nadir, land glint, ocean glint, and target. Following Liu et al. (2017), we only use land nadir observations. The super observations have more uniform spatial coverage and are more comparable to the spatial representation of ACOS-GOSAT observations and the transport model (see Fig. B5 in Liu et al., 2017).

We directly use observational uncertainty provided with ACOS-GOSAT b7.3 to represent the observation error statistics, R, in Eq. (1). The uncertainty of the OCO-2 super observations is the sum of the variability of X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ used to generate each individual super observation and the mean uncertainty provided in the original OCO-2 retrievals. Kulawik et al. (2019) showed that both OCO-2 and ACOS-GOSAT bias-corrected retrievals have a mean bias of −0.1 ppm when compared with X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ from the Total Carbon Column Observing Network (TCCON) (Wunch et al., 2011), indicating consistency between ACOS-GOSAT and OCO-2 retrievals. O'Dell et al. (2018) showed that the OCO-2 X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ land nadir retrievals have an RMSE of ∼1.1 ppm when compared to TCCON retrievals; the differences between OCO-2 X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ retrievals and surface-CO₂-constrained model simulations are well within 1.0 ppm over most of the locations in the Northern Hemisphere (NH), where most of the surface CO₂ observations are located.

The magnitude of observation errors used in R is generally above 1.0 ppm, larger than the sum of random error and biases in the observations. The ACOS-GOSAT b7.3 observations from July 2009–June 2015 are used to optimize fluxes between 2010 and 2014, and the OCO-2 X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ observations from September 2014–June 2019 are used to optimize fluxes between 2015 and 2018.

The observational coverage of ACOS-GOSAT and OCO-2 is spatiotemporally dependent, with more coverage during summer than winter over the NH and more observations over midlatitudes than over the tropics (Fig. B3). The variability (i.e., standard deviation) of the annual total number of observations from 2010–2014 is within 4 % of the annual mean number for ACOS-GOSAT. Except for a data gap in 2017 caused by a malfunction of the OCO-2 instrument, the variability of the annual total number of observations between 2015 and 2018 is within 8 % of the annual mean number for OCO-2.

2.4 Uncertainty quantification

The posterior flux error covariance is the inverse Hessian, which incorporates the transport, measurement, and background errors at the 4D-Var solution (Eq. 13 in Bowman et al., 2017). Posterior flux uncertainty projected to regions can be estimated analytically based on the methods described in Fisher and Courtier (1995) and Meirink et al. (2008), using either flux singular vectors or flux increments obtained during the iterative optimization (e.g., Niwa and Fujii, 2020). In this study, we rely on a Monte Carlo approach to quantify posterior flux uncertainties following Chevallier et al. (2010) and Liu et al. (2014), which is simpler and widely used. In this approach, an ensemble of flux inversions is carried out with an ensemble of priors and simulated observations to sample the uncertainties of prior fluxes (i.e., B in Eq. 1) and observations (R in Eq. 1), respectively. The magnitude of posterior flux uncertainties is a function of assumed uncertainties in prior fluxes and observations, as well as the density of observations. Since the densities of GOSAT and OCO-2 observations are stable (Sect. 2.3) within their respective data record, we characterize the posterior flux uncertainties for 2010 and 2015 only and assume the flux uncertainties for 2011–2014 are the same as 2010 and flux uncertainties for 2016–2018 are the same as 2015.

2.5 Evaluation of posterior fluxes

Direct NBE estimates from flux towers only provide a spatial representation of roughly 1–3 km (Running et al., 1999), which is not appropriate to evaluate regional NBE from top-down flux inversions. Thus, we use two methods to indirectly evaluate the posterior NBE and its uncertainties. One is to compare annual NBE anomalies and seasonal cycle to a gross primary production (GPP) product. The other is to compare posterior CO₂ mole fractions to independent (i.e., not assimilated in the inversion) aircraft and the NOAA MBL reference observations. The second method has been broadly used to indirectly evaluate posterior fluxes from top-down flux inversions (e.g., Stephens et al., 2007; Liu and Bowman, 2016; Chevallier et al., 2019; Crowell et al., 2019). In addition to these two methods, we also compare the NBE seasonal cycles to three publicly available top-down NBE estimates that are constrained by surface CO₂ observations (Tables 3 and 7).

2.5.1 Evaluation against the independent gross primary production (GPP) product

NBE is a small residual difference between two large terms: total ecosystem respiration (TER) and GPP, plus fire. A positive NBE anomaly (i.e., less uptake from the atmosphere) has been shown to correspond to reduced GPP caused by climate anomalies (e.g., Bastos et al., 2018), and the magnitude of net uptake is proportional to GPP in most biomes observed by flux tower observations (e.g., Falk et al., 2008). Since NBE is related not only to GPP, the comparison to GPP only serves as a qualitative measure of the NBE quality. For example, we would expect the posterior NBE seasonality to be anti-correlated with GPP in the temperate and high latitudes. In this study, we use FLUXSAT GPP (Joiner et al., 2018), which is an upscaled GPP product based on flux tower GPP observations and satellite-based geometry-adjusted reflectance from the Moderate Resolution Imaging Spectroradiometer (MODIS) and solar-induced chlorophyll fluorescence observations from GOME-2 (Joiner et al., 2013). Joiner et al. (2018) show that the agreement between FLUXSAT GPP and GPP from flux towers is better than other available upscaled GPP products.

Figure 2The spatial and temporal distributions of aircraft observations used in evaluation of posterior NBE. (a) The total number of aircraft observations between 1–5 km between 2010–2018 at each 4^∘ × 5^∘ grid point. The rectangle boxes show the range of the nine sub-regions. (b) The total number of monthly aircraft observations at each longitude as a function of time.

2.5.2 Evaluation against aircraft and the NOAA marine boundary layer (MBL) reference CO₂ observations

The aircraft observations used in this study include those published in OCO-2 MIP ObsPack in August 2019 (CarbonTracker team, 2019), which include regular vertical profiles from flask samples collected on light aircraft by NOAA (Sweeney et al., 2015) and other laboratories, regular (2- to 4-weekly) vertical profiles from the Instituto de Pesquisas Espaciais (INPE) over tropical South America (SA) (Gatti et al., 2014), and vertical profiles from the Atmospheric Tomography (ATom, Wofsy et al., 2018), HIAPER Pole-to-Pole Observations (HIPPO, Wofsy et al., 2011), the ${\mathrm{O}}_{\mathrm{2}}/{\mathrm{N}}_{\mathrm{2}}$ Ratio and CO₂ airborne Southern Ocean Study (ORCAS) (Stephens et al., 2017), and Atmospheric Carbon and Transport – America (ACT-America, Davis et al., 2018) aircraft campaigns (Table 3). Figure 2 shows the aircraft observation coverage and density between 2010 and 2018. Most of the aircraft observations are concentrated over NA. ATom had four (1–4) campaigns between August 2016 and May 2018, spanning four seasons over the Pacific and Atlantic Ocean. HIPPO had five (HIPPO 1–5) campaigns over the Pacific, but only HIPPO 3–5 occurred between 2010 and 2011. HIPPO 1–2 occurred in 2009. Based on the spatial distribution of aircraft observations, we divide the comparison into nine regions: Alaska, midlatitude NA, Europe, East Asia, South Asia, Africa, Australia, Southern Ocean, and South America (Table 4 and Fig. 2).

Table 4Latitude and longitude ranges for seven sub-regions.

Download Print Version | Download XLSX

We calculate several quantities to evaluate the posterior fluxes and their uncertainty with aircraft observations. One is the monthly mean differences between posterior and aircraft CO₂ mole fractions. The second is the monthly root mean square errors (RMSEs) over each of the nine sub-regions, which is defined as

$$\begin{array}{}\text{(2)}& \text{RMSE}={\left({\displaystyle \frac{\mathrm{1}}{n}}\sum _{i=\mathrm{1}}^n(y_{\mathrm{aircraft}}^{\mathrm{o}}-y_{\mathrm{aircraft}}^b{)}_i^{\mathrm{2}}\right)}^{\frac{\mathrm{1}}{\mathrm{2}}},\end{array}$$

where $y_{\mathrm{aircraft}}^{\mathrm{o}}$ is the ith aircraft observation, $y_{\mathrm{aircraft}}^b$ is the corresponding posterior CO₂ mole fraction sampled at the ith aircraft location, and n is the number of aircraft observations over each region. The RMSE is computed over the n aircraft observations within one of the nine sub-regions. The mean differences indicate the magnitude of the mean posterior CO₂ bias, while the RMSE includes both random and systematic errors in posterior CO₂. The bias and RMSE could be due to errors in posterior fluxes, transport, and initial CO₂ concentrations. When errors in transport and initial CO₂ concentrations are smaller than the errors in the posterior fluxes, the magnitude of biases and RMSE indicates the accuracy of the posterior fluxes.

To evaluate the magnitude of posterior flux uncertainty estimates, we compare RMSE against the standard deviation of ensemble simulated aircraft observations (Eq. 3) from the Monte Carlo method (RMSE_MC). The quantity RMSE_MC can be written as

$$\begin{array}{}\text{(3)}& {\text{RMSE}}_{\mathrm{MC}}={\left[{\displaystyle \frac{\mathrm{1}}{n\mathrm{ens}}}\sum _{i\mathrm{ens}=\mathrm{1}}^{n\mathrm{ens}}\left(\right(y_{\mathrm{aircraft}}^{b\left(\text{MC}\right)}{)}_{i\mathrm{ens}}-{\overline{y}}_{\mathrm{aircraft}}^{b\left(\mathrm{MC}\right)}{)}^{\mathrm{2}}\right]}^{\frac{\mathrm{1}}{\mathrm{2}}}.\end{array}$$

The variable $(y_{\mathrm{aircraft}}^{b\left(\mathrm{MC}\right)}{)}_{i\mathrm{ens}}$ is the ith ensemble member of simulated aircraft observations from the Monte Carlo ensemble simulations, ${\overline{y}}_{\mathrm{aircraft}}^{b\left(\mathrm{MC}\right)}$ is the mean, and nens is the total number of ensemble members. For simplicity, in Eq. (3), we drop the indices for the aircraft observations used in Eq. (2). In the absence of errors in transport and initial CO₂ concentrations, when the estimated posterior flux uncertainty reflects the true posterior flux uncertainty, we show in the Appendix that

$$\begin{array}{}\text{(4)}& {\text{RMSE}}^{\mathrm{2}}={\displaystyle \frac{\mathrm{1}}{n}}\sum _{i=\mathrm{1}}^n{R}_{i,i}+{\text{RMSE}}_{\mathrm{MC}}^{\mathrm{2}},\end{array}$$

where R_aircraft is the aircraft observation error variance, which could be neglected on a regional scale.

We further calculate the ratio r between RMSE and RMSE_MC:

$$\begin{array}{}\text{(5)}& r={\displaystyle \frac{\text{RMSE}}{{\text{RMSE}}_{\mathrm{MC}}}}.\end{array}$$

A ratio close to one indicates that the posterior flux uncertainty reflects the true uncertainty in the posterior fluxes when the transport errors are small.

The presence of transport errors will make the comparison between RMSE and RMSE_MC potentially difficult to interpret. Even when RMSE_MC represents the actual uncertainty in posterior fluxes, the RMSE could be larger than RMSE_MC, since the differences between aircraft observations and model-simulated posterior mole fraction RMSE could be due to errors in both transport and the posterior fluxes, while RMSE_MC only reflects the impact of posterior flux uncertainty on simulated aircraft observations. In this study, we assume the primary sources of RMSE come from errors in posterior fluxes.

The RMSE and RMSE_MC comparison only shows differences in CO₂ space. We further calculate the sensitivity of the RMSE to the posterior flux using the GEOS-Chem adjoint. We first define a cost function J as

$$\begin{array}{}\text{(6)}& J={\text{RMSE}}^{\mathrm{2}}.\end{array}$$

The sensitivity of the mean-square error to a flux, x, at location i and month j is

$$\begin{array}{}\text{(7)}& w_{i,j}={\displaystyle \frac{\partial J}{\partial x_{i,j}}}\times x_{i,j}.\end{array}$$

This sensitivity is normalized by the flux magnitude. Equation (7) can be interpreted as the sensitivity of the RMSE² to a fractional change in the fluxes. We can estimate the time-integrated magnitude of the sensitivity over the entire assimilation window by calculating

$$\begin{array}{}\text{(8)}& S_i={\displaystyle \frac{\sum _{j=\mathrm{1}}^M\left|w_{i,j}\right|}{\sum _{k=\mathrm{1}}^P\sum _{j=\mathrm{1}}^M\left|w_{k,j}\right|}},\end{array}$$

where P is the total number of grid points and M is the total number of months from the time of the aircraft data to the beginning of the inversion. The numerator of Eq. (8) quantifies the absolute total sensitivity of the RMSE² to the fluxes at the ith grid. Normalized by the total absolute sensitivity across the globe, the quantity S_i indicates the relative sensitivity of RMSE² to fluxes at the ith grid point. Note that S_i is unitless, and it only quantifies sensitivity, not the contribution of fluxes at each grid to RMSE².

We use the NOAA MBL reference dataset (Table 7) to evaluate the CO₂ seasonal cycle over four latitude bands: 90–60^∘ N, 60–20^∘ N, 20^∘ N–20^∘ S, and 20–90^∘ S. The MBL reference is based on a subset of sites from the NOAA Cooperative Global Air Sampling Network. Only measurements that are representative of a large of volume air over a broad region are considered. In the comparison, we first remove the global mean CO₂ (https://www.esrl.noaa.gov/gmd/ccgg/trends/global.html, last access: 19 May 2020) from both the NOAA MBL reference and the posterior CO₂.

Figure 3Three types of regional masks used in calculating regional fluxes. (a) The mask is based on a combination of seven condensed MODIS IGBP plant functional types, TransCom 3 regions (Gurney et al., 2004), and continents. (b) The mask is based on latitude and continents. (c) The TransCom region mask.

2.6 Regional masks

We provide posterior NBE from 2010–2018 using three sets of regional masks (Fig. 3), in addition to the gridded product. The regional mask in Fig. 3a is based on a combination of seven plant function types condensed from MODIS IGBP and the TransCom 3 regions (Gurney et al., 2004), which is referred to as Region Mask 1 (RM1) in the following description. There are 28 regions in Fig. 3a: six in NA, four in SA, five in Eurasia (north of 40^∘ N), three in tropical Asia, three in Australia, and seven in Africa. The regional mask in Fig. 3b is based on latitude and continents with 13 regions in total, which is referred to as Region Mask 2 (RM2) in later description. Figure 3c is the TransCom regional mask with 11 regions on land.

4.1 Global fluxes

The annual atmospheric CO₂ growth rate, which is the net difference between fossil fuel emissions and total annual sink over land and ocean, is well observed by the NOAA surface CO₂ observing network (https://www.esrl.noaa.gov/gmd/ccgg/ggrn.php, last access: 12 March 2020). We compare the global total flux estimates constrained by GOSAT and OCO-2 with the NOAA CO₂ growth rate from 2010–2018 and discuss the mean carbon sink over land and ocean. Over these 9 years, the satellite-constrained atmospheric CO₂ growth rate agrees with the NOAA observed CO₂ growth rate within the uncertainty of the posterior fluxes (Fig. 4). The mean annual global surface CO₂ fluxes (in Gt C/yr) are derived from the NOAA observed CO₂ growth rate (in ppm/yr) using a conversion factor of 2.124 Gt C/ppm (Le Quéré et al., 2018). The estimated growth rate has the largest discrepancy with the NOAA observed growth rate in 2014, which may be due to a failure of one of the two solar paddles of GOSAT in May 2014 (Kuze et al., 2016). Over the 9 years, the estimated total accumulated carbon in the atmosphere is 41.5±2.4 Gt C, which is slightly lower than the accumulated carbon based on the NOAA CO₂ growth rate (45.2±0.4 Gt C). On average, we estimate that the NBE is 2.0±0.7 Gt C, $\sim \mathrm{20}\pm \mathrm{8}$ % of fossil fuel emissions, and the ocean sink is 3.0±0.1 Gt C, $\sim \mathrm{30}\pm \mathrm{1}$ % of fossil fuel emissions (Fig. 4). These numbers are within the ranges of the corresponding GCB estimates from Freidlingstein et al. (2019) (referred to as GCB-2019 hereafter). The mean NBE and ocean sink from GCB-2019 are 2.0±1.0 Gt C and 2.5±0.5 Gt C, respectively, which are 21±10 % and 26±5 % of fossil fuel emissions, respectively, between 2010–2018. The GCB does not report NBE directly, we calculate NBE from GCB-2019 as the residual differences between fossil fuel emissions, ocean net carbon sink, and atmospheric CO₂ growth rate. It is also equivalent to ($S_{\mathrm{LAND}}+B_{\mathrm{IM}}-E_{\mathrm{LUC}}$) reported by Freidlingstein et al. (2019), where S_LAND is terrestrial sink, B_IM is a budget imbalance, and E_LUC is land use change. Over these 9 years, we estimate that NBE ranges from 3.6 Gt C (∼37 % of fossil fuel emissions) in 2011 (a La Niña year) to only 0.5 Gt C (∼5 % of fossil fuel emissions) in 2015 (an El Niño year), consistent with 3.3 Gt C (35 % of fossil fuel) in 2011 to 0.9 Gt C (7 % of fossil fuel) in 2015 estimated from GCB-2019. We estimate that the ocean sinks range from 3.5 Gt C in 2015 to 2.3 Gt C in 2012, which is larger than the estimated ocean flux ranges of 2.7 in 2016 to 2.5 in 2012 reported by Freidlingstein et al. (2019).

Figure 5Mean annual regional NBE (a, b, c), uncertainty (d, e, f), and variability between 2010–2018 (g, h, i) with the three types of regional masks (Fig. 3). The first column uses a region mask based on plant functional types (PFTs) and continents (RM1). The second column uses a region mask based on latitude and continents (RM2), and the third column uses a TransCom mask.

4.2 Mean regional fluxes and uncertainties

Figure 5 shows the 9-year mean regional annual fluxes, uncertainty, and its variability between 2010–2018. Table 6 shows an example of the dataset corresponding to Fig. 5a, d, and g. It shows that large net carbon uptake occurs over Eurasia, NA, and the Southern Hemisphere (SH) midlatitudes. The largest net carbon uptake is over the eastern US ($-\mathrm{0.4}\pm \mathrm{0.1}$ Gt C (1σ uncertainty)) and high-latitude Eurasia ($-\mathrm{0.5}\pm \mathrm{0.1}$ Gt C) (Fig. 5a, b). We estimate a net land carbon sink of 2.5±0.3 Gt C/yr between 2010–2013 over the NH mid- to high latitudes, which agrees with 2.4±0.6 Gt C estimates over the same time periods based on a two-box model (Ciais et al., 2019). Net uptake in the tropics ranges from close to neutral in tropical South America (0.1±0.1 Gt C) to a net source in northern Africa (0.6±0.2 Gt C) (Fig. 5a, b). The tropics exhibit both large uncertainty and large variability. The NBE interannual variabilities over northern Africa and tropical SA are 0.5 and 0.3 Gt C, respectively, which are larger than the 0.2 and 0.1 Gt C uncertainty (Fig. 5d, e). We also find collocation of regions with large NBE and FLUXSAT GPP interannual variability (Fig. B4). The availability of flux estimates over the broadly used TransCom regions makes it easy to compare to previous studies. For example, we estimate that the annual net carbon uptake over North America is 0.7±0.1 Gt C/yr with 0.2 Gt C variability between 2010 and 2018, which agrees with 0.7±0.5 Gt C/yr estimates based on surface CO₂ observations between 1996–2007 (Peylin et al., 2013).

Figure 6The NBE interannual variability over the globe (black), the tropics (20^∘ S–20^∘ N), SH midlatitudes (60–20^∘ S), and NH midlatitudes (20–90^∘ N). For reference, the residual net land carbon sink from GCB-2019 (Friedlingstein et al., 2019) and its uncertainty is also shown (magenta).

Download

Figure 7The NBE interannual variability over six selected regions. Blue: annual NBE anomaly and its uncertainties. Green: annual GPP anomaly based on FLUXSAT.

Download

4.3 Interannual variabilities and uncertainties

Here we present hemispheric and regional NBE interannual variabilities and corresponding uncertainties (Figs. 6 and 7, and corresponding tabular data files). In Fig. 6, we further divide the globe into three large latitude bands: tropics (20^∘ S–20^∘ N), NH extra-tropics (20–85^∘ N), and SH extra-tropics (60–20^∘ S). The tropical NBE contributes 90 % to the global NBE interannual variability (IAV). The IAV of NBE over the extra-tropics is only about one-third of that over the tropics. The dominant role of tropical NBE in the global IAV of NBE agrees with Fig. 4 in Sellers et al. (2018). The top-down global annual NBE anomaly is within the 1.0 Gt C/yr uncertainty of residual NBE (i.e., fossil fuel–atmospheric growth–ocean sink) calculated from GCB-2019 (Friedlinston et al., 2019) (Fig. 6).

Figure 7 shows the annual NBE anomalies and uncertainties over a few selected regions based on RM1. Positive NBE indicates reduced net uptake relative to the 2010–2018 mean, and vice versa. Also shown in Fig. 7 are GPP anomalies estimated from FLUXSAT. Positive GPP indicates increased productivity, and vice versa. GPP drives NBE in years where anomalies are inversely correlated (e.g., positive NBE and negative GPP), and TER drives NBE in years where anomalies of GPP and NBE have the same sign or are weakly correlated. Over tropical SA evergreen broadleaf forest, the largest positive NBE anomalies occur during the 2015–2016 El Niño, corresponding to large reductions in productivity, consistent with Liu et al. (2017). In 2017, the region sees increased net uptake and increased productivity, implying a recovery from the 2015–2016 El Niño event. The variability in GPP explains 80 % of NBE variability over this region over the 9-year period. In Australian shrubland, our inversion captures the increased net uptake in 2010 and 2011 due to increased precipitation (Poulter et al., 2014) and increased productivity. The variability in GPP explains 70 % of the interannual variability in NBE. Over tropical South America savanna, the NBE interannual variability also shows strong negative correlations with GPP, with GPP explaining 40 % of NBE interannual variability. Over the midlatitude regions where the IAV is small, the R² between GPP and NBE is also small (0.0–0.5) as expected. But the increased net uptake generally corresponds to increased productivity. We also do not expect perfect negative correlation between NBE anomalies and GPP anomalies, as discussed in Sect. 2.5. The comparison between NBE and GPP provides insight into when and where net fluxes are likely dominated by productivity.

Table 5List of the data products.

n/a: not applicable.

Download Print Version | Download XLSX

Table 6The 9-year mean regional annual fluxes, uncertainties, and variability. Regions are based on the mask shown in Fig. 5a (Figure5.csv). Unit: Gt C/yr.

Download Print Version | Download XLSX

4.4 Seasonal cycle

We provide the regional mean NBE seasonal cycle, its variability, and uncertainty based on the three regional masks (Table 5). Here we briefly describe the characteristics of the NBE seasonal cycle over the 11 TransCom regions and its comparison to three independent top-down inversion results based on surface CO₂, which are CT-Europe (e.g., van der Laan-Luijkx et al., 2017), CAMS (Chevallier et al., 2005), and Jena CarbonScope (Rödenbeck et al., 2003). CMS-Flux NBE differs the most from surface-CO₂-based inversions over the South American tropical, northern Africa, tropical Asia, and NH boreal regions. The CMS-Flux NBE has a larger seasonal cycle amplitude over tropical Asia and northern Africa, where the surface CO₂ constraint is weak, while it has a smaller seasonal cycle amplitude over the boreal region; this may be due to the sparse satellite observations over the high latitudes and weaker seasonal amplitude of the prior CARDAMOM fluxes. The comparison to FLUXSAT GPP can only qualitatively evaluate the NBE seasonal cycle but cannot differentiate among different estimates. In general, the months that have larger productivity correspond to months with a net uptake of carbon from the atmosphere, especially over the NH (Fig. 8). More research is still needed to understand the seasonal cycles of NBE, including its phase (i.e., transition from source to sink) and amplitude (peak-to-trough difference), as well as its relationships with GPP and respiration.

Figure 8The NBE climatological seasonality over TransCom regions. The seasonal cycle is calculated over 2010–2017 since CT-Europe only covers until 2017. Black: CMS-Flux NBE and its uncertainty; blue shaded: mean NBE seasonality based on surface CO₂ inversion results from CAMS, CT-Europe, and Jena CarbonScope; red: CAMS; magenta: CT-Europe; green: Jena CarbonScope. The names of each region are shown on individual subplots.

Download

Figure 9Comparison between posterior CO₂ mole fraction and aircraft observations. Left panel: detrended posterior CO₂ (y axis) vs. detrended aircraft CO₂ (x axis) over nine regions. The dashed line is the one-to-one line; right panel: the differences between posterior CO₂ and aircraft CO₂ as a function of time (black), RMSE (blue; unit: ppm), and RMSE_MC (red).

Download

5.1 Comparison to aircraft observations over nine sub-regions

In this section, we evaluate posterior CO₂ against aircraft observations over the nine sub-regions listed in Table 4 and Fig. 2. We compare the posterior CO₂ to aircraft CO₂ mole fractions above the planetary boundary layer and up to the mid-troposphere (1–5 km) at the locations and time of aircraft observations, and then we calculate the monthly mean error statistics between 1–5 km. The aircraft observations between 1–5 km are more sensitive to regional fluxes (Liu et al., 2015; Liu and Bowman, 2016). Scatter plots in the left column of Fig. 9 show regional monthly mean detrended aircraft CO₂ observations (x axis) versus the simulated detrended posterior CO₂ (y axis). We used the NOAA global CO₂ trend to detrend both the observations and model-simulated mole fractions (ftp://aftp.cmdl.noaa.gov/products/trends/co2/co2_trend_gl.txt, last access: 19 May 2020). Over the NH regions (Fig. 9a, b, c, d) and Africa (Fig. 9f), the R² is greater than or equal to 0.9, which indicates that the posterior CO₂ captures the observed seasonality. The low R² (0.7) value in South Asia is caused by one outlier. Over the Southern Ocean, Australia, and SA, the R² is between 0.2 and 0.4, reflecting weaker CO₂ seasonality over these regions and possible bias in ocean flux estimates (see discussions later).

The right panel of Fig. 9 shows the monthly mean differences between posterior CO₂ and aircraft observations (black), RMSE (Eq. 2) (blue line), and RMSE_MC (Eq. 3) (red line). The magnitude of the mean differences between the posterior CO₂ and aircraft observations is less than 0.5 ppm, except over the Southern Ocean, which has a −0.8 ppm bias. The mean differences between posterior CO₂ and aircraft observations are primarily caused by errors in transport and biases in assimilated satellite observations, while RMSE_MC is the internal flux error projected into mole fraction space. With the exception of the Southern Ocean, for all regions mean bias is significantly less than RMSE_MC, which suggests that transport and data bias in satellite observations may be much smaller than the internal flux errors. Note that RMSE_MC is smaller than RMSE over the first ∼ 6 months of simulation, which may indicate a dominant impact of errors in transport and initial CO₂ concentration on posterior CO₂ RMSE.

As demonstrated in Sect. 2.5, comparing RMSE and RMSE_MC is a test of the accuracy of posterior flux uncertainty estimate. Over all the regions, the differences between RMSE and RMSE_MC are smaller than 0.3 ppm, which indicates a comparable magnitude between empirical posterior flux uncertainty estimates from the Monte Carlo method and the actual posterior flux uncertainty over the regions that these aircraft observations are sensitive to. These aircraft observations are sensitive to NBE over a broad region as shown in Fig. B5. Note that Figs. B5 and B8–B10 are calculated using Eq. (8).

Figure 10Left column: the mean differences between posterior CO₂ and aircraft observations from ATom 1–4 aircraft campaigns between 1–5 km (a–d). Middle column: the root mean square errors (RMSEs) between aircraft observations and posterior CO₂ between 1–5 km. The color bar is the same as in the left column. Right column: the ratio between RMSE and RMSE_MC based on ensemble CO₂ from the Monte Carlo uncertainty estimation method.

Figure 11Left column: the mean differences between posterior CO₂ and aircraft observations from HIPPO 3–5 aircraft campaigns between 1–5 km (a–c) (unit: ppm). The time frame of each campaign is in the figure. Middle column: the root mean square errors (RMSEs) between aircraft observations and posterior CO₂ between 1–5 km (unit: ppm). The color bar is the same as in the left column. Right column: the ratio between RMSE and RMSE_MC based on ensemble CO₂ from the Monte Carlo method.

5.2 Comparison to aircraft observations from ATom and HIPPO aircraft campaigns

Figures 10 and 11 show comparisons to aircraft CO₂ from the ATom 1–4 campaigns spanning four seasons and HIPPO 3–5 over the Pacific Ocean between 1–5 km. The vertical curtain comparisons are shown in Figs. B6 and B7. The mean differences between posterior CO₂ and aircraft CO₂ are quite uniform (within 0.5 ppm) throughout the column, except over the Atlantic Ocean during ATom 1–2 and the Southern Ocean during ATom 1 (Figs. S6 and S7). Also shown in Figs. 10 and 11 are RMSEs of each aircraft campaign (middle column) and the ratio between RMSE and RMSE_MC (right column). A ratio larger than one between RMSE and RMSE_MC indicates errors in either transport or underestimation of the posterior flux uncertainty (Sect. 2.5).

Over most of the flight tracks during ATom 1–4, the posterior CO₂ errors are between −0.5 and 0.5 ppm, the RMSE is smaller than 0.5 ppm, and the ratio between RMSE and RMSE_MC is smaller than or equal to 1. However, off the coast of Africa during ATom 1 and 2 and over the Southern Ocean during ATom 1, the mean differences between posterior CO₂ and aircraft observations are larger than 0.5 ppm. During ATom 1 (29 July–23 August 2016), the mean differences between posterior CO₂ and aircraft CO₂ show large negative biases, while during ATom 2 (26 January–21 February 2017) it has large positive biases off the coast of Africa. The ratio between RMSE and RMSE_MC is significantly larger than one over these regions, which indicates an underestimation of posterior flux uncertainty or a large magnitude of transport errors during that time period.

We further run adjoint sensitivity analyses over the three regions with ratios significantly larger than one to identify the posterior fluxes that could contribute to the large differences between posterior CO₂ and aircraft observations during ATom 1–2. We run the adjoint model backward for 3 months from the observation time and calculate S_i as defined in Eq. (7). The adjoint sensitivity analysis indicates that the large mismatch between aircraft observations and model simulations during ATom 1 and -2 off the coast of Africa could be potentially driven by errors in posterior fluxes over tropical Africa (Fig. B8). The large posterior CO₂ errors and large ratio between RMSE and RMSE_MC over the Southern Ocean during ATom 1 are driven by flux errors in oceanic fluxes around 30^∘ S and over Australia (Fig. B9), which also contribute to the large errors in comparison to aircraft observations over the Southern Ocean shown in Fig. 9h.

During the HIPPO aircraft campaigns, the absolute errors in posterior CO₂ across the Pacific are less than 0.5 ppm, except over the Arctic Ocean and over Alaska in summer (Fig. 11), consistent with Fig. 10a. The large errors over the Arctic Ocean may be related to both transport errors and the accuracy of high-latitude fluxes. Byrne et al. (2020) provide a brief summary of the challenges in simulating CO₂ over high latitudes using a transport model with 4^∘ × 5^∘ resolution. Increasing the resolution of the transport model may reduce transport errors over high latitudes.

We run adjoint sensitivity analyses over the high-latitude regions where the differences between posterior CO₂ and aircraft observations are large (Fig. 11). The adjoint sensitivity analysis (Fig. B10) shows that the large errors over these regions could be driven by errors in fluxes over Alaska as well as broad NH midlatitude regions.

Figure 12Comparison between posterior CO₂ and the NOAA marine boundary layer (MBL) reference sites. Black: observations averaged over each latitude band; blue and shaded area: posterior CO₂ and its uncertainty. The global mean CO₂ (https://www.esrl.noaa.gov/gmd/ccgg/trends/global.html, last access: 29 May 2020) was subtracted from both the NOAA MBL reference and posterior CO₂ before the comparison.

Download

5.3 Comparison to MBL reference sites

Since MBL reference sites sample air over broad regions, the comparison to detrended MBL observations indirectly evaluates the NBE over large regions. Figure 12 shows the comparison over four latitude bands. The uncertainty of posterior CO₂ concentration is from the MC method. Except over 90–20^∘ S, the differences between observations and posterior CO₂ are within posterior CO₂ uncertainty estimates. The posterior CO₂ concentrations have the smallest bias and random errors over the tropical latitude band. The R² is above 0.9 over NH mid- to high latitudes, consistent with Fig. 9. Over 90–20^∘ S, the posterior CO₂ has positive bias in 2013 and 2014 and negative bias and much weaker seasonality between January 2015–December 2018 compared to observations, which indicates possible biases in Southern Ocean flux estimates (Fig. B11). The low bias over the Southern Ocean is consistent with aircraft comparison during the OCO-2 period (Figs. 9–10, B9). The changes of performance after 2013 over 90–20^∘ S are most likely due to the prior ocean carbon fluxes. Evaluation of ocean carbon fluxes is out of the scope of this study. Note that since we only assimilate land nadir X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ observations in this study due to known issues with the OCO-2 v9 ocean glint observations (O'Dell et all., 2018), the constraint of top-down inversion on air–sea CO₂ exchanges is weak (not shown). The ocean glint observations of OCO-2 v10 observations have been improved compared to v9 (Osterman et al., 2020). We expect to have a better estimate of ocean carbon fluxes over the Southern Ocean when assimilating both land and ocean X${}_{{\mathrm{CO}}_{\mathrm{2}}}$ observations from GOSAT and OCO-2 in the future.